v0.14.0
Public Member Functions | Public Attributes | List of all members
FractureMechanics::GriffithForceElement::OpConstrainsDelta Struct Reference

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

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

Public Member Functions

 OpConstrainsDelta (MoFEM::Interface &m_field, int tag, BlockData &block_data, CommonData &common_data, const std::string &lambda_field_name, SmartPetscObj< Vec > &delta_vec)
 
MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 
- Public Member Functions inherited from FractureMechanics::GriffithForceElement::AuxOp
 AuxOp (int tag, BlockData &block_data, CommonData &common_data)
 
MoFEMErrorCode setIndices (DataForcesAndSourcesCore::EntData &data)
 
MoFEMErrorCode setVariables (FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, DataForcesAndSourcesCore::EntData &data)
 
MoFEMErrorCode setLambdaNodes (FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, const std::string &lambda_field_name)
 
MoFEMErrorCode setLambdaIndices (FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, const std::string &lambda_field_name)
 

Public Attributes

MoFEM::InterfacemField
 
const std::string lambdaFieldName
 
SmartPetscObj< Vec > & deltaVec
 
AuxFunctions< doubleauxFun
 
- Public Attributes inherited from FractureMechanics::GriffithForceElement::AuxOp
int tAg
 
BlockDatablockData
 
CommonDatacommonData
 
ublas::vector< int > rowIndices
 
ublas::vector< int > rowIndicesLocal
 
VectorInt3 rowLambdaIndices
 
VectorInt3 rowLambdaIndicesLocal
 
VectorDouble3 lambdaAtNodes
 
VectorDouble activeVariables
 

Detailed Description

Definition at line 854 of file GriffithForceElement.hpp.

Constructor & Destructor Documentation

◆ OpConstrainsDelta()

FractureMechanics::GriffithForceElement::OpConstrainsDelta::OpConstrainsDelta ( MoFEM::Interface m_field,
int  tag,
BlockData block_data,
CommonData common_data,
const std::string &  lambda_field_name,
SmartPetscObj< Vec > &  delta_vec 
)
inline

Definition at line 862 of file GriffithForceElement.hpp.

867  "MESH_NODE_POSITIONS", UserDataOperator::OPROW),
868  AuxOp(tag, block_data, common_data), mField(m_field),
869  lambdaFieldName(lambda_field_name), deltaVec(delta_vec) {}

Member Function Documentation

◆ doWork()

MoFEMErrorCode FractureMechanics::GriffithForceElement::OpConstrainsDelta::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData data 
)
inline

Definition at line 873 of file GriffithForceElement.hpp.

874  {
876  if (type != MBVERTEX) {
878  }
879  // get indices of material DOFs
880  CHKERR setIndices(data);
881  // get indices of Lagrange multiplier DOFs
883  // do calculations
884  auxFun.currentCoords.resize(9, false);
885  for (int dd = 0; dd != 9; ++dd)
886  auxFun.currentCoords[dd] = data.getFieldData()[dd];
887 
888  CHKERR auxFun.calculateGriffithForce(1, getGaussPts()(2, 0) * 0.5,
889  data.getDiffN(0));
890 
891  // assemble area change
892  double delta[] = {0, 0, 0};
893  for (int nn = 0; nn != 3; ++nn) {
894  if (rowLambdaIndices[nn] == -1) {
895  if (rowIndices[3 * nn + 0] != -1 || rowIndices[3 * nn + 1] != -1 ||
896  rowIndices[3 * nn + 2] != -1) {
897  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
898  "Element has no nodes at crack front");
899  }
900  continue;
901  }
902  for (int dd = 0; dd != 3; ++dd) {
903  int idx = 3 * nn + dd;
904  delta[nn] += auxFun.griffithForce[idx] *
905  (data.getFieldData()[idx] - getCoords()[idx]);
906  }
907  }
908  // add nodes to vector of lambda indices
909  CHKERR VecSetValues(deltaVec, 3, &rowLambdaIndices[0], delta, ADD_VALUES);
910 
912  }

Member Data Documentation

◆ auxFun

AuxFunctions<double> FractureMechanics::GriffithForceElement::OpConstrainsDelta::auxFun

Definition at line 871 of file GriffithForceElement.hpp.

◆ deltaVec

SmartPetscObj<Vec>& FractureMechanics::GriffithForceElement::OpConstrainsDelta::deltaVec

Definition at line 860 of file GriffithForceElement.hpp.

◆ lambdaFieldName

const std::string FractureMechanics::GriffithForceElement::OpConstrainsDelta::lambdaFieldName

Definition at line 859 of file GriffithForceElement.hpp.

◆ mField

MoFEM::Interface& FractureMechanics::GriffithForceElement::OpConstrainsDelta::mField

Definition at line 858 of file GriffithForceElement.hpp.


The documentation for this struct was generated from the following file:
FractureMechanics::GriffithForceElement::AuxFunctions::calculateGriffithForce
MoFEMErrorCode calculateGriffithForce(const double gc, const double beta, const MatrixDouble &diff_n)
Definition: GriffithForceElement.hpp:201
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
FractureMechanics::GriffithForceElement::AuxOp::AuxOp
AuxOp(int tag, BlockData &block_data, CommonData &common_data)
Definition: GriffithForceElement.hpp:221
FractureMechanics::GriffithForceElement::AuxOp::setIndices
MoFEMErrorCode setIndices(DataForcesAndSourcesCore::EntData &data)
Definition: GriffithForceElement.hpp:232
FractureMechanics::GriffithForceElement::OpConstrainsDelta::lambdaFieldName
const std::string lambdaFieldName
Definition: GriffithForceElement.hpp:859
MoFEM::VecSetValues
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
Definition: EntitiesFieldData.hpp:1589
FractureMechanics::GriffithForceElement::AuxFunctions::currentCoords
ublas::vector< TYPE > currentCoords
Definition: GriffithForceElement.hpp:126
FractureMechanics::GriffithForceElement::OpConstrainsDelta::deltaVec
SmartPetscObj< Vec > & deltaVec
Definition: GriffithForceElement.hpp:860
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
delta
static constexpr double delta
Definition: prism_elements_from_surface.cpp:18
FractureMechanics::GriffithForceElement::OpConstrainsDelta::mField
MoFEM::Interface & mField
Definition: GriffithForceElement.hpp:858
FractureMechanics::GriffithForceElement::OpConstrainsDelta::auxFun
AuxFunctions< double > auxFun
Definition: GriffithForceElement.hpp:871
convert.type
type
Definition: convert.py:64
FractureMechanics::GriffithForceElement::AuxOp::setLambdaIndices
MoFEMErrorCode setLambdaIndices(FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, const std::string &lambda_field_name)
Definition: GriffithForceElement.hpp:301
FractureMechanics::GriffithForceElement::AuxOp::rowLambdaIndices
VectorInt3 rowLambdaIndices
Definition: GriffithForceElement.hpp:227
FTensor::dd
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
FractureMechanics::GriffithForceElement::AuxFunctions::griffithForce
ublas::vector< TYPE > griffithForce
Definition: GriffithForceElement.hpp:127
UserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Definition: HookeElement.hpp:75
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
FractureMechanics::GriffithForceElement::AuxOp::rowIndices
ublas::vector< int > rowIndices
Definition: GriffithForceElement.hpp:222
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROW
@ OPROW
operator doWork function is executed on FE rows
Definition: ForcesAndSourcesCore.hpp:567