834 {
835
837
838 using CachePhi = boost::tuple<int, int, MatrixDouble>;
839
840 CGGUserPolynomialBase(boost::shared_ptr<CachePhi> cache_phi = nullptr);
841 ~CGGUserPolynomialBase() = default;
842
843 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
844 BaseFunctionUnknownInterface **iface) const;
846 boost::shared_ptr<BaseFunctionCtx> ctx_ptr);
847
848private:
849 MatrixDouble shapeFun;
850 boost::shared_ptr<CachePhi> cachePhiPtr;
851
853};
854
855struct ContactTree;
856
864};
866
867using MatrixPtr = boost::shared_ptr<MatrixDouble>;
868using VectorPtr = boost::shared_ptr<VectorDouble>;
869
870using EntData = EntitiesFieldData::EntData;
873using FaceUserDataOperator = FaceElementForcesAndSourcesCore::UserDataOperator;
874using EleOnSide = PipelineManager::ElementsAndOpsByDim<SPACE_DIM>::FaceSideEle;
875using SideEleOp = EleOnSide::UserDataOperator;
876
877struct PhysicalEquations;
878
879struct AnalyticalExprPython;
880
882 : public boost::enable_shared_from_this<DataAtIntegrationPts> {
883
884 MatrixDouble approxPAtPts;
885 MatrixDouble approxSigmaAtPts;
886 MatrixDouble divPAtPts;
887 MatrixDouble divSigmaAtPts;
888 MatrixDouble wL2AtPts;
889 MatrixDouble wL2DotAtPts;
890 MatrixDouble wL2DotDotAtPts;
891 MatrixDouble logStretchTensorAtPts;
892 MatrixDouble stretchTensorAtPts;
893 VectorDouble jacobianAtPts;
894 MatrixDouble tractionAtPts;
895
896 MatrixDouble diffStretchTensorAtPts;
897 VectorDouble detStretchTensorAtPts;
898 MatrixDouble detStretchTensorAtPts_du;
899 MatrixDouble logStretchDotTensorAtPts;
900 MatrixDouble gradLogStretchDotTensorAtPts;
901 MatrixDouble rotAxisAtPts;
902 MatrixDouble rotAxisDotAtPts;
903 MatrixDouble rotAxisGradDotAtPts;
904 MatrixDouble rotAxis0AtPts;
905 MatrixDouble WAtPts;
906 MatrixDouble W0AtPts;
907 MatrixDouble GAtPts;
908 MatrixDouble G0AtPts;
909 MatrixDouble wH1AtPts;
910 MatrixDouble XH1AtPts;
911 MatrixDouble contactL2AtPts;
912 MatrixDouble wGradH1AtPts;
913 MatrixDouble logStretch2H1AtPts;
914 MatrixDouble logStretchTotalTensorAtPts;
915
916 MatrixDouble hAtPts;
917 MatrixDouble hdOmegaAtPts;
918 MatrixDouble hdLogStretchAtPts;
919 MatrixDouble leviKirchhoffAtPts;
920 MatrixDouble leviKirchhoffdOmegaAtPts;
921 MatrixDouble leviKirchhoffdLogStreatchAtPts;
922 MatrixDouble leviKirchhoffPAtPts;
923 MatrixDouble adjointPdstretchAtPts;
924 MatrixDouble adjointPdUAtPts;
925 MatrixDouble adjointPdUdOmegaAtPts;
926 MatrixDouble adjointPdUdPAtPts;
927 MatrixDouble rotMatAtPts;
928 MatrixDouble PAtPts;
929 MatrixDouble SigmaAtPts;
930 VectorDouble energyAtPts;
931 MatrixDouble flowL2AtPts;
932
933 MatrixDouble varRotAxis;
934 MatrixDouble varLogStreach;
935 MatrixDouble varPiola;
936 MatrixDouble varDivPiola;
937 MatrixDouble varWL2;
938
939 MatrixDouble P_du;
940
941 MatrixDouble eigenVals;
942 MatrixDouble eigenVecs;
943 VectorInt nbUniq;
944 MatrixDouble eigenValsC;
945 MatrixDouble eigenVecsC;
946 VectorInt nbUniqC;
947
948 MatrixDouble matD;
949 MatrixDouble matAxiatorD;
950 MatrixDouble matDeviatorD;
951 MatrixDouble matInvD;
952
953 MatrixDouble internalStressAtPts;
954
955 MatrixDouble facePiolaAtPts;
956 MatrixDouble gradHybridDispAtPts;
957 MatrixDouble faceMaterialForceAtPts;
958 VectorDouble normalPressureAtPts;
959
962 double piolaScale = 1.;
963 double faceEnergy = 0.;
964
965 inline auto getPiolaScalePtr() {
966 return boost::shared_ptr<double>(shared_from_this(), &piolaScale);
967 }
968
970 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
971 &approxSigmaAtPts);
972 }
974 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &approxPAtPts);
975 }
976
978 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divPAtPts);
979 }
980
982 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divSigmaAtPts);
983 }
984
986 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2AtPts);
987 }
988
990 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2DotAtPts);
991 }
992
993 inline MatrixPtr getSmallWL2DotDotAtPts() {
994 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2DotDotAtPts);
995 }
996
997 inline MatrixPtr getLogStretchTensorAtPts() {
998 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
999 &logStretchTensorAtPts);
1000 }
1001
1002 inline MatrixPtr getStretchTensorAtPts() {
1003 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1004 &stretchTensorAtPts);
1005 }
1006
1007 inline MatrixPtr getLogStretchDotTensorAtPts() {
1008 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1009 &logStretchDotTensorAtPts);
1010 }
1011
1012 inline MatrixPtr getGradLogStretchDotTensorAtPts() {
1013 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1014 &gradLogStretchDotTensorAtPts);
1015 }
1016
1018 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxisAtPts);
1019 }
1020
1022 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxis0AtPts);
1023 }
1024
1026 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1027 &rotAxisDotAtPts);
1028 }
1029
1030 inline MatrixPtr getRotAxisGradDotAtPts() {
1031 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1032 &rotAxisGradDotAtPts);
1033 }
1034
1036 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
1037 }
1038
1040 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
1041 }
1042
1044 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matD);
1045 }
1046
1048 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matInvD);
1049 }
1050
1052 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matAxiatorD);
1053 }
1054
1056 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matDeviatorD);
1057 }
1058
1060 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wH1AtPts);
1061 }
1062
1064 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &XH1AtPts);
1065 }
1066
1068 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &contactL2AtPts);
1069 }
1070
1071 inline MatrixPtr getSmallWGradH1AtPts() {
1072 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wGradH1AtPts);
1073 }
1074
1076 return boost::shared_ptr<VectorDouble>(shared_from_this(), &jacobianAtPts);
1077 }
1078
1079 inline MatrixPtr getLeviKirchhoffAtPts() {
1080 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1081 &leviKirchhoffAtPts);
1082 };
1083
1085 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varRotAxis);
1086 }
1087
1088 inline MatrixPtr getVarLogStreachPts() {
1089 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varLogStreach);
1090 }
1091
1093 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varPiola);
1094 }
1095
1097 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varDivPiola);
1098 }
1099
1101 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varWL2);
1102 }
1103
1104 inline MatrixPtr getInternalStressAtPts() {
1105 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1106 &internalStressAtPts);
1107 }
1108
1110 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &facePiolaAtPts);
1111 }
1112
1113 inline MatrixPtr getGradHybridDispAtPts() {
1114 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1115 &gradHybridDispAtPts);
1116 }
1117
1119 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &eigenVals);
1120 }
1121
1122 boost::shared_ptr<PhysicalEquations> physicsPtr;
1123};
1124
1126
1127
1128struct ExternalStrain;
1130struct PhysicalEquations {
1131
1138
1140
1143
1144 PhysicalEquations() = delete;
1145 PhysicalEquations(const int size_active, const int size_dependent)
1146 : activeVariables(size_active, 0),
1147 dependentVariablesPiola(size_dependent, 0),
1148 dependentVariablesPiolaDirevatives(size_dependent * size_active, 0) {}
1149 virtual ~PhysicalEquations() = default;
1150
1151 virtual MoFEMErrorCode recordTape(
const int tag, DTensor2Ptr *t_h) = 0;
1152
1154 returnOpJacobian(const int tag, const bool eval_rhs, const bool eval_lhs,
1155 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1156 boost::shared_ptr<PhysicalEquations> physics_ptr);
1157
1159 returnOpSpatialPhysical(
const std::string &
field_name,
1160 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1161 const double alpha_u);
1162
1164 returnOpSpatialPhysicalExternalStrain(
1166 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1167 boost::shared_ptr<ExternalStrainVec> external_strain_vec_ptr,
1168 std::map<std::string, boost::shared_ptr<ScalingMethod>> smv);
1169
1171 std::string row_field, std::string col_field,
1172 boost::shared_ptr<DataAtIntegrationPts> data_ptr, const double alpha);
1173
1175 returnOpCalculateEnergy(boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1176 boost::shared_ptr<double> total_energy_ptr);
1177
1179 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1180 boost::shared_ptr<PhysicalEquations> physics_ptr);
1181
1183 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1184 boost::shared_ptr<PhysicalEquations> physics_ptr);
1185
1187 returnOpSetScale(boost::shared_ptr<double> scale_ptr,
1188 boost::shared_ptr<PhysicalEquations> physics_ptr);
1189
1190 std::vector<double> activeVariables;
1191 std::vector<double> dependentVariablesPiola;
1192 std::vector<double> dependentVariablesPiolaDirevatives;
1193
1194
1195
1196
1197
1198 template <int S>
1199 inline static DTensor2Ptr get_VecTensor2(std::vector<double> &
v) {
1200 return DTensor2Ptr(&
v[S + 0], &
v[S + 1], &
v[S + 2], &
v[S + 3], &
v[S + 4],
1201 &
v[S + 5], &
v[S + 6], &
v[S + 7], &
v[S + 8]);
1202 }
1203
1204 template <int S>
1205 inline static DTensor0Ptr get_VecTensor0(std::vector<double> &
v) {
1206 return DTensor0Ptr(&
v[S + 0]);
1207 }
1208
1209 template <int S0>
1210 inline static DTensor3Ptr get_vecTensor3(std::vector<double> &
v,
1211 const int nba) {
1212
1213 const int A00 = nba * 0 + S0;
1214 const int A01 = nba * 1 + S0;
1215 const int A02 = nba * 2 + S0;
1216 const int A10 = nba * 3 + S0;
1217 const int A11 = nba * 4 + S0;
1218 const int A12 = nba * 5 + S0;
1219 const int A20 = nba * 6 + S0;
1220 const int A21 = nba * 7 + S0;
1221 const int A22 = nba * 8 + S0;
1222
1223 return DTensor3Ptr(
1224
1225 &
v[A00 + 0], &
v[A00 + 1], &
v[A00 + 2], &
v[A01 + 0], &
v[A01 + 1],
1226 &
v[A01 + 2], &
v[A02 + 0], &
v[A02 + 1], &
v[A02 + 2],
1227
1228 &
v[A10 + 0], &
v[A10 + 1], &
v[A10 + 2], &
v[A11 + 0], &
v[A11 + 1],
1229 &
v[A11 + 2], &
v[A12 + 0], &
v[A12 + 1], &
v[A12 + 2],
1230
1231 &
v[A20 + 0], &
v[A20 + 1], &
v[A20 + 2], &
v[A21 + 0], &
v[A21 + 1],
1232 &
v[A21 + 2], &
v[A22 + 0], &
v[A22 + 1], &
v[A22 + 2]
1233
1234 );
1235 }
1236
1237
1238
1239
1240
1241
1242
1243 inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
1244
1245
1246
1247
1248
1249
1250
1251 inline DTensor2Ptr get_P() {
1252 return get_VecTensor2<0>(dependentVariablesPiola);
1253 }
1254
1255
1256
1257
1258
1259
1260
1261 inline DTensor3Ptr get_P_dh0() {
1262 return get_vecTensor3<0>(dependentVariablesPiolaDirevatives,
1263 activeVariables.size());
1264 }
1265 inline DTensor3Ptr get_P_dh1() {
1266 return get_vecTensor3<3>(dependentVariablesPiolaDirevatives,
1267 activeVariables.size());
1268 }
1269 inline DTensor3Ptr get_P_dh2() {
1270 return get_vecTensor3<6>(dependentVariablesPiolaDirevatives,
1271 activeVariables.size());
1272 }
1273
1274
1275};
1276
1277struct BcDisp {
1278 BcDisp(std::string name, std::vector<double> attr,
Range faces);
1279 std::string blockName;
1281 VectorDouble3 vals;
1282 VectorInt3 flags;
1283};
1285
1286struct BcRot {
1287 BcRot(std::string name, std::vector<double> attr,
Range faces);
1288 std::string blockName;
1290 VectorDouble vals;
1291 double theta;
1292};
1293using BcRotVec = std::vector<BcRot>;
1294
1296
1297struct TractionBc {
1298 TractionBc(std::string name, std::vector<double> attr,
Range faces);
1299 std::string blockName;
1301 VectorDouble3 vals;
1302 VectorInt3 flags;
1303};
1305
1306struct NormalDisplacementBc {
1307 NormalDisplacementBc(std::string name, std::vector<double> attr,
Range faces);
1308 std::string blockName;
1310 double val;
1311};
1313
1314struct AnalyticalDisplacementBc {
1315 AnalyticalDisplacementBc(std::string name, std::vector<double> attr,
1317 std::string blockName;
1319 VectorInt3 flags;
1320};
1322
1323struct AnalyticalTractionBc {
1324 AnalyticalTractionBc(std::string name, std::vector<double> attr,
Range faces);
1325 std::string blockName;
1327 VectorInt3 flags;
1328};
1330
1331struct PressureBc {
1332 PressureBc(std::string name, std::vector<double> attr,
Range faces);
1333 std::string blockName;
1335 double val;
1336};
1338
1339struct ExternalStrain {
1340 ExternalStrain(std::string name, std::vector<double> attr,
Range ents);
1341 std::string blockName;
1343 double val;
1344 double bulkModulusK;
1345};
1346
1347 #ifdef ENABLE_PYTHON_BINDING
1348struct AnalyticalExprPython {
1349 AnalyticalExprPython() = default;
1350 virtual ~AnalyticalExprPython() = default;
1351
1353 MoFEMErrorCode evalAnalyticalDisp(
double delta_t,
double t, np::ndarray x,
1354 np::ndarray y, np::ndarray z,
1355 np::ndarray nx, np::ndarray ny,
1356 np::ndarray nz,
1357 const std::string &block_name,
1358 np::ndarray &analytical_expr);
1359 MoFEMErrorCode evalAnalyticalTraction(
double delta_t,
double t, np::ndarray x,
1360 np::ndarray y, np::ndarray z,
1361 np::ndarray nx, np::ndarray ny,
1362 np::ndarray nz,
1363 const std::string &block_name,
1364 np::ndarray &analytical_expr);
1365
1366 template <typename T>
1367 inline std::vector<T>
1368 py_list_to_std_vector(const boost::python::object &iterable) {
1369 return std::vector<T>(boost::python::stl_input_iterator<T>(iterable),
1370 boost::python::stl_input_iterator<T>());
1371 }
1372
1373private:
1374 bp::object mainNamespace;
1375 bp::object analyticalDispFun;
1376 bp::object analyticalTractionFun;
1377};
1378
1379extern boost::weak_ptr<AnalyticalExprPython> AnalyticalExprPythonWeakPtr;
1380 #endif
1381
1384
1385}
1386
1387#endif
1388
ForcesAndSourcesCore::UserDataOperator UserDataOperator
const double v
phase velocity of light in medium (cm/ns)
boost::shared_ptr< VectorDouble > VectorPtr
std::vector< AnalyticalTractionBc > AnalyticalTractionBcVec
std::vector< BcDisp > BcDispVec
std::vector< NormalDisplacementBc > NormalDisplacementBcVec
std::vector< BcRot > BcRotVec
std::vector< ExternalStrain > ExternalStrainVec
std::vector< AnalyticalDisplacementBc > AnalyticalDisplacementBcVec
std::vector< Range > TractionFreeBc
std::vector< PressureBc > PressureBcVec
std::vector< TractionBc > TractionBcVec
boost::shared_ptr< MatrixDouble > MatrixPtr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
constexpr double t
plate stiffness
constexpr auto field_name
PipelineManager::ElementsAndOpsByDim< SPACE_DIM >::FaceSideEle EleOnSide
Data on single entity (This is passed as argument to DataOperator::doWork)