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

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

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

Public Member Functions

 OpConstrainsLhs (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 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

MoFEM::InterfacemField
 
SmartPetscObj< Vec > & deltaVec
 
MatrixDouble nG_dX
 
MatrixDouble nF_dLambda
 
MatrixDouble nF_dX
 
MatrixDouble jacDelta
 
ublas::vector< double * > jacDeltaPtr
 
VectorDouble d_dA
 
VectorDouble d_Delta
 
ublas::vector< adoublea_Delta
 
AuxFunctions< adoublea_auxFun
 
ublas::vector< adoublea_dA
 
ublas::vector< adoublea_lambdaAtNodes
 
- 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 1048 of file GriffithForceElement.hpp.

Constructor & Destructor Documentation

◆ OpConstrainsLhs()

FractureMechanics::GriffithForceElement::OpConstrainsLhs::OpConstrainsLhs ( 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 1055 of file GriffithForceElement.hpp.

1059 : FaceElementForcesAndSourcesCore::UserDataOperator(
1060 "MESH_NODE_POSITIONS", lambda_field_name,
1061 UserDataOperator::OPROWCOL,
1062 false // not symmetric operator
1063 ),
1064 AuxOp(tag, block_data, common_data), mField(m_field),
1065 deltaVec(delta_vec) {}
AuxOp(int tag, BlockData &block_data, CommonData &common_data)

Member Function Documentation

◆ doWork()

MoFEMErrorCode FractureMechanics::GriffithForceElement::OpConstrainsLhs::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 1081 of file GriffithForceElement.hpp.

1084 {
1086
1087 if (row_type != MBVERTEX || col_type != MBVERTEX)
1089
1090 CHKERR setIndices(row_data);
1091 CHKERR setLambdaNodes(this, colFieldName);
1092 CHKERR setLambdaIndices(this, colFieldName);
1093
1094 activeVariables.resize(18, false);
1095 a_auxFun.referenceCoords.resize(9, false);
1096 a_auxFun.currentCoords.resize(9, false);
1097
1098 a_Delta.resize(3, false);
1099 d_Delta.resize(3, false);
1100 trace_on(tAg);
1101 for (int dd = 0; dd != 9; dd++) {
1102 double val = getCoords()[dd];
1103 a_auxFun.referenceCoords[dd] <<= val;
1104 activeVariables[dd] = val;
1105 }
1106 for (int dd = 0; dd != 9; dd++) {
1107 double val = row_data.getFieldData()[dd];
1108 a_auxFun.currentCoords[dd] <<= val;
1109 activeVariables[9 + dd] = val;
1110 }
1111
1112 CHKERR a_auxFun.calculateGriffithForce(1, getGaussPts()(2, 0) * 0.5,
1113 row_data.getDiffN(0));
1114
1115 // calculate area change
1116 a_Delta.clear();
1117 for (int nn = 0; nn != 3; nn++) {
1118 if (col_data.getIndices()[nn] == -1)
1119 continue;
1120 if (col_data.getIndices()[nn] != rowLambdaIndices[nn]) {
1121 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1122 "data inconsistency");
1123 }
1124 for (int dd = 0; dd != 3; dd++) {
1125 int idx = 3 * nn + dd;
1126 if (row_data.getFieldDofs()[idx]->getEnt() !=
1127 col_data.getFieldDofs()[nn]->getEnt()) {
1128 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1129 "data inconsistency");
1130 }
1131 a_Delta[nn] +=
1132 a_auxFun.griffithForce[idx] *
1133 (a_auxFun.currentCoords[idx] - a_auxFun.referenceCoords[idx]);
1134 }
1135 }
1136 for (int nn = 0; nn != 3; nn++) {
1137 a_Delta[nn] >>= d_Delta[nn];
1138 }
1139 trace_off();
1140
1141 jacDelta.resize(3, 18, false);
1142 jacDeltaPtr.resize(3, false);
1143 for (int nn = 0; nn != 3; nn++) {
1144 jacDeltaPtr[nn] = &jacDelta(nn, 0);
1145 }
1146
1147 {
1148 int r;
1149 // play recorder for jacobian
1150 r = ::jacobian(tAg, 3, 18, &activeVariables[0], &jacDeltaPtr[0]);
1151 if (r < 3) { // function is locally analytic
1152 SETERRQ1(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
1153 "ADOL-C function evaluation with error r = %d", r);
1154 }
1155 }
1156
1157 const double *delta;
1158 Vec delta_vec_local;
1159 CHKERR VecGhostGetLocalForm(deltaVec, &delta_vec_local);
1160 int vec_size;
1161 CHKERR VecGetSize(delta_vec_local, &vec_size);
1162 CHKERR VecGetArrayRead(delta_vec_local, &delta);
1163 nG_dX.resize(3, 9, false);
1164 nG_dX.clear();
1165 for (int nn = 0; nn != 3; ++nn) {
1166 int local_idx = col_data.getLocalIndices()[nn];
1167 if (local_idx != rowLambdaIndicesLocal[nn]) {
1168 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1169 "data inconsistency");
1170 }
1171 if (local_idx == -1)
1172 continue;
1173 if (vec_size < local_idx) {
1174 int g_vec_size;
1175 CHKERR VecGetLocalSize(deltaVec, &g_vec_size);
1176 SETERRQ3(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1177 "data inconsistency %d < %d (%d)", vec_size, local_idx,
1178 g_vec_size);
1179 }
1180 double W = col_data.getFieldData()[nn] - blockData.E * delta[local_idx];
1181 double diff_constrain = -diffCalMax_a(W, 0, blockData.r);
1182 for (int dd = 0; dd != 9; dd++) {
1183 double diffW = -blockData.E * jacDelta(nn, 9 + dd);
1184 nG_dX(nn, dd) = diff_constrain * diffW;
1185 }
1186 }
1187 CHKERR VecRestoreArrayRead(delta_vec_local, &delta);
1188 CHKERR VecGhostRestoreLocalForm(deltaVec, &delta_vec_local);
1189 CHKERR MatSetValues(getFEMethod()->snes_B, 3,
1190 &*col_data.getIndices().data().begin(), 9,
1191 &*row_data.getIndices().data().begin(),
1192 &*nG_dX.data().begin(), ADD_VALUES);
1193
1194 // Dervatives of nF
1195
1196 // Linearisation with adol-C
1197 activeVariables.resize(21, false);
1198 d_dA.resize(9, false);
1199 a_dA.resize(9, false);
1200 a_lambdaAtNodes.resize(3, false);
1201 a_auxFun.referenceCoords.resize(9, false);
1202
1203 trace_on(tAg);
1204 for (int dd = 0; dd != 9; dd++) {
1205 double val = getCoords()[dd];
1206 a_auxFun.referenceCoords[dd] <<= val;
1207 activeVariables[dd] = val;
1208 }
1209 for (int dd = 0; dd != 9; dd++) {
1210 double val = row_data.getFieldData()[dd];
1211 a_auxFun.currentCoords[dd] <<= val;
1212 activeVariables[9 + dd] = val;
1213 }
1214 for (int nn = 0; nn != 3; nn++) {
1215 double val = col_data.getFieldData()[nn];
1216 a_lambdaAtNodes[nn] <<= val;
1217 activeVariables[18 + nn] = val;
1218 }
1219
1220 CHKERR a_auxFun.calculateGriffithForce(1, getGaussPts()(2, 0) * 0.5,
1221 row_data.getDiffN(0));
1222
1223 a_dA.clear();
1224 for (int nn = 0; nn != 3; nn++) {
1225 if (col_data.getIndices()[nn] == -1)
1226 continue;
1227 for (int dd = 0; dd != 3; dd++) {
1228 int idx = 3 * nn + dd;
1229 a_dA[idx] -= a_lambdaAtNodes[nn] * a_auxFun.griffithForce[idx];
1230 }
1231 }
1232 for (int dd = 0; dd != 9; dd++) {
1233 a_dA[dd] >>= d_dA[dd];
1234 }
1235 trace_off();
1236
1237 commonData.tangentGriffithForce.resize(9, 21, false);
1239 for (int dd = 0; dd != 9; ++dd) {
1242 }
1243
1244 {
1245 int r;
1246 // play recorder for jacobian
1247 r = ::jacobian(tAg, 9, 21, &activeVariables[0],
1249 if (r < 3) { // function is locally analytic
1250 SETERRQ1(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
1251 "ADOL-C function evaluation with error r = %d", r);
1252 }
1253 }
1254
1255 nF_dLambda.resize(9, 3);
1256 for (int rr = 0; rr != 9; rr++) {
1257 for (int cc = 0; cc != 3; cc++) {
1258 nF_dLambda(rr, cc) = commonData.tangentGriffithForce(rr, 18 + cc);
1259 }
1260 }
1261 CHKERR MatSetValues(getFEMethod()->snes_B, 9, &*rowIndices.data().begin(),
1262 3, &*col_data.getIndices().data().begin(),
1263 &*nF_dLambda.data().begin(), ADD_VALUES);
1264
1265 nF_dX.resize(9, 9, false);
1266 for (int rr = 0; rr != 9; rr++) {
1267 for (int cc = 0; cc != 9; cc++) {
1268 nF_dX(rr, cc) = commonData.tangentGriffithForce(rr, 9 + cc);
1269 }
1270 }
1271 CHKERR MatSetValues(getFEMethod()->snes_B, 9, &*rowIndices.data().begin(),
1272 9, &*row_data.getIndices().data().begin(),
1273 &*nF_dX.data().begin(), ADD_VALUES);
1274
1276 }
#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
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition: definitions.h:34
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#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
constexpr double W
const FTensor::Tensor2< T, Dim, Dim > Vec
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
static double diffCalMax_a(double a, double b, double r)
const double r
rate factor
static constexpr double delta
MoFEMErrorCode setLambdaNodes(FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, const std::string &lambda_field_name)
MoFEMErrorCode setLambdaIndices(FaceElementForcesAndSourcesCore::UserDataOperator *fe_ptr, const std::string &lambda_field_name)
MoFEMErrorCode setIndices(DataForcesAndSourcesCore::EntData &data)
vector< double * > tangentGriffithForceRowPtr
Pointers to rows of tangent matrix.

Member Data Documentation

◆ a_auxFun

AuxFunctions<adouble> FractureMechanics::GriffithForceElement::OpConstrainsLhs::a_auxFun

Definition at line 1077 of file GriffithForceElement.hpp.

◆ a_dA

ublas::vector<adouble> FractureMechanics::GriffithForceElement::OpConstrainsLhs::a_dA

Definition at line 1078 of file GriffithForceElement.hpp.

◆ a_Delta

ublas::vector<adouble> FractureMechanics::GriffithForceElement::OpConstrainsLhs::a_Delta

Definition at line 1076 of file GriffithForceElement.hpp.

◆ a_lambdaAtNodes

ublas::vector<adouble> FractureMechanics::GriffithForceElement::OpConstrainsLhs::a_lambdaAtNodes

Definition at line 1079 of file GriffithForceElement.hpp.

◆ d_dA

VectorDouble FractureMechanics::GriffithForceElement::OpConstrainsLhs::d_dA

Definition at line 1073 of file GriffithForceElement.hpp.

◆ d_Delta

VectorDouble FractureMechanics::GriffithForceElement::OpConstrainsLhs::d_Delta

Definition at line 1074 of file GriffithForceElement.hpp.

◆ deltaVec

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

Definition at line 1053 of file GriffithForceElement.hpp.

◆ jacDelta

MatrixDouble FractureMechanics::GriffithForceElement::OpConstrainsLhs::jacDelta

Definition at line 1071 of file GriffithForceElement.hpp.

◆ jacDeltaPtr

ublas::vector<double *> FractureMechanics::GriffithForceElement::OpConstrainsLhs::jacDeltaPtr

Definition at line 1072 of file GriffithForceElement.hpp.

◆ mField

MoFEM::Interface& FractureMechanics::GriffithForceElement::OpConstrainsLhs::mField

Definition at line 1052 of file GriffithForceElement.hpp.

◆ nF_dLambda

MatrixDouble FractureMechanics::GriffithForceElement::OpConstrainsLhs::nF_dLambda

Definition at line 1068 of file GriffithForceElement.hpp.

◆ nF_dX

MatrixDouble FractureMechanics::GriffithForceElement::OpConstrainsLhs::nF_dX

Definition at line 1069 of file GriffithForceElement.hpp.

◆ nG_dX

MatrixDouble FractureMechanics::GriffithForceElement::OpConstrainsLhs::nG_dX

Definition at line 1067 of file GriffithForceElement.hpp.


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