v0.15.4
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Public Attributes | List of all members
EshelbianMonitor Struct Reference
Inheritance diagram for EshelbianMonitor:
[legend]
Collaboration diagram for EshelbianMonitor:
[legend]

Public Types

using Ele = ForcesAndSourcesCore
 
using VolEle = VolumeElementForcesAndSourcesCore
 
using VolOp = VolumeElementForcesAndSourcesCore::UserDataOperator
 
using SetPtsData = FieldEvaluatorInterface::SetPtsData
 
using PtsHashMap = std::map< std::string, std::array< double, 3 > >
 
using FaceSideEle = PipelineManager::ElementsAndOpsByDim< SPACE_DIM >::FaceSideEle
 

Public Member Functions

 EshelbianMonitor (EshelbianCore &ep)
 
MoFEMErrorCode preProcess ()
 
MoFEMErrorCode operator() ()
 
MoFEMErrorCode postProcess ()
 

Public Attributes

EshelbianCoreeP
 
boost::shared_ptr< SetPtsDatadataFieldEval
 
boost::shared_ptr< VolElevolPostProcEnergy
 
boost::shared_ptr< doublegEnergy
 
PtsHashMap ptsHashMap
 
PetscBool writeRestart = PETSC_FALSE
 

Detailed Description

Definition at line 6 of file EshelbianMonitor.cpp.

Member Typedef Documentation

◆ Ele

Definition at line 8 of file EshelbianMonitor.cpp.

◆ FaceSideEle

using EshelbianMonitor::FaceSideEle = PipelineManager::ElementsAndOpsByDim<SPACE_DIM>::FaceSideEle

Definition at line 13 of file EshelbianMonitor.cpp.

◆ PtsHashMap

using EshelbianMonitor::PtsHashMap = std::map<std::string, std::array<double, 3> >

Definition at line 12 of file EshelbianMonitor.cpp.

◆ SetPtsData

using EshelbianMonitor::SetPtsData = FieldEvaluatorInterface::SetPtsData

Definition at line 11 of file EshelbianMonitor.cpp.

◆ VolEle

using EshelbianMonitor::VolEle = VolumeElementForcesAndSourcesCore

Definition at line 9 of file EshelbianMonitor.cpp.

◆ VolOp

using EshelbianMonitor::VolOp = VolumeElementForcesAndSourcesCore::UserDataOperator

Definition at line 10 of file EshelbianMonitor.cpp.

Constructor & Destructor Documentation

◆ EshelbianMonitor()

EshelbianMonitor::EshelbianMonitor ( EshelbianCore ep)
inline

Definition at line 23 of file EshelbianMonitor.cpp.

24 : eP(ep), dataFieldEval(ep.mField.getInterface<FieldEvaluatorInterface>()
25 ->getData<VolEle>()),
26 volPostProcEnergy(new VolEle(ep.mField)), gEnergy(new double) {
28 ep.mField.getInterface<FieldEvaluatorInterface>()->buildTree<SPACE_DIM>(
29 dataFieldEval, "EP"),
30 "build field evaluator tree");
31
32 auto no_rule = [](int, int, int) { return -1; };
33
34 auto set_element_for_field_eval = [&]() {
36 boost::shared_ptr<Ele> vol_ele(dataFieldEval->feMethodPtr.lock());
37 vol_ele->getRuleHook = no_rule;
38 vol_ele->getUserPolynomialBase() =
39 boost::make_shared<CGGUserPolynomialBase>();
40 EshelbianPlasticity::AddHOOps<SPACE_DIM, SPACE_DIM, SPACE_DIM>::add(
41 vol_ele->getOpPtrVector(), {HDIV, H1, L2}, eP.materialH1Positions,
42 ep.frontAdjEdges);
43
44 auto piola_scale_ptr = boost::make_shared<double>(1.0);
45 vol_ele->getOpPtrVector().push_back(
46 eP.physicalEquations->returnOpSetScale(piola_scale_ptr,
48 vol_ele->getOpPtrVector().push_back(new OpCalculateHVecTensorField<3, 3>(
49 eP.piolaStress, eP.dataAtPts->getApproxPAtPts(), piola_scale_ptr));
50 vol_ele->getOpPtrVector().push_back(
51 new OpCalculateHTensorTensorField<3, 3>(
52 eP.bubbleField, eP.dataAtPts->getApproxPAtPts(), piola_scale_ptr,
53 SmartPetscObj<Vec>(), MBMAXTYPE));
55 vol_ele->getOpPtrVector().push_back(
56 eP.physicalEquations->returnOpCalculateStretchFromStress(
58 } else {
59 vol_ele->getOpPtrVector().push_back(
60 new OpCalculateTensor2SymmetricFieldValues<3>(
61 eP.stretchTensor, eP.dataAtPts->getLogStretchTensorAtPts(),
62 MBTET));
63 }
64 vol_ele->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
65 eP.rotAxis, eP.dataAtPts->getRotAxisAtPts(), MBTET));
66 vol_ele->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
67 eP.rotAxis, eP.dataAtPts->getRotAxis0AtPts(), eP.solTSStep, MBTET));
68 vol_ele->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
69 eP.spatialL2Disp, eP.dataAtPts->getSmallWL2AtPts(), MBTET));
70
71 // H1 displacements
72 vol_ele->getOpPtrVector().push_back(new OpCalculateVectorFieldValues<3>(
73 eP.spatialH1Disp, eP.dataAtPts->getSmallWH1AtPts()));
74 vol_ele->getOpPtrVector().push_back(
75 new OpCalculateVectorFieldGradient<3, 3>(
76 eP.spatialH1Disp, eP.dataAtPts->getSmallWGradH1AtPts()));
77
78 vol_ele->getOpPtrVector().push_back(
81 };
82
83 auto set_element_for_post_process_energy = [&]() {
85
87 1, false, false, false, SmartPetscObj<Vec>(), volPostProcEnergy);
88
89 if (auto op = eP.physicalEquations->returnOpCalculateEnergy(eP.dataAtPts,
90 gEnergy)) {
91
92 if (eP.noStretch ||
94 // We have to use actual strains to evaluate J integral and energy,
95 // in this case. Note actual stresses, and actual energy can only
96 // drive crack growth
97
98 // Note: Calling this bellow we overwrite the stretch calculation
99 // from setBaseVolumeElementOps
100
101 volPostProcEnergy->getOpPtrVector().push_back(
102 eP.physicalEquations->returnOpCalculateStretchFromStress(
104 }
105
106 volPostProcEnergy->getOpPtrVector().push_back(op);
107 }
108
110 };
111
112 auto reads_post_proc_data = [](PtsHashMap &pts_hash_map) {
114 std::ifstream file(
115 "points.txt"); // Open the file with the name "data.txt"
116
117 if (!file.is_open()) {
119 }
120
121 std::string line;
122
123 while (std::getline(file, line)) {
124 std::istringstream iss(line);
125 std::string col1;
126 double col2, col3, col4;
127
128 if (iss >> col1 >> col2 >> col3 >> col4) {
129 MOFEM_LOG("EP", Sev::verbose) << "Read: " << col1 << ", " << col2
130 << ", " << col3 << ", " << col4;
131 pts_hash_map[col1] = {col2, col3, col4};
132 } else {
133 MOFEM_LOG("EP", Sev::error) << "Error parsing line: " << line;
134 }
135 }
136
137 file.close(); // Close the file
139 };
140
141 CHK_THROW_MESSAGE(set_element_for_field_eval(), "set element for field");
142 CHK_THROW_MESSAGE(set_element_for_post_process_energy(),
143 "set element for post energy");
144
145 PetscBool test_cook_flg = PETSC_FALSE;
146 CHK_THROW_MESSAGE(PetscOptionsGetBool(PETSC_NULLPTR, "", "-test_cook_pts",
147 &test_cook_flg, PETSC_NULLPTR),
148 "get post proc points");
149 if (test_cook_flg) {
150 ptsHashMap["Point A"] = {48., 60., 4.999};
151 ptsHashMap["Point B"] = {48. / 2., 44. + (60. - 44.) / 2., 0.};
152 ptsHashMap["Point C"] = {48. / 2., (44. - 0.) / 2., 0.};
153 }
154
155 PetscInt atom_test = 0;
156
157 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "", "-atom_test", &atom_test,
158 PETSC_NULLPTR);
159 if (atom_test == 14) {
160 // Points for atom test 14: check external strain
161 ptsHashMap["Point (2.5, 0., 0.)"] = {2.5, 0, 0.};
162 }
163 CHK_MOAB_THROW(reads_post_proc_data(ptsHashMap), "read post proc points");
164
165 CHK_THROW_MESSAGE(PetscOptionsGetBool(PETSC_NULLPTR, "", "-write_restart",
166 &writeRestart, PETSC_NULLPTR),
167 "get write restart option");
168
169 MOFEM_LOG("EP", Sev::inform)
170 << "EshelbianMonitor writeRestart " << writeRestart
171 ? " true"
172 : " false";
173 }
constexpr int SPACE_DIM
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MOFEM_LOG(channel, severity)
Log.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
boost::shared_ptr< Range > frontAdjEdges
MoFEM::Interface & mField
const std::string spatialL2Disp
const std::string materialH1Positions
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, const bool calc_rates, SmartPetscObj< Vec > ver_vec, boost::shared_ptr< VolumeElementForcesAndSourcesCore > fe)
const std::string spatialH1Disp
static enum MaterialModel materialModel
const std::string piolaStress
const std::string bubbleField
static PetscBool noStretch
boost::shared_ptr< PhysicalEquations > physicalEquations
const std::string rotAxis
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
SmartPetscObj< Vec > solTSStep
const std::string stretchTensor
boost::shared_ptr< double > gEnergy
boost::shared_ptr< VolEle > volPostProcEnergy
boost::shared_ptr< SetPtsData > dataFieldEval
EshelbianCore & eP
std::map< std::string, std::array< double, 3 > > PtsHashMap
VolumeElementForcesAndSourcesCore VolEle
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
int atom_test
Atom test.
Definition plastic.cpp:121

Member Function Documentation

◆ operator()()

MoFEMErrorCode EshelbianMonitor::operator() ( )
inline

Definition at line 177 of file EshelbianMonitor.cpp.

177{ return 0; }

◆ postProcess()

MoFEMErrorCode EshelbianMonitor::postProcess ( )
inline

Definition at line 179 of file EshelbianMonitor.cpp.

179 {
181
182 MOFEM_LOG("EP", Sev::inform) << "Monitor postProcess";
183
184 auto get_step = [](auto ts_step) {
185 std::ostringstream ss;
186 ss << boost::str(boost::format("%d") % static_cast<int>(ts_step));
187 std::string s = ss.str();
188 return s;
189 };
190
191 auto calculate_reaction_forces = [&]() {
193
194 // Get boundary faces marked in blocks name "SPATIAL_DISP_...".
195 auto get_ents_on_mesh_skin = [&]() {
196 std::map<std::string, Range> boundary_entities_vec;
197
198 auto get_block_vec_impl = [&](auto block_name) {
199 return eP.mField.getInterface<MeshsetsManager>()->getCubitMeshsetPtr(
200 std::regex(
201
202 (boost::format("%s(.*)") % block_name).str()
203
204 ));
205 };
206
207 auto add_blockset_ents_impl = [&](auto &&vec) {
209 for (auto it : vec) {
210 Range boundary_entities;
211 CHKERR it->getMeshsetIdEntitiesByDimension(eP.mField.get_moab(), 2,
212 boundary_entities, true);
213 std::string meshset_name = it->getName();
214 boundary_entities_vec[meshset_name] = boundary_entities;
215 boundary_entities.clear();
216 }
218 };
219
220 for (auto block_name : {"SPATIAL_DISP", "FIX", "CONTACT",
221 "SPATIAL_ROTATION", "NORMAL_DISPLACEMENT"}) {
223 add_blockset_ents_impl(get_block_vec_impl(block_name)),
224 "add blockset entities");
225 }
226
227 return boundary_entities_vec;
228 };
229
230 auto boundary_entities_vec = get_ents_on_mesh_skin();
231
232 std::vector<std::tuple<std::string, Range, std::array<double, 6>>>
233 reactionForces;
234 for (const auto &pair : boundary_entities_vec) {
235 reactionForces.push_back(std::make_tuple(
236 pair.first, pair.second,
237 std::array<double, 6>{0.0, 0.0, 0.0, 0.0, 0.0, 0.0}));
238 }
239
240 auto integration_rule_face = [](int, int, int approx_order) {
241 return 2 * approx_order + 1;
242 };
243 auto face_fe =
244 boost::make_shared<FaceElementForcesAndSourcesCore>(eP.mField);
245 auto no_rule = [](int, int, int) { return -1; };
246 face_fe->getRuleHook = integration_rule_face;
247 CHKERR
248 EshelbianPlasticity::AddHOOps<SPACE_DIM - 1, SPACE_DIM, SPACE_DIM>::add(
249 face_fe->getOpPtrVector(), {L2}, eP.materialH1Positions,
251
252 auto op_side =
253 new OpLoopSide<FaceSideEle>(eP.mField, "EP", SPACE_DIM, Sev::noisy);
254 face_fe->getOpPtrVector().push_back(op_side);
255 auto side_fe_ptr = op_side->getSideFEPtr();
256 auto base_ptr =
257 boost::make_shared<EshelbianPlasticity::CGGUserPolynomialBase>();
258 side_fe_ptr->getUserPolynomialBase() = base_ptr;
259 CHKERR
260 EshelbianPlasticity::AddHOOps<SPACE_DIM, SPACE_DIM, SPACE_DIM>::add(
261 side_fe_ptr->getOpPtrVector(), {H1, HDIV, L2}, eP.materialH1Positions,
263 auto piola_scale_ptr = boost::make_shared<double>(1.0);
264 side_fe_ptr->getOpPtrVector().push_back(
265 eP.physicalEquations->returnOpSetScale(piola_scale_ptr,
267 side_fe_ptr->getOpPtrVector().push_back(
268 new OpCalculateHVecTensorField<SPACE_DIM, SPACE_DIM>(
269 eP.piolaStress, eP.dataAtPts->getApproxPAtPts(),
270 piola_scale_ptr));
271 side_fe_ptr->getOpPtrVector().push_back(
272 new OpCalculateHTensorTensorField<3, 3>(
273 eP.bubbleField, eP.dataAtPts->getApproxPAtPts(), MBMAXTYPE));
274 side_fe_ptr->getOpPtrVector().push_back(
276 side_fe_ptr->getOpPtrVector().push_back(
277 new OpCalculateVectorFieldValues<3>(
278 eP.spatialL2Disp, eP.dataAtPts->getSmallWL2AtPts(), MBTET));
279 for (auto &[name, ents, reaction_vec] : reactionForces) {
280 face_fe->getOpPtrVector().push_back(new OpCalculateReactionForces(
281 eP.dataAtPts, name, ents, reaction_vec));
282 }
283
284 CHKERR eP.mField.loop_finite_elements(problemPtr->getName(),
285 eP.skinElement, *face_fe);
286
287 for (auto &[name, ents, reaction_vec] : reactionForces) {
288 std::array<double, 6> block_reaction_force{0.0, 0.0, 0.0,
289 0.0, 0.0, 0.0};
290
291 MPI_Allreduce(reaction_vec.data(), &block_reaction_force, 6, MPI_DOUBLE,
292 MPI_SUM, eP.mField.get_comm());
293
294 for (auto &force : block_reaction_force) {
295 if (std::abs(force) < 1e-12) {
296 force = 0.0;
297 }
298 }
299 MOFEM_LOG_C("EP", Sev::inform,
300 "Step %d time %3.4g Block %s Reaction force [%3.6e, "
301 "%3.6e, %3.6e]",
302 ts_step, ts_t, name.c_str(), block_reaction_force[0],
303 block_reaction_force[1], block_reaction_force[2]);
304 MOFEM_LOG_C("EP", Sev::inform,
305 "Step %d time %3.4g Block %s Moment [%3.6e, %3.6e, %3.6e]",
306 ts_step, ts_t, name.c_str(), block_reaction_force[3],
307 block_reaction_force[4], block_reaction_force[5]);
308 }
310 };
311
312 auto save_restart_file = [&]() {
314
315 if (writeRestart && ts_u) {
316 PetscViewer viewer;
317 CHKERR PetscViewerBinaryOpen(
318 PETSC_COMM_WORLD, ("restart_" + get_step(ts_step) + ".dat").c_str(),
319 FILE_MODE_WRITE, &viewer);
320 CHKERR PetscViewerBinarySetSkipInfo(viewer, PETSC_TRUE);
321 CHKERR VecView(ts_u, viewer);
322 CHKERR PetscViewerDestroy(&viewer);
323 }
325 };
326
327 auto post_process_skin = [&]() {
329 CHKERR eP.postProcessResults(1, "out_sol_elastic_" + get_step(ts_step) +
330 ".h5m");
332 };
333
334 auto post_process_material_forces = [&]() {
336 // Function to get material force tags
339 };
340
341 auto post_process_skeleton_results = [&]() {
343 auto get_material_force_tags = [&]() {
344 auto &moab = eP.mField.get_moab();
345 std::vector<Tag> tag(2);
346 CHK_MOAB_THROW(moab.tag_get_handle("MaterialForce", tag[0]),
347 "can't get tag");
348 CHK_MOAB_THROW(moab.tag_get_handle("FacePressure", tag[1]),
349 "can't get tag");
350 return tag;
351 };
352 // Post-process skeleton elements
353 bool post_process_skeleton = true;
354#ifndef NDEBUG
355 post_process_skeleton = true;
356#endif
357 if (post_process_skeleton) {
359 1, "out_skeleton_" + get_step(ts_step) + ".h5m", PETSC_NULLPTR,
360 get_material_force_tags());
361 }
363 };
364
365 auto calculate_energy = [&]() {
367 // Loop boundary elements with traction boundary conditions
368 *gEnergy = 0;
369 TetPolynomialBase::switchCacheBaseOn<HDIV>({volPostProcEnergy.get()});
370 CHKERR eP.mField.loop_finite_elements(problemPtr->getName(), "EP",
372 TetPolynomialBase::switchCacheBaseOff<HDIV>({volPostProcEnergy.get()});
373
374 double body_energy = 0;
375 MPI_Allreduce(gEnergy.get(), &body_energy, 1, MPI_DOUBLE, MPI_SUM,
376 eP.mField.get_comm());
378 auto crack_area_ptr = boost::make_shared<double>(0.0);
379 CHKERR eP.calculateCrackArea(crack_area_ptr);
380 MOFEM_LOG_C("EP", Sev::inform,
381 "Step %d time %3.4g strain energy %3.6e crack "
382 "area %3.6e",
383 ts_step, ts_t, body_energy, *crack_area_ptr);
384 if (eP.mField.get_comm_rank() == 0) {
385 auto crack_faces = *eP.crackFaces;
386 if (crack_faces.empty()) {
387 crack_faces = *eP.frontEdges;
388 }
391 (boost::format("crack_faces_step_%d.vtk") % ts_step).str(),
392 crack_faces);
395 (boost::format("front_edges_step_%d.vtk") % ts_step).str(),
396 *eP.frontEdges);
397 }
398
399 } else {
400 MOFEM_LOG_C("EP", Sev::inform, "Step %d time %3.4g strain energy %3.6e",
401 ts_step, ts_t, body_energy);
402 }
404 };
405
406 auto post_proc_at_points = [&](std::array<double, 3> point,
407 std::string str) {
409
410 dataFieldEval->setEvalPoints(point.data(), point.size() / 3);
411
412 struct OpPrint : public VolOp {
413
414 EshelbianCore &eP;
415 std::array<double, 3> point;
416 std::string str;
417
418 OpPrint(EshelbianCore &ep, std::array<double, 3> &point,
419 std::string &str)
420 : VolOp(ep.spatialL2Disp, VolOp::OPROW), eP(ep), point(point),
421 str(str) {}
422
423 MoFEMErrorCode doWork(int side, EntityType type,
424 EntitiesFieldData::EntData &data) {
426 if (type == MBTET) {
427 if (getGaussPts().size2()) {
428
429 auto t_h = getFTensor2FromMat<3, 3>(eP.dataAtPts->hAtPts);
430 auto t_approx_P =
431 getFTensor2FromMat<3, 3>(eP.dataAtPts->approxPAtPts);
432
433 FTensor::Index<'i', 3> i;
434 FTensor::Index<'j', 3> j;
435 FTensor::Index<'k', 3> k;
436 const double jac = determinantTensor3by3(t_h);
438 t_cauchy(i, j) = (1. / jac) * (t_approx_P(i, k) * t_h(j, k));
439
440 auto add = [&]() {
441 std::ostringstream s;
442 s << str << " elem " << getFEEntityHandle() << " ";
443 return s.str();
444 };
445
446 auto print_tensor = [](auto &t) {
447 std::ostringstream s;
448 s << t;
449 return s.str();
450 };
451
452 std::ostringstream print;
453 MOFEM_LOG("EPSYNC", Sev::inform)
454 << add() << "comm rank " << eP.mField.get_comm_rank();
455 MOFEM_LOG("EPSYNC", Sev::inform)
456 << add() << "point " << getVectorAdaptor(point.data(), 3);
457 MOFEM_LOG("EPSYNC", Sev::inform)
458 << add() << "coords at gauss pts " << getCoordsAtGaussPts();
459 MOFEM_LOG("EPSYNC", Sev::inform)
460 << add() << "w " << eP.dataAtPts->wL2AtPts;
461 MOFEM_LOG("EPSYNC", Sev::inform)
462 << add() << "Piola " << eP.dataAtPts->approxPAtPts;
463 MOFEM_LOG("EPSYNC", Sev::inform)
464 << add() << "Cauchy " << print_tensor(t_cauchy);
465 }
466 }
468 }
469 };
470
471 if (auto fe_ptr = dataFieldEval->feMethodPtr.lock()) {
472
473 fe_ptr->getOpPtrVector().push_back(new OpPrint(eP, point, str));
474 CHKERR eP.mField.getInterface<FieldEvaluatorInterface>()
475 ->evalFEAtThePoint<SPACE_DIM>(
476 point.data(), 1e-12, problemPtr->getName(), "EP", dataFieldEval,
478 MF_EXIST, QUIET);
479 fe_ptr->getOpPtrVector().pop_back();
480 }
481
483 };
484
485 auto check_external_strain = [&](std::array<double, 3> point,
486 std::string str, PetscInt atom_test) {
488
489 dataFieldEval->setEvalPoints(point.data(), point.size() / 3);
490 auto vec = createVectorMPI(eP.mField.get_comm(), 1, 1);
491 eP.dataAtPts->wL2AtPts.resize(0, 0, false);
492
493 if (auto fe_ptr = dataFieldEval->feMethodPtr.lock()) {
494 CHKERR eP.mField.getInterface<FieldEvaluatorInterface>()
495 ->evalFEAtThePoint<SPACE_DIM>(
496 point.data(), 1e-12, problemPtr->getName(), "EP", dataFieldEval,
498 MF_EXIST, QUIET);
499 }
500 if (eP.dataAtPts->wL2AtPts.size2() == 0) {
501 VecSetValue(vec, 0, 0.0, ADD_VALUES);
502 } else if (eP.dataAtPts->wL2AtPts.size2() == 1) {
503 auto add = [&]() {
504 std::ostringstream s;
505 s << str << " elem " << getFEEntityHandle() << " ";
506 return s.str();
507 };
508 MOFEM_LOG("EPSYNC", Sev::inform)
509 << add() << "comm rank " << eP.mField.get_comm_rank();
510 MOFEM_LOG("EPSYNC", Sev::inform)
511 << add() << "point " << getVectorAdaptor(point.data(), 3);
512 MOFEM_LOG("EPSYNC", Sev::inform)
513 << add() << "w " << eP.dataAtPts->wL2AtPts;
514 double disp_at_point = eP.dataAtPts->wL2AtPts(0);
515 VecSetValue(vec, 0, disp_at_point, ADD_VALUES);
516 }
517
519 if (ts_t > 0.0) {
520 CHKERR VecAssemblyBegin(vec);
521 CHKERR VecAssemblyEnd(vec);
522 double error;
523 PetscInt idx = 0;
524 if (eP.mField.get_comm_rank() == 0) {
525 CHKERR VecGetValues(vec, 1, &idx, &error);
526 }
527 MPI_Bcast(&error, 1, MPI_DOUBLE, 0, PETSC_COMM_WORLD);
528 // Check if the displacement is correct for the applied external
529 // strain For the bar problem, we expect a displacement of 0.25 at
530 // point A
531 if (std::abs(error - 0.25) > 1e-5) {
532 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
533 "Atom test %d failed: wrong displacement.", atom_test);
534 }
535 }
537 };
538
539 CHKERR save_restart_file();
540 CHKERR calculate_reaction_forces();
541 CHKERR calculate_energy();
542
543 CHKERR post_process_material_forces();
544 CHKERR post_process_skin();
545 CHKERR post_process_skeleton_results();
546
547 PetscInt atom_test = 0;
548 CHKERR PetscOptionsGetInt(PETSC_NULLPTR, "", "-atom_test", &atom_test,
549 PETSC_NULLPTR);
550
551 if (atom_test == 14) {
552 // Points for external strain
553 for (auto &pts : ptsHashMap) {
554 CHKERR check_external_strain(pts.second, pts.first, atom_test);
555 }
556 } else {
557 // Points for Cook beam
558 for (auto &pts : ptsHashMap) {
559 CHKERR post_proc_at_points(pts.second, pts.first);
561 }
562 }
563
565 }
#define MOFEM_LOG_SEVERITY_SYNC(comm, severity)
Synchronise "SYNC" on curtain severity level.
#define MOFEM_LOG_C(channel, severity, format,...)
@ QUIET
@ MF_EXIST
@ MOFEM_ATOM_TEST_INVALID
Definition definitions.h:40
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition Templates.hpp:31
auto createVectorMPI(MPI_Comm comm, PetscInt n, PetscInt N)
Create MPI Vector.
MoFEMErrorCode determinantTensor3by3(T1 &t, T2 &det)
Calculate determinant 3 by 3.
constexpr double t
plate stiffness
Definition plate.cpp:58
static constexpr int approx_order
static PetscBool crackingOn
MoFEMErrorCode calculateFaceMaterialForce(const int tag, TS ts)
MoFEMErrorCode calculateCrackArea(boost::shared_ptr< double > area_ptr)
MoFEMErrorCode postProcessSkeletonResults(const int tag, const std::string file, Vec f_residual=PETSC_NULLPTR, std::vector< Tag > tags_to_transfer={})
const std::string skinElement
MoFEMErrorCode postProcessResults(const int tag, const std::string file, Vec f_residual=PETSC_NULLPTR, Vec var_vec=PETSC_NULLPTR, std::vector< Tag > tags_to_transfer={})
boost::shared_ptr< Range > crackFaces
boost::shared_ptr< Range > frontEdges
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
auto save_range

◆ preProcess()

MoFEMErrorCode EshelbianMonitor::preProcess ( )
inline

Definition at line 175 of file EshelbianMonitor.cpp.

175{ return 0; }

Member Data Documentation

◆ dataFieldEval

boost::shared_ptr<SetPtsData> EshelbianMonitor::dataFieldEval

Definition at line 17 of file EshelbianMonitor.cpp.

◆ eP

EshelbianCore& EshelbianMonitor::eP

Definition at line 16 of file EshelbianMonitor.cpp.

◆ gEnergy

boost::shared_ptr<double> EshelbianMonitor::gEnergy

Definition at line 19 of file EshelbianMonitor.cpp.

◆ ptsHashMap

PtsHashMap EshelbianMonitor::ptsHashMap

Definition at line 20 of file EshelbianMonitor.cpp.

◆ volPostProcEnergy

boost::shared_ptr<VolEle> EshelbianMonitor::volPostProcEnergy

Definition at line 18 of file EshelbianMonitor.cpp.

◆ writeRestart

PetscBool EshelbianMonitor::writeRestart = PETSC_FALSE

Definition at line 21 of file EshelbianMonitor.cpp.


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