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

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

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

Public Member Functions

 OpLhs (int tag, BlockData &block_data, CommonData &common_data, const double beta_gc=0)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_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

MatrixDouble mat
 
const double betaGC
 
- 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 599 of file GriffithForceElement.hpp.

Constructor & Destructor Documentation

◆ OpLhs()

FractureMechanics::GriffithForceElement::OpLhs::OpLhs ( int  tag,
BlockData block_data,
CommonData common_data,
const double  beta_gc = 0 
)
inline

Definition at line 602 of file GriffithForceElement.hpp.

604 : FaceElementForcesAndSourcesCore::UserDataOperator(
605 "MESH_NODE_POSITIONS", "MESH_NODE_POSITIONS",
606 UserDataOperator::OPROWCOL),
607 betaGC(beta_gc), AuxOp(tag, block_data, common_data) {}
AuxOp(int tag, BlockData &block_data, CommonData &common_data)

Member Function Documentation

◆ doWork()

MoFEMErrorCode FractureMechanics::GriffithForceElement::OpLhs::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData &  row_data,
DataForcesAndSourcesCore::EntData &  col_data 
)
inline

Definition at line 612 of file GriffithForceElement.hpp.

615 {
617
618 if (row_type != MBVERTEX) {
620 }
621
629
630 AuxFunctions<double> auxFun;
631 auxFun.currentCoords.resize(9, false);
632
633 int nb_dofs = row_data.getFieldData().size();
634 for (int dd = 0; dd != nb_dofs; dd++)
635 auxFun.currentCoords[dd] = row_data.getFieldData()[dd];
636
637 CHKERR setIndices(row_data);
638 CHKERR setVariables(this, row_data);
639 int nb_gauss_pts = row_data.getN().size1();
640 int nb_base_row = row_data.getFieldData().size() / 3;
641 int nb_base_col = col_data.getFieldData().size() / 3;
642 int row_nb_dofs = row_data.getIndices().size();
643
644 mat.resize(9, 9, false);
645 mat.clear();
646 auto calculate_derivative = [&](FTensor::Tensor1<double *, 2> &t_diff) {
648 t_d_n(i, j) = FTensor::levi_civita(i, j, k) * auxFun.t_Tangent1_ksi(k) *
649 t_diff(N1);
650 t_d_n(i, j) -= FTensor::levi_civita(i, j, k) *
651 auxFun.t_Tangent2_eta(k) * t_diff(N0);
652
653 return t_d_n;
654 };
655
656 CHKERR auxFun.calculateAreaAndNormal(col_data.getDiffN());
657 auto t_normal = auxFun.t_Normal;
658 const double coefficient_1 = 0.5 * pow((t_normal(i) * t_normal(i)), -0.5);
659 const double coefficient_2 = 0.5 * pow((t_normal(i) * t_normal(i)), -1.5);
660
661 // for homogeneous case, this is a bit ugly
662 auto &rho_v = commonData.densityRho;
663 if (rho_v.empty() || rho_v.size() != nb_base_row) {
664 rho_v.resize(nb_base_row, false);
665 }
666 auto rho = getFTensor0FromVec(rho_v);
667 for (int gg = 0; gg != nb_gauss_pts; gg++) {
668
669 double val = getGaussPts()(2, gg) * 0.5;
670 auto t_row_diff = row_data.getFTensor1DiffN<2>(gg, 0);
671 for (int rrr = 0; rrr != nb_base_row; rrr++) {
672
674 &mat(3 * rrr + 0, 0), &mat(3 * rrr + 0, 1), &mat(3 * rrr + 0, 2),
675 &mat(3 * rrr + 1, 0), &mat(3 * rrr + 1, 1), &mat(3 * rrr + 1, 2),
676 &mat(3 * rrr + 2, 0), &mat(3 * rrr + 2, 1), &mat(3 * rrr + 2, 2),
677 3);
678
679 auto t_tan_row = calculate_derivative(t_row_diff);
680 auto t_col_diff = col_data.getFTensor1DiffN<2>(gg, 0);
681 const double coef = blockData.gc * val * pow(rho, betaGC);
682 for (int ccc = 0; ccc != nb_base_col; ccc++) {
683
684 auto t_tan_col = calculate_derivative(t_col_diff);
685
686 t_mat(i, j) -= coefficient_1 *
687 (t_tan_row(i, l) * t_tan_col(l, j) +
688 FTensor::levi_civita(i, j, k) * t_col_diff(N0) *
689 t_row_diff(N1) * t_normal(k) -
690 FTensor::levi_civita(i, j, k) * t_col_diff(N1) *
691 t_row_diff(N0) * t_normal(k));
692
693 t_mat(i, j) += coefficient_2 * (t_tan_row(i, k) * t_normal(k)) *
694 (t_tan_col(l, j) * t_normal(l));
695
696 t_mat(i, j) *= coef;
697 ++t_col_diff;
698 ++t_mat;
699 }
700 ++t_row_diff;
701 ++rho;
702 }
703 }
704
705 int col_nb_dofs = col_data.getIndices().size();
706
707 CHKERR MatSetValues(getFEMethod()->snes_B, row_nb_dofs,
708 &*rowIndices.data().begin(), col_nb_dofs,
709 &*col_data.getIndices().data().begin(),
710 &*mat.data().begin(), ADD_VALUES);
711
713 }
static Number< 1 > N1
static Number< 0 > N0
#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
FTensor::Index< 'm', SPACE_DIM > m
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
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
double rho
Definition: plastic.cpp:98
MoFEMErrorCode setVariables(FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, DataForcesAndSourcesCore::EntData &data)
MoFEMErrorCode setIndices(DataForcesAndSourcesCore::EntData &data)

Member Data Documentation

◆ betaGC

const double FractureMechanics::GriffithForceElement::OpLhs::betaGC

Definition at line 610 of file GriffithForceElement.hpp.

◆ mat

MatrixDouble FractureMechanics::GriffithForceElement::OpLhs::mat

Definition at line 609 of file GriffithForceElement.hpp.


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