v0.13.1
Loading...
Searching...
No Matches
Classes | Public Member Functions | Public Attributes | List of all members
FractureMechanics::GriffithForceElement::FrontArcLengthControl Struct Reference

arc-length element More...

#include <users_modules/fracture_mechanics/src/GriffithForceElement.hpp>

Inheritance diagram for FractureMechanics::GriffithForceElement::FrontArcLengthControl:
[legend]
Collaboration diagram for FractureMechanics::GriffithForceElement::FrontArcLengthControl:
[legend]

Classes

struct  OpAssemble_db
 assemble internal residual derivative More...
 
struct  OpIntLambda
 assemble internal residual More...
 

Public Member Functions

 FrontArcLengthControl (int tag, GriffithForceElement::BlockData &block_data, GriffithForceElement::CommonData &common_data, const std::string &lambda_field_name, boost::shared_ptr< ArcLengthCtx > &arc_ptr)
 
virtual ~FrontArcLengthControl ()
 
MoFEMErrorCode preProcess ()
 
double calculateLambdaInt ()
 Calculate internal lambda. More...
 
MoFEMErrorCode operator() ()
 
MoFEMErrorCode postProcess ()
 
MoFEMErrorCode calculateDb ()
 Calculate db. More...
 
MoFEMErrorCode calculateDxAndDlambda (Vec x)
 

Public Attributes

int tAg
 
boost::shared_ptr< ArcLengthCtxarcPtr
 
BlockDatablockData
 
CommonDatacommonData
 
std::string lambdaFieldName
 
Vec ghostIntLambda
 
double intLambdaArray [2]
 
doubleintLambdaDelta
 
doubleintLambdaLambda
 
double intLambda
 
boost::shared_ptr< MyTriangleFEfeIntLambda
 

Detailed Description

arc-length element

Calualte residual for arc length method amd vector db, which is derivative of residual.

Definition at line 1286 of file GriffithForceElement.hpp.

Constructor & Destructor Documentation

◆ FrontArcLengthControl()

FractureMechanics::GriffithForceElement::FrontArcLengthControl::FrontArcLengthControl ( int  tag,
GriffithForceElement::BlockData block_data,
GriffithForceElement::CommonData common_data,
const std::string &  lambda_field_name,
boost::shared_ptr< ArcLengthCtx > &  arc_ptr 
)
inline

Definition at line 1495 of file GriffithForceElement.hpp.

1499 : tAg(tag), arcPtr(arc_ptr), blockData(block_data),
1500 commonData(common_data), lambdaFieldName(lambda_field_name),
1503 int two[] = {0, 1};
1504 int nb_local = arcPtr->mField.get_comm_rank() == 0 ? 2 : 0;
1505 int nb_ghost = nb_local ? 0 : 2;
1506 ierr = VecCreateGhostWithArray(PETSC_COMM_WORLD, nb_local, 2, nb_ghost,
1508 CHKERRABORT(PETSC_COMM_WORLD, ierr);
1509 feIntLambda = boost::shared_ptr<MyTriangleFE>(
1510 new GriffithForceElement::MyTriangleFE(arcPtr->mField));
1511 }
static PetscErrorCode ierr

◆ ~FrontArcLengthControl()

virtual FractureMechanics::GriffithForceElement::FrontArcLengthControl::~FrontArcLengthControl ( )
inlinevirtual

Definition at line 1512 of file GriffithForceElement.hpp.

1512 {
1513 ierr = VecDestroy(&ghostIntLambda);
1514 CHKERRABORT(PETSC_COMM_WORLD, ierr);
1515 }

Member Function Documentation

◆ calculateDb()

MoFEMErrorCode FractureMechanics::GriffithForceElement::FrontArcLengthControl::calculateDb ( )
inline

Calculate db.

Definition at line 1586 of file GriffithForceElement.hpp.

1586 {
1588
1589 MOFEM_LOG_CHANNEL("WORLD");
1591 MOFEM_LOG_TAG("WORLD", "Arc");
1592
1593 double alpha = arcPtr->alpha;
1594 double beta = arcPtr->beta;
1595 double d_lambda = arcPtr->dLambda;
1596
1597 CHKERR VecZeroEntries(ghostIntLambda);
1598 CHKERR VecGhostUpdateBegin(ghostIntLambda, INSERT_VALUES,
1599 SCATTER_FORWARD);
1600 CHKERR VecGhostUpdateEnd(ghostIntLambda, INSERT_VALUES, SCATTER_FORWARD);
1601 CHKERR VecZeroEntries(arcPtr->db);
1602 CHKERR VecGhostUpdateBegin(arcPtr->db, INSERT_VALUES, SCATTER_FORWARD);
1603 CHKERR VecGhostUpdateEnd(arcPtr->db, INSERT_VALUES, SCATTER_FORWARD);
1604
1605 // calculate internal lambda
1606 feIntLambda->getOpPtrVector().clear();
1607 feIntLambda->getOpPtrVector().push_back(new OpIntLambda(
1608 arcPtr->mField, blockData, commonData, lambdaFieldName, *this));
1609 feIntLambda->getOpPtrVector().push_back(new OpAssemble_db(
1610 tAg, arcPtr->mField, blockData, commonData, lambdaFieldName, *this));
1611 CHKERR arcPtr->mField.loop_finite_elements(
1612 problemPtr->getName(), "CRACKFRONT_AREA_ELEM", (*feIntLambda));
1613 const double *dx;
1614 const auto bit_field_number = getFieldBitNumber(lambdaFieldName);
1615 CHKERR VecGetArrayRead(arcPtr->dx, &dx);
1617 bit_field_number, dit)) {
1618 if (static_cast<int>(dit->get()->getPart()) !=
1619 arcPtr->mField.get_comm_rank())
1620 continue;
1621 int idx = dit->get()->getPetscLocalDofIdx();
1622 double val = dx[idx];
1623 CHKERR VecSetValue(ghostIntLambda, 1, val, ADD_VALUES);
1624 }
1625 CHKERR VecRestoreArrayRead(arcPtr->dx, &dx);
1626
1627 CHKERR VecAssemblyBegin(ghostIntLambda);
1628 CHKERR VecAssemblyEnd(ghostIntLambda);
1629 CHKERR VecGhostUpdateBegin(ghostIntLambda, ADD_VALUES, SCATTER_REVERSE);
1630 CHKERR VecGhostUpdateEnd(ghostIntLambda, ADD_VALUES, SCATTER_REVERSE);
1631 CHKERR VecGhostUpdateBegin(ghostIntLambda, INSERT_VALUES,
1632 SCATTER_FORWARD);
1633 CHKERR VecGhostUpdateEnd(ghostIntLambda, INSERT_VALUES, SCATTER_FORWARD);
1634
1635 CHKERR VecAssemblyBegin(arcPtr->db);
1636 CHKERR VecAssemblyEnd(arcPtr->db);
1637 CHKERR VecGhostUpdateBegin(arcPtr->db, ADD_VALUES, SCATTER_REVERSE);
1638 CHKERR VecGhostUpdateEnd(arcPtr->db, ADD_VALUES, SCATTER_REVERSE);
1639 CHKERR VecGhostUpdateBegin(arcPtr->db, INSERT_VALUES, SCATTER_FORWARD);
1640 CHKERR VecGhostUpdateEnd(arcPtr->db, INSERT_VALUES, SCATTER_FORWARD);
1641
1643 alpha * beta * beta * d_lambda * d_lambda * arcPtr->F_lambda2;
1644
1646 "WORLD", Sev::noisy,
1647 "\tintLambdaLambda = %9.8e intLambdaDelta = %9.8e intLambda = %9.8e",
1649
1650 double nrm;
1651 CHKERR VecNorm(arcPtr->db, NORM_2, &nrm);
1652 MOFEM_LOG_C("WORLD", Sev::noisy, "\tdb nrm = %9.8e", nrm);
1653
1655 }
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:304
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:277
#define MOFEM_LOG_FUNCTION()
Set scope.
Definition: LogManager.hpp:318
#define _IT_NUMEREDDOF_ROW_BY_BITNUMBER_FOR_LOOP_(PROBLEMPTR, FIELD_BIT_NUMBER, IT)

◆ calculateDxAndDlambda()

MoFEMErrorCode FractureMechanics::GriffithForceElement::FrontArcLengthControl::calculateDxAndDlambda ( Vec  x)
inline

Definition at line 1657 of file GriffithForceElement.hpp.

1657 {
1659 MOFEM_LOG_CHANNEL("WORLD");
1661 MOFEM_LOG_TAG("WORLD", "Arc");
1662
1663 // Calculate dx
1664 CHKERR VecGhostUpdateBegin(arcPtr->x0, INSERT_VALUES, SCATTER_FORWARD);
1665 CHKERR VecGhostUpdateEnd(arcPtr->x0, INSERT_VALUES, SCATTER_FORWARD);
1666
1667 Vec l_x, l_x0, l_dx;
1668 CHKERR VecGhostGetLocalForm(x, &l_x);
1669 CHKERR VecGhostGetLocalForm(arcPtr->x0, &l_x0);
1670 CHKERR VecGhostGetLocalForm(arcPtr->dx, &l_dx);
1671 {
1672 double *x_array, *x0_array, *dx_array;
1673 CHKERR VecGetArray(l_x, &x_array);
1674 CHKERR VecGetArray(l_x0, &x0_array);
1675 CHKERR VecGetArray(l_dx, &dx_array);
1676 int size =
1677 problemPtr->getNbLocalDofsRow() + problemPtr->getNbGhostDofsRow();
1678 for (int i = 0; i != size; ++i) {
1679 dx_array[i] = x_array[i] - x0_array[i];
1680 }
1681 CHKERR VecRestoreArray(l_x, &x_array);
1682 CHKERR VecRestoreArray(l_x0, &x0_array);
1683 CHKERR VecRestoreArray(l_dx, &dx_array);
1684 }
1685 CHKERR VecGhostRestoreLocalForm(x, &l_x);
1686 CHKERR VecGhostRestoreLocalForm(arcPtr->x0, &l_x0);
1687 CHKERR VecGhostRestoreLocalForm(arcPtr->dx, &l_dx);
1688
1689 // Calculate dlambda
1690 if (arcPtr->getPetscLocalDofIdx() != -1) {
1691 double *array;
1692 CHKERR VecGetArray(arcPtr->dx, &array);
1693 arcPtr->dLambda = array[arcPtr->getPetscLocalDofIdx()];
1694 array[arcPtr->getPetscLocalDofIdx()] = 0;
1695 CHKERR VecRestoreArray(arcPtr->dx, &array);
1696 }
1697 CHKERR VecGhostUpdateBegin(arcPtr->ghosTdLambda, INSERT_VALUES,
1698 SCATTER_FORWARD);
1699 CHKERR VecGhostUpdateEnd(arcPtr->ghosTdLambda, INSERT_VALUES,
1700 SCATTER_FORWARD);
1701
1702 // Calculate dx2
1703 double x_nrm, x0_nrm;
1704 CHKERR VecNorm(x, NORM_2, &x_nrm);
1705 CHKERR VecNorm(arcPtr->x0, NORM_2, &x0_nrm);
1706 CHKERR VecDot(arcPtr->dx, arcPtr->dx, &arcPtr->dx2);
1707 MOFEM_LOG_C("WORLD", Sev::noisy,
1708 "\tx norm = %6.4e x0 norm = %6.4e dx2 = %6.4e", x_nrm, x0_nrm,
1709 arcPtr->dx2);
1711 }
FTensor::Index< 'i', SPACE_DIM > i
const FTensor::Tensor2< T, Dim, Dim > Vec

◆ calculateLambdaInt()

double FractureMechanics::GriffithForceElement::FrontArcLengthControl::calculateLambdaInt ( )
inline

Calculate internal lambda.

Definition at line 1534 of file GriffithForceElement.hpp.

1534 {
1536 return intLambda;
1538 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440

◆ operator()()

MoFEMErrorCode FractureMechanics::GriffithForceElement::FrontArcLengthControl::operator() ( )
inline

Definition at line 1540 of file GriffithForceElement.hpp.

1540 {
1542 MOFEM_LOG_CHANNEL("WORLD");
1544 MOFEM_LOG_TAG("WORLD", "Arc");
1545
1546 switch (snes_ctx) {
1547 case CTX_SNESSETFUNCTION: {
1548 arcPtr->res_lambda = calculateLambdaInt() - arcPtr->alpha * arcPtr->s;
1549 CHKERR VecSetValue(snes_f, arcPtr->getPetscGlobalDofIdx(),
1550 arcPtr->res_lambda, ADD_VALUES);
1551 } break;
1552 case CTX_SNESSETJACOBIAN: {
1553 arcPtr->dIag = 2 * arcPtr->dLambda * arcPtr->alpha *
1554 pow(arcPtr->beta, 2) * arcPtr->F_lambda2;
1555 MOFEM_LOG_C("WORLD", Sev::noisy, "diag = %9.8e", arcPtr->dIag);
1556 CHKERR MatSetValue(snes_B, arcPtr->getPetscGlobalDofIdx(),
1557 arcPtr->getPetscGlobalDofIdx(), 1, ADD_VALUES);
1558 } break;
1559 default:
1560 break;
1561 }
1563 }

◆ postProcess()

MoFEMErrorCode FractureMechanics::GriffithForceElement::FrontArcLengthControl::postProcess ( )
inline

Definition at line 1565 of file GriffithForceElement.hpp.

1565 {
1567 switch (snes_ctx) {
1568 case CTX_SNESSETFUNCTION: {
1569 CHKERR VecAssemblyBegin(snes_f);
1570 CHKERR VecAssemblyEnd(snes_f);
1571 } break;
1572 case CTX_SNESSETJACOBIAN: {
1573 CHKERR VecGhostUpdateBegin(arcPtr->ghostDiag, INSERT_VALUES,
1574 SCATTER_FORWARD);
1575 CHKERR VecGhostUpdateEnd(arcPtr->ghostDiag, INSERT_VALUES,
1576 SCATTER_FORWARD);
1577 } break;
1578 default:
1579 break;
1580 }
1582 }

◆ preProcess()

MoFEMErrorCode FractureMechanics::GriffithForceElement::FrontArcLengthControl::preProcess ( )
inline

Definition at line 1517 of file GriffithForceElement.hpp.

1517 {
1519 switch (snes_ctx) {
1520 case CTX_SNESSETFUNCTION: {
1521 CHKERR VecAssemblyBegin(snes_f);
1522 CHKERR VecAssemblyEnd(snes_f);
1525 } break;
1526 default:
1527 break;
1528 }
1530 }

Member Data Documentation

◆ arcPtr

boost::shared_ptr<ArcLengthCtx> FractureMechanics::GriffithForceElement::FrontArcLengthControl::arcPtr

Definition at line 1289 of file GriffithForceElement.hpp.

◆ blockData

BlockData& FractureMechanics::GriffithForceElement::FrontArcLengthControl::blockData

Definition at line 1290 of file GriffithForceElement.hpp.

◆ commonData

CommonData& FractureMechanics::GriffithForceElement::FrontArcLengthControl::commonData

Definition at line 1291 of file GriffithForceElement.hpp.

◆ feIntLambda

boost::shared_ptr<MyTriangleFE> FractureMechanics::GriffithForceElement::FrontArcLengthControl::feIntLambda

Definition at line 1300 of file GriffithForceElement.hpp.

◆ ghostIntLambda

Vec FractureMechanics::GriffithForceElement::FrontArcLengthControl::ghostIntLambda

Definition at line 1294 of file GriffithForceElement.hpp.

◆ intLambda

double FractureMechanics::GriffithForceElement::FrontArcLengthControl::intLambda

Definition at line 1298 of file GriffithForceElement.hpp.

◆ intLambdaArray

double FractureMechanics::GriffithForceElement::FrontArcLengthControl::intLambdaArray[2]

Definition at line 1295 of file GriffithForceElement.hpp.

◆ intLambdaDelta

double& FractureMechanics::GriffithForceElement::FrontArcLengthControl::intLambdaDelta

Definition at line 1296 of file GriffithForceElement.hpp.

◆ intLambdaLambda

double& FractureMechanics::GriffithForceElement::FrontArcLengthControl::intLambdaLambda

Definition at line 1297 of file GriffithForceElement.hpp.

◆ lambdaFieldName

std::string FractureMechanics::GriffithForceElement::FrontArcLengthControl::lambdaFieldName

Definition at line 1292 of file GriffithForceElement.hpp.

◆ tAg

int FractureMechanics::GriffithForceElement::FrontArcLengthControl::tAg

Definition at line 1288 of file GriffithForceElement.hpp.


The documentation for this struct was generated from the following file: