848 {
849
851
852 using CachePhi = boost::tuple<int, int, MatrixDouble>;
853
854 CGGUserPolynomialBase(boost::shared_ptr<CachePhi> cache_phi = nullptr);
855 ~CGGUserPolynomialBase() = default;
856
857 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
858 BaseFunctionUnknownInterface **iface) const;
860 boost::shared_ptr<BaseFunctionCtx> ctx_ptr);
861
862private:
863 MatrixDouble shapeFun;
864 boost::shared_ptr<CachePhi> cachePhiPtr;
865
867};
868
869struct ContactTree;
870
878};
880
881using MatrixPtr = boost::shared_ptr<MatrixDouble>;
882using VectorPtr = boost::shared_ptr<VectorDouble>;
883
884using EntData = EntitiesFieldData::EntData;
887using FaceUserDataOperator = FaceElementForcesAndSourcesCore::UserDataOperator;
888using EleOnSide = PipelineManager::ElementsAndOpsByDim<SPACE_DIM>::FaceSideEle;
889using SideEleOp = EleOnSide::UserDataOperator;
890
891struct PhysicalEquations;
892
893struct AnalyticalExprPython;
894
896 : public boost::enable_shared_from_this<DataAtIntegrationPts> {
897
898 MatrixDouble approxPAtPts;
899 MatrixDouble approxSigmaAtPts;
900 MatrixDouble divPAtPts;
901 MatrixDouble divSigmaAtPts;
902 MatrixDouble wL2AtPts;
903 MatrixDouble wL2DotAtPts;
904 MatrixDouble wL2DotDotAtPts;
905 MatrixDouble logStretchTensorAtPts;
906 MatrixDouble stretchTensorAtPts;
907 VectorDouble jacobianAtPts;
908 MatrixDouble tractionAtPts;
909
910 MatrixDouble diffStretchTensorAtPts;
911 VectorDouble detStretchTensorAtPts;
912 MatrixDouble detStretchTensorAtPts_du;
913 MatrixDouble logStretchDotTensorAtPts;
914 MatrixDouble gradLogStretchDotTensorAtPts;
915 MatrixDouble rotAxisAtPts;
916 MatrixDouble rotAxisDotAtPts;
917 MatrixDouble rotAxisGradDotAtPts;
918 MatrixDouble rotAxis0AtPts;
919 MatrixDouble WAtPts;
920 MatrixDouble W0AtPts;
921 MatrixDouble GAtPts;
922 MatrixDouble G0AtPts;
923 MatrixDouble wH1AtPts;
924 MatrixDouble XH1AtPts;
925 MatrixDouble contactL2AtPts;
926 MatrixDouble wGradH1AtPts;
927 MatrixDouble logStretch2H1AtPts;
928 MatrixDouble logStretchTotalTensorAtPts;
929
930 MatrixDouble hAtPts;
931 MatrixDouble hdOmegaAtPts;
932 MatrixDouble hdLogStretchAtPts;
933 MatrixDouble leviKirchhoffAtPts;
934 MatrixDouble leviKirchhoffdOmegaAtPts;
935 MatrixDouble leviKirchhoffdLogStreatchAtPts;
936 MatrixDouble leviKirchhoffPAtPts;
937 MatrixDouble adjointPdstretchAtPts;
938 MatrixDouble adjointPdUAtPts;
939 MatrixDouble adjointPdUdOmegaAtPts;
940 MatrixDouble adjointPdUdPAtPts;
941 MatrixDouble rotMatAtPts;
942 MatrixDouble PAtPts;
943 MatrixDouble SigmaAtPts;
944 VectorDouble energyAtPts;
945 MatrixDouble flowL2AtPts;
946
947 MatrixDouble varRotAxis;
948 MatrixDouble varLogStreach;
949 MatrixDouble varPiola;
950 MatrixDouble varDivPiola;
951 MatrixDouble varWL2;
952
953 MatrixDouble P_du;
954
955 MatrixDouble eigenVals;
956 MatrixDouble eigenVecs;
957 VectorInt nbUniq;
958 MatrixDouble eigenValsC;
959 MatrixDouble eigenVecsC;
960 VectorInt nbUniqC;
961
962 MatrixDouble matD;
963 MatrixDouble matAxiatorD;
964 MatrixDouble matDeviatorD;
965 MatrixDouble matInvD;
966
967 MatrixDouble internalStressAtPts;
968
969 MatrixDouble facePiolaAtPts;
970 MatrixDouble gradHybridDispAtPts;
971 MatrixDouble faceMaterialForceAtPts;
972 VectorDouble normalPressureAtPts;
973
976 double piolaScale = 1.;
977 double faceEnergy = 0.;
978
979 inline auto getPiolaScalePtr() {
980 return boost::shared_ptr<double>(shared_from_this(), &piolaScale);
981 }
982
984 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
985 &approxSigmaAtPts);
986 }
988 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &approxPAtPts);
989 }
990
992 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divPAtPts);
993 }
994
996 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divSigmaAtPts);
997 }
998
1000 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2AtPts);
1001 }
1002
1003 inline MatrixPtr getSmallWL2DotAtPts() {
1004 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2DotAtPts);
1005 }
1006
1007 inline MatrixPtr getSmallWL2DotDotAtPts() {
1008 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wL2DotDotAtPts);
1009 }
1010
1011 inline MatrixPtr getLogStretchTensorAtPts() {
1012 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1013 &logStretchTensorAtPts);
1014 }
1015
1016 inline MatrixPtr getStretchTensorAtPts() {
1017 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1018 &stretchTensorAtPts);
1019 }
1020
1021 inline MatrixPtr getLogStretchDotTensorAtPts() {
1022 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1023 &logStretchDotTensorAtPts);
1024 }
1025
1026 inline MatrixPtr getGradLogStretchDotTensorAtPts() {
1027 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1028 &gradLogStretchDotTensorAtPts);
1029 }
1030
1032 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxisAtPts);
1033 }
1034
1036 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxis0AtPts);
1037 }
1038
1040 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1041 &rotAxisDotAtPts);
1042 }
1043
1044 inline MatrixPtr getRotAxisGradDotAtPts() {
1045 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1046 &rotAxisGradDotAtPts);
1047 }
1048
1050 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
1051 }
1052
1054 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
1055 }
1056
1058 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matD);
1059 }
1060
1062 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matInvD);
1063 }
1064
1066 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matAxiatorD);
1067 }
1068
1070 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matDeviatorD);
1071 }
1072
1074 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wH1AtPts);
1075 }
1076
1078 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &XH1AtPts);
1079 }
1080
1082 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &contactL2AtPts);
1083 }
1084
1085 inline MatrixPtr getSmallWGradH1AtPts() {
1086 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wGradH1AtPts);
1087 }
1088
1090 return boost::shared_ptr<VectorDouble>(shared_from_this(), &jacobianAtPts);
1091 }
1092
1093 inline MatrixPtr getLeviKirchhoffAtPts() {
1094 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1095 &leviKirchhoffAtPts);
1096 };
1097
1099 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varRotAxis);
1100 }
1101
1102 inline MatrixPtr getVarLogStreachPts() {
1103 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varLogStreach);
1104 }
1105
1107 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varPiola);
1108 }
1109
1111 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varDivPiola);
1112 }
1113
1115 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &varWL2);
1116 }
1117
1118 inline MatrixPtr getInternalStressAtPts() {
1119 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1120 &internalStressAtPts);
1121 }
1122
1124 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &facePiolaAtPts);
1125 }
1126
1127 inline MatrixPtr getGradHybridDispAtPts() {
1128 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
1129 &gradHybridDispAtPts);
1130 }
1131
1133 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &eigenVals);
1134 }
1135
1136 boost::shared_ptr<PhysicalEquations> physicsPtr;
1137};
1138
1140
1141
1142struct ExternalStrain;
1144struct PhysicalEquations {
1145
1152
1154
1157
1158 PhysicalEquations() = delete;
1159 PhysicalEquations(const int size_active, const int size_dependent)
1160 : activeVariables(size_active, 0),
1161 dependentVariablesPiola(size_dependent, 0),
1162 dependentVariablesPiolaDirevatives(size_dependent * size_active, 0) {}
1163 virtual ~PhysicalEquations() = default;
1164
1165 virtual MoFEMErrorCode recordTape(
const int tag, DTensor2Ptr *t_h) = 0;
1166
1168 returnOpJacobian(const int tag, const bool eval_rhs, const bool eval_lhs,
1169 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1170 boost::shared_ptr<PhysicalEquations> physics_ptr);
1171
1173 returnOpSpatialPhysical(
const std::string &
field_name,
1174 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1175 const double alpha_u);
1176
1178 returnOpSpatialPhysicalExternalStrain(
1180 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1181 boost::shared_ptr<ExternalStrainVec> external_strain_vec_ptr,
1182 std::map<std::string, boost::shared_ptr<ScalingMethod>> smv);
1183
1185 std::string row_field, std::string col_field,
1186 boost::shared_ptr<DataAtIntegrationPts> data_ptr, const double alpha);
1187
1189 returnOpCalculateEnergy(boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1190 boost::shared_ptr<double> total_energy_ptr);
1191
1193 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1194 boost::shared_ptr<PhysicalEquations> physics_ptr);
1195
1197 boost::shared_ptr<DataAtIntegrationPts> data_ptr,
1198 boost::shared_ptr<PhysicalEquations> physics_ptr);
1199
1201 returnOpSetScale(boost::shared_ptr<double> scale_ptr,
1202 boost::shared_ptr<PhysicalEquations> physics_ptr);
1203
1204 std::vector<double> activeVariables;
1205 std::vector<double> dependentVariablesPiola;
1206 std::vector<double> dependentVariablesPiolaDirevatives;
1207
1208
1209
1210
1211
1212 template <int S>
1213 inline static DTensor2Ptr get_VecTensor2(std::vector<double> &
v) {
1214 return DTensor2Ptr(&
v[S + 0], &
v[S + 1], &
v[S + 2], &
v[S + 3], &
v[S + 4],
1215 &
v[S + 5], &
v[S + 6], &
v[S + 7], &
v[S + 8]);
1216 }
1217
1218 template <int S>
1219 inline static DTensor0Ptr get_VecTensor0(std::vector<double> &
v) {
1220 return DTensor0Ptr(&
v[S + 0]);
1221 }
1222
1223 template <int S0>
1224 inline static DTensor3Ptr get_vecTensor3(std::vector<double> &
v,
1225 const int nba) {
1226
1227 const int A00 = nba * 0 + S0;
1228 const int A01 = nba * 1 + S0;
1229 const int A02 = nba * 2 + S0;
1230 const int A10 = nba * 3 + S0;
1231 const int A11 = nba * 4 + S0;
1232 const int A12 = nba * 5 + S0;
1233 const int A20 = nba * 6 + S0;
1234 const int A21 = nba * 7 + S0;
1235 const int A22 = nba * 8 + S0;
1236
1237 return DTensor3Ptr(
1238
1239 &
v[A00 + 0], &
v[A00 + 1], &
v[A00 + 2], &
v[A01 + 0], &
v[A01 + 1],
1240 &
v[A01 + 2], &
v[A02 + 0], &
v[A02 + 1], &
v[A02 + 2],
1241
1242 &
v[A10 + 0], &
v[A10 + 1], &
v[A10 + 2], &
v[A11 + 0], &
v[A11 + 1],
1243 &
v[A11 + 2], &
v[A12 + 0], &
v[A12 + 1], &
v[A12 + 2],
1244
1245 &
v[A20 + 0], &
v[A20 + 1], &
v[A20 + 2], &
v[A21 + 0], &
v[A21 + 1],
1246 &
v[A21 + 2], &
v[A22 + 0], &
v[A22 + 1], &
v[A22 + 2]
1247
1248 );
1249 }
1250
1251
1252
1253
1254
1255
1256
1257 inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
1258
1259
1260
1261
1262
1263
1264
1265 inline DTensor2Ptr get_P() {
1266 return get_VecTensor2<0>(dependentVariablesPiola);
1267 }
1268
1269
1270
1271
1272
1273
1274
1275 inline DTensor3Ptr get_P_dh0() {
1276 return get_vecTensor3<0>(dependentVariablesPiolaDirevatives,
1277 activeVariables.size());
1278 }
1279 inline DTensor3Ptr get_P_dh1() {
1280 return get_vecTensor3<3>(dependentVariablesPiolaDirevatives,
1281 activeVariables.size());
1282 }
1283 inline DTensor3Ptr get_P_dh2() {
1284 return get_vecTensor3<6>(dependentVariablesPiolaDirevatives,
1285 activeVariables.size());
1286 }
1287
1288
1289};
1290
1291struct BcDisp {
1292 BcDisp(std::string name, std::vector<double> attr,
Range faces);
1293 std::string blockName;
1295 VectorDouble3 vals;
1296 VectorInt3 flags;
1297};
1299
1300struct BcRot {
1301 BcRot(std::string name, std::vector<double> attr,
Range faces);
1302 std::string blockName;
1304 VectorDouble vals;
1305 double theta;
1306};
1307using BcRotVec = std::vector<BcRot>;
1308
1310
1311struct TractionBc {
1312 TractionBc(std::string name, std::vector<double> attr,
Range faces);
1313 std::string blockName;
1315 VectorDouble3 vals;
1316 VectorInt3 flags;
1317};
1319
1320struct NormalDisplacementBc {
1321 NormalDisplacementBc(std::string name, std::vector<double> attr,
Range faces);
1322 std::string blockName;
1324 double val;
1325};
1327
1328struct AnalyticalDisplacementBc {
1329 AnalyticalDisplacementBc(std::string name, std::vector<double> attr,
1331 std::string blockName;
1333 VectorInt3 flags;
1334};
1336
1337struct AnalyticalTractionBc {
1338 AnalyticalTractionBc(std::string name, std::vector<double> attr,
Range faces);
1339 std::string blockName;
1341 VectorInt3 flags;
1342};
1344
1345struct PressureBc {
1346 PressureBc(std::string name, std::vector<double> attr,
Range faces);
1347 std::string blockName;
1349 double val;
1350};
1352
1353struct ExternalStrain {
1354 ExternalStrain(std::string name, std::vector<double> attr,
Range ents);
1355 std::string blockName;
1357 double val;
1358 double bulkModulusK;
1359};
1360
1361 #ifdef ENABLE_PYTHON_BINDING
1362struct AnalyticalExprPython {
1363 AnalyticalExprPython() = default;
1364 virtual ~AnalyticalExprPython() = default;
1365
1367 MoFEMErrorCode evalAnalyticalDisp(
double delta_t,
double t, np::ndarray x,
1368 np::ndarray y, np::ndarray z,
1369 np::ndarray nx, np::ndarray ny,
1370 np::ndarray nz,
1371 const std::string &block_name,
1372 np::ndarray &analytical_expr);
1373 MoFEMErrorCode evalAnalyticalTraction(
double delta_t,
double t, np::ndarray x,
1374 np::ndarray y, np::ndarray z,
1375 np::ndarray nx, np::ndarray ny,
1376 np::ndarray nz,
1377 const std::string &block_name,
1378 np::ndarray &analytical_expr);
1379
1380 template <typename T>
1381 inline std::vector<T>
1382 py_list_to_std_vector(const boost::python::object &iterable) {
1383 return std::vector<T>(boost::python::stl_input_iterator<T>(iterable),
1384 boost::python::stl_input_iterator<T>());
1385 }
1386
1387private:
1388 bp::object mainNamespace;
1389 bp::object analyticalDispFun;
1390 bp::object analyticalTractionFun;
1391};
1392
1393extern boost::weak_ptr<AnalyticalExprPython> AnalyticalExprPythonWeakPtr;
1394 #endif
1395
1398
1399}
1400
1401#endif
1402
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)