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

calculate griffith force vector More...

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

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

Public Member Functions

 OpCalculateGriffithForce (int tag, BlockData &block_data, CommonData &common_data)
 
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

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

calculate griffith force vector

\[ \mathbf{f}_{\mathrm{m}}^{\mathrm{h}}=\frac{1}{2}\left(\tilde{\mathbf{A}}_{\Gamma}^{\mathrm{h}}\right)^{\mathrm{T}} \mathbf{g}_{c}(\rho(\mathbf{X})) ]\f The matrix $\mathbf{A}_{\Gamma}^h$ comprising of direction vectors along the crack front that are normal to the crack front and tangent to the crack surface is defined as follows: \f[ A^h_{\Gamma} = \| \mathbf{N}(\tilde{\mathbf{X}}) \| = \left\| \epsilon_{ijk} \frac{\partial \Phi^\alpha_p}{\partial \xi_i} \frac{\partial \Phi^\beta_r}{\partial \xi_j} \tilde{X}^\alpha_p \tilde{X}^\beta_r \right\| \]

Definition at line 378 of file GriffithForceElement.hpp.

Constructor & Destructor Documentation

◆ OpCalculateGriffithForce()

FractureMechanics::GriffithForceElement::OpCalculateGriffithForce::OpCalculateGriffithForce ( int  tag,
BlockData block_data,
CommonData common_data 
)
inline

Definition at line 381 of file GriffithForceElement.hpp.

383 : FaceElementForcesAndSourcesCore::UserDataOperator(
384 "MESH_NODE_POSITIONS", UserDataOperator::OPROW),
385 AuxOp(tag, block_data, common_data) {}
AuxOp(int tag, BlockData &block_data, CommonData &common_data)

Member Function Documentation

◆ doWork()

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

Definition at line 387 of file GriffithForceElement.hpp.

388 {
390
391 if (type != MBVERTEX)
393 CHKERR setIndices(data);
394 CHKERR setVariables(this, data);
395 int nb_dofs = data.getIndices().size();
396 int nb_gauss_pts = data.getN().size1();
397
398 auxFun.referenceCoords.resize(nb_dofs, false);
399 auxFun.currentCoords.resize(nb_dofs, false);
400
401 for (int dd = 0; dd != nb_dofs; dd++) {
402 auxFun.referenceCoords[dd] = getCoords()[dd];
403 auxFun.currentCoords[dd] = data.getFieldData()[dd];
404 }
405
406 for (int gg = 0; gg != nb_gauss_pts; gg++) {
407 double val = getGaussPts()(2, gg) * 0.5;
409 data.getDiffN(gg));
410 }
411 commonData.griffithForce.resize(nb_dofs, false);
413
414 std::copy(auxFun.griffithForce.begin(), auxFun.griffithForce.end(),
415 commonData.griffithForce.begin());
416
418 }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#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
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
MoFEMErrorCode calculateGriffithForce(const double gc, const double beta, const MatrixDouble &diff_n)
MoFEMErrorCode setVariables(FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, DataForcesAndSourcesCore::EntData &data)
MoFEMErrorCode setIndices(DataForcesAndSourcesCore::EntData &data)
VectorDouble griffithForce
Griffith force at nodes.

Member Data Documentation

◆ auxFun

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

Definition at line 386 of file GriffithForceElement.hpp.


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