v0.14.0
Loading...
Searching...
No Matches
CrackPropagation.hpp
Go to the documentation of this file.
1/** \file CrackPropagation.hpp
2 \brief Main class for crack propagation
3*/
4
5/* This file is part of MoFEM.
6 * MoFEM is free software: you can redistribute it and/or modify it under
7 * the terms of the GNU Lesser General Public License as published by the
8 * Free Software Foundation, either version 3 of the License, or (at your
9 * option) any later version.
10 *
11 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
12 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14 * License for more details.
15 *
16 * You should have received a copy of the GNU Lesser General Public
17 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
18
19#ifndef __CRACK_PROPAGATION_HPP__
20#define __CRACK_PROPAGATION_HPP__
21
22namespace FractureMechanics {
23
24MoFEMErrorCode clean_pcomms(moab::Interface &moab,
25 boost::shared_ptr<WrapMPIComm> moab_comm_wrap);
26
27MoFEMErrorCode broadcast_entities(moab::Interface &moab,
28 moab::Interface &moab_tmp,
29 const int from_proc, Range &entities,
30 const bool adjacencies, const bool tags);
31
32/**
33 * @brief Tapes numbers used by ADOL-C.
34 *
35 */
48};
49
50/**
51 * @brief Names of tangent matrices tests
52 *
53 */
59};
60/**
61 * Crack propagation data and methods
62 * \nosubgrouping
63 *
64 * - Functions strarting with build, for example buildElastic,
65 * add and build fields.
66 *
67 * - Functions starting with declare, for example declareElasticFE,
68 * add and build finite elements.
69 *
70 * - Functions starting with assemble, for example assembleElasticDM,
71 * create finite element instances and operators.
72 *
73 * - Functions starting with create, for example createElasticDM, create
74 * problems.
75 *
76 */
78
79 static bool debug;
80
81 static bool parallelReadAndBroadcast; ///< That is unstable development, for
82 ///< some meshses (propably generated by
83 ///< cubit) this is not working. Error
84 ///< can be attributed to bug in MOAB.
85
86 Version fileVersion;
87
88 MoFEMErrorCode getInterfaceVersion(Version &version) const {
90 version = Version(FM_VERSION_MAJOR, FM_VERSION_MINOR, FM_VERSION_BUILD);
92 }
93
94 /**
95 * \brief Getting interface of core database
96 * @param uuid unique ID of interface that can be either for
97 * CrackPropagation, CPSolvers or CPMeshCut interface
98 * @param iface returned pointer to interface
99 * @return error code
100 */
101 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
102 UnknownInterface **iface) const;
103
104 boost::shared_ptr<WrapMPIComm> moabCommWorld;
105
110
111 double gC; ///< Griffith energy
112 double griffithE; ///< Griffith stability parameter
113 double griffithR; ///< Griffith regularisation parameter
114
115 double betaGc; ///< heterogeneous Griffith energy exponent
116 PetscBool isGcHeterogeneous; ///< flag for heterogeneous gc
117
118 PetscBool isPartitioned;
119 PetscBool propagateCrack; ///< If true crack propagation is calculated
122 PetscBool addSingularity;
124
125 PetscBool isPressureAle; ///< If true surface pressure is considered in ALE
126 PetscBool areSpringsAle; ///< If true surface spring is considered in ALE
127 PetscBool
128 isSurfaceForceAle; ///< If true surface pressure is considered in ALE
129 PetscBool
130 ignoreMaterialForce; ///< If true surface pressure is considered in ALE
131
132 std::string mwlsApproxFile; ///< Name of file with internal stresses
133 std::string mwlsStressTagName; ///< Name of tag with internal stresses
134 std::string mwlsEigenStressTagName; ///< Name of tag with eigen stresses
135 std::string mwlsRhoTagName; ///< Name of tag with density
136 int residualStressBlock; ///< Block on which residual stress is applied
137 int densityMapBlock; ///< Block on which density is applied (Note: if used
138 ///< with internal stress, it has to be the same ID!)
139 std::string configFile;
141
142 double
143 partitioningWeightPower; ///< Weight controling load balance. Elements at
144 ///< the crack front are higher order, also have
145 ///< more fildes, associated with ALE
146 ///< formulation. Higher number put more weight
147 ///< on those elements while partitioning.
148
149 PetscBool resetMWLSCoeffsEveryPropagationStep; ///< If true, MWLS coefficients
150 ///< are recalulated every
151 ///< propagaton step.
152
154
156
157 );
158
159 PetscBool
160 addAnalyticalInternalStressOperators; ///< Add analytical internal stress
161 ///< operators for MWLS test
162
163 PetscBool solveEigenStressProblem; ///< Solve eigen problem
164 PetscBool useEigenPositionsSimpleContact; ///< Use eigen positions for
165 ///< matching meshes contact
166
167 int postProcLevel; ///< level of postprocessing (amount of output files)
169
172
173 /**
174 * Map associating string literals to bit refinement levels
175 * The following convention is currently used:
176 * "mesh_cut" : level 0
177 * "spatial_domain" : level 1
178 * "material_domain" : level 2
179 * */
180 std::map<std::string, BitRefLevel> mapBitLevel;
181
182 /**
183 * \brief Constructor of Crack Propagation, prints fracture module version and
184 registers the three interfaces
185 * used to solve fracture problem, i.e. CrackPropagation, CPSolvers and
186 * CPMeshCut
187 * @param m_field interface to MoAB database
188 * @param approx_order order of approximation space, default is 2
189 * @param geometry_order order of geometry, default is 1 (can be up to 2)
190 */
191 CrackPropagation(MoFEM::Interface &m_field, const int approx_order = 2,
192 const int geometry_order = 1);
193
194 virtual ~CrackPropagation();
195
196 /**
197 * \brief Get options form command line
198 * @return error code
199 */
200 MoFEMErrorCode getOptions();
201
202 /**
203 * \brief This is run with ctest
204 * @return error code
205 */
206 MoFEMErrorCode tetsingReleaseEnergyCalculation();
207
208 inline int &getNbLoadSteps() { return nbLoadSteps; }
209 inline int getNbLoadSteps() const { return nbLoadSteps; }
210 inline double &getLoadScale() { return loadScale; }
211 inline double getLoadScale() const { return loadScale; }
212 inline int &getNbCutSteps() { return nbCutSteps; }
213 inline int getNbCutSteps() const { return nbCutSteps; }
214
215 /// \brief read mesh file
216 MoFEMErrorCode readMedFile();
217
218 /**
219 * @brief Save entities on ech processor.
220 *
221 * Usually used for debugging, to check if meshes are consistent on all
222 * processors.
223 *
224 * @param prefix file name prefix
225 * @param ents entities to save
226 * @return MoFEMErrorCode error code
227 */
228 MoFEMErrorCode saveEachPart(const std::string prefix, const Range &ents);
229
230 /** \name Resolve shared entities and partition mesh */
231
232 /**@{*/
233
234 /// \brief partotion mesh
235 MoFEMErrorCode partitionMesh(BitRefLevel bit1, BitRefLevel bit2,
236 int verb = QUIET, const bool debug = false);
237
238 /// \brief resolve shared entities
239 MoFEMErrorCode resolveShared(const Range &tets, Range &proc_ents,
240 const int verb = QUIET,
241 const bool debug = false);
242
243 /// \brief resole shared entities by bit level
244 MoFEMErrorCode resolveSharedBitRefLevel(const BitRefLevel bit,
245 const int verb = QUIET,
246 const bool debug = false);
247
248 /**@}*/
249
250 /** \name Set bit level for material problem */
251
252 /**@{*/
253
254 /**
255 * @brief Set bit ref level for entities adjacent to crack front
256 *
257 * @param from_bit
258 * @param bit
259 * @param nb_levels adjacency bridge levels
260 * @return MoFEMErrorCode
261 */
262 MoFEMErrorCode setCrackFrontBitLevel(BitRefLevel from_bit, BitRefLevel bit,
263 const int nb_levels = 2,
264 const bool debug = false);
265
266 /**@}*/
267
268 /** \name Build problem */
269
270 /**@{*/
271
272 /**
273 * \brief Build problem fields
274 * @param bit1 bit level for computational domain
275 * @param mask1 bit mask
276 * @param bit2 bit level for elements near crack front
277 * @param mask2 bit mask 2
278 * @param verb verbosity level
279 * @param debug
280 * @return error code
281 *
282 */
283 MoFEMErrorCode buildProblemFields(const BitRefLevel &bit1,
284 const BitRefLevel &mask1,
285 const BitRefLevel &bit2,
286 const int verb = QUIET,
287 const bool debug = false);
288
289 /**
290 * \brief Build problem finite elements
291 * @param bit1 bit level for computational domain
292 * @param bit2 bit level for elements near crack front
293 * @param surface_ids meshsets Ids on body skin
294 * @param verb verbosity level
295 * @param debug
296 * @return error code
297 *
298 */
299 MoFEMErrorCode buildProblemFiniteElements(BitRefLevel bit1, BitRefLevel bit2,
300 std::vector<int> surface_ids,
301 const int verb = QUIET,
302 const bool debug = false);
303
304 /**
305 * \brief Construct problem data structures
306 * @param surface_ids surface_ids meshsets Ids on body skin
307 * @param verb verbosity level
308 * @param debug
309 * @return error code
310 */
311 MoFEMErrorCode createProblemDataStructures(const std::vector<int> surface_ids,
312 const int verb = QUIET,
313 const bool debug = false);
314
315 /**
316 * \brief Crate DMs for all problems and sub problems
317 * @param dm_elastic elastic DM
318 * @param dm_eiegn_elastci eigen_elastic DM to solve proble wih eigen
319 * stresses
320 * @param dm_material material DM
321 * @param dm_crack_propagation composition of elastic DM and material DM
322 * @param dm_material_forces used to calculate material forces, sub-dm of
323 * dm_crack_propagation
324 * @param dm_surface_projection dm to project material forces on surfaces,
325 * sub-dm of dm_crack_propagation
326 * @param dm_crack_srf_area dm to project material forces on crack
327 * surfaces, sub-dm of dm_crack_propagation
328 * @param surface_ids IDs of surface meshsets
329 * @param fe_surf_proj_list name of elements on surface elements
330 * @return error code
331 */
332 MoFEMErrorCode createDMs(SmartPetscObj<DM> &dm_elastic,
333 SmartPetscObj<DM> &dm_eigen_elastic,
334 SmartPetscObj<DM> &dm_material,
335 SmartPetscObj<DM> &dm_crack_propagation,
336 SmartPetscObj<DM> &dm_material_forces,
337 SmartPetscObj<DM> &dm_surface_projection,
338 SmartPetscObj<DM> &dm_crack_srf_area,
339 std::vector<int> surface_ids,
340 std::vector<std::string> fe_surf_proj_list);
341
342 MoFEMErrorCode deleteEntities(const int verb = QUIET,
343 const bool debug = false);
344
345 /**@}*/
346
347 /** \name Build fields */
348
349 /**@{*/
350
351 /**
352 * \brief Declate fields for elastic analysis
353 * @param bit bit refinement level
354 * @return error code
355 */
356 MoFEMErrorCode buildElasticFields(
357 const BitRefLevel bit, const BitRefLevel mask = BitRefLevel().set(),
358 const bool proc_only = true, const bool build_fields = true,
359 const int verb = QUIET, const bool debug = false);
360
361 /**
362 * \brief Declate field for arc-length
363 * @param bit bit refinement level
364 * @return error code
365 */
366 MoFEMErrorCode buildArcLengthField(const BitRefLevel bit,
367 const bool build_fields = true,
368 const int verb = QUIET);
369
370 /// \brief build fields with Lagrange multipliers to constrain surfaces
371 MoFEMErrorCode buildSurfaceFields(const BitRefLevel bit,
372 const bool proc_only = true,
373 const bool build_fields = true,
374 const int verb = QUIET,
375 const bool debug = false);
376
377 /// \brief build fields with Lagrange multipliers to constrain edges
378 MoFEMErrorCode buildEdgeFields(const BitRefLevel bit,
379 const bool proc_only = true,
380 const bool build_fields = true,
381 const int verb = QUIET,
382 const bool debug = false);
383
384 /// \brief declare crack surface files
385 MoFEMErrorCode buildCrackSurfaceFieldId(const BitRefLevel bit,
386 const bool proc_only = true,
387 const bool build_fields = true,
388 const int verb = QUIET,
389 const bool debug = false);
390
391 /// \brief declare crack surface files
392 MoFEMErrorCode buildCrackFrontFieldId(const BitRefLevel bit,
393 const bool build_fields = true,
394 const int verb = QUIET,
395 const bool debug = false);
396
397 /**
398 * \brief Lagrange multipliers field which constrains material displacements
399 *
400 * Nodes on both sides are constrainded to have the same material
401 * displacements
402 *
403 */
404 MoFEMErrorCode buildBothSidesFieldId(const BitRefLevel bit_spatial,
405 const BitRefLevel bit_material,
406 const bool proc_only = false,
407 const bool build_fields = true,
408 const int verb = QUIET,
409 const bool debug = false);
410
411 /** \brief Zero fields with lagrange multipliers
412 */
413 MoFEMErrorCode zeroLambdaFields();
414
415 /**@}*/
416
417 /** \name Build finite elements */
418
419 /**@{*/
420
421 /// \brief declare elastic finite elements
422 MoFEMErrorCode
423 declareElasticFE(const BitRefLevel bit1, const BitRefLevel mask1,
424 const BitRefLevel bit2, const BitRefLevel mask2,
425 const bool add_forces = true, const bool proc_only = true,
426 const int verb = QUIET);
427
428 /// \brief create arc-length element entity and declare elemets
429 MoFEMErrorCode declareArcLengthFE(const BitRefLevel bits,
430 const int verb = QUIET);
431
432 MoFEMErrorCode
433 declareExternalForcesFE(const BitRefLevel bit,
434 const BitRefLevel mask = BitRefLevel().set(),
435 const bool proc_only = true);
436
437 /** \brief Declare FE for pressure BC in ALE formulation (in material domain)
438 * @param bit material domain bit ref level
439 * @param mask bit ref level mask for the problem
440 * @return error code
441 */
442 MoFEMErrorCode
443 declareSurfaceForceAleFE(const BitRefLevel bit,
444 const BitRefLevel mask = BitRefLevel().set(),
445 const bool proc_only = true);
446
447 /** \brief Declare FE for pressure BC in ALE formulation (in material domain)
448 * @param bit material domain bit ref level
449 * @param mask bit ref level mask for the problem
450 * @return error code
451 */
452 MoFEMErrorCode
453 declarePressureAleFE(const BitRefLevel bit,
454 const BitRefLevel mask = BitRefLevel().set(),
455 const bool proc_only = true);
456
457 /** \brief Declare FE for spring BC in ALE formulation (in material domain)
458 * @param bit material domain bit ref level
459 * @param mask bit ref level mask for the problem
460 * @return error code
461 */
462 MoFEMErrorCode
463 declareSpringsAleFE(const BitRefLevel bit,
464 const BitRefLevel mask = BitRefLevel().set(),
465 const bool proc_only = true);
466
467 /** \brief Declare FE for pressure BC in ALE formulation (in material domain)
468 * @param bit material domain bit ref level
469 * @param mask bit ref level mask for the problem
470 * @return error code
471 */
472 MoFEMErrorCode
473 declareSimpleContactAleFE(const BitRefLevel bit,
474 const BitRefLevel mask = BitRefLevel().set(),
475 const bool proc_only = true);
476
477 /// \brief declare material finite elements
478 MoFEMErrorCode declareMaterialFE(const BitRefLevel bit,
479 const BitRefLevel mask = BitRefLevel().set(),
480 const bool proc_only = true,
481 const bool verb = QUIET);
482
483 MoFEMErrorCode declareFrontFE(const BitRefLevel bit,
484 const BitRefLevel mask = BitRefLevel().set(),
485 const bool proc_only = true,
486 const bool verb = QUIET);
487
488 MoFEMErrorCode
489 declareBothSidesFE(const BitRefLevel bit_spatial,
490 const BitRefLevel bit_material,
491 const BitRefLevel mask = BitRefLevel().set(),
492 const bool proc_only = true, const bool verb = QUIET);
493
494 /// \brief declare mesh smoothing finite elements
495 MoFEMErrorCode
496 declareSmoothingFE(const BitRefLevel bit,
497 const BitRefLevel mask = BitRefLevel().set(),
498 const bool proc_only = true, const bool verb = QUIET);
499
500 /// \brief declare surface sliding elements
501 MoFEMErrorCode declareSurfaceFE(std::string fe_name, const BitRefLevel bit,
502 const BitRefLevel mask,
503 const std::vector<int> &ids,
504 const bool proc_only = true,
505 const int verb = QUIET,
506 const bool debug = false);
507
508 // \brief declare edge sliding elements
509 MoFEMErrorCode declareEdgeFE(std::string fe_name, const BitRefLevel bit,
510 const BitRefLevel mask,
511 const bool proc_only = true,
512 const int verb = QUIET,
513 const bool debug = false);
514
515 MoFEMErrorCode
516 declareCrackSurfaceFE(std::string fe_name, const BitRefLevel bit,
517 const BitRefLevel mask, const bool proc_only = true,
518 const int verb = QUIET, const bool debug = false);
519
520 /**@}*/
521
522 /** \name Create DMs & problems */
523
524 /**@{*/
525
526 /** \brief Create DM by composition of elastic DM and material DM
527 * @param dm composition of elastic DM and material DM
528 * @param prb_name problem name
529 * @param dm_elastic elastic DM
530 * @param dm_material material DM
531 * @param bit domain bit ref level
532 * @param mask mask ref level for problem
533 * @return error code
534 */
535 MoFEMErrorCode
536 createCrackPropagationDM(SmartPetscObj<DM> &dm, const std::string prb_name,
537 SmartPetscObj<DM> dm_elastic,
538 SmartPetscObj<DM> dm_material, const BitRefLevel bit,
539 const BitRefLevel mask,
540 const std::vector<std::string> fe_list);
541
542 /** \brief Create elastic problem DM
543 * @param dm elastic DM
544 * @param prb_name problem name
545 * @param bit domain bit ref level
546 * @param mask mask ref level for problem
547 * @return error code
548 */
549 MoFEMErrorCode createElasticDM(SmartPetscObj<DM> &dm,
550 const std::string prb_name,
551 const BitRefLevel bit,
552 const BitRefLevel mask = BitRefLevel().set());
553
554 /** \brief Create elastic problem DM
555 * @param dm elastic DM
556 * @param prb_name problem name
557 * @param bit domain bit ref level
558 * @param mask mask ref level for problem
559 * @return error code
560 */
561 MoFEMErrorCode
562 createEigenElasticDM(SmartPetscObj<DM> &dm, const std::string prb_name,
563 const BitRefLevel bit,
564 const BitRefLevel mask = BitRefLevel().set());
565
566 /** \brief Create problem to calculate boundary conditions
567 * @param dm newly created boundary condition DM
568 * @param prb_name problem name
569 * @param bit auxiliary bit ref level
570 * @param mask mask ref level for problem
571 * @return error code
572 * \note BcDM is only useed for testing, to enforce boundary conditions
573staisifing anlylitical solution to which numerical solution is compared.
574 */
575 MoFEMErrorCode createBcDM(SmartPetscObj<DM> &dm, const std::string prb_name,
576 const BitRefLevel bit,
577 const BitRefLevel mask = BitRefLevel().set());
578
579 /** \brief Create DM fto calculate material problem
580 * @param dm used to calculate material forces, sub-dm of
581 * dm_crack_propagation
582 * @param prb_name problem name
583 * @param bit material domain
584 * @param mask dof mask ref level for problem
585 * @param fe_list name of elements on surface elements
586 * @param debug flag for debugging
587 * @return error code
588 */
589 MoFEMErrorCode createMaterialDM(SmartPetscObj<DM> &dm,
590 const std::string prb_name,
591 const BitRefLevel bit, const BitRefLevel mask,
592 const std::vector<std::string> fe_list,
593 const bool debug = false);
594
595 /** \brief Create DM for calculation of material forces (sub DM of DM
596 * material)
597 * @param dm used to calculate material forces, sub-dm of
598 * dm_crack_propagation
599 * @param prb_name problem name
600 * @param verb compilation parameter determining the amount
601 * of information printed
602 * @return error code
603 */
604 MoFEMErrorCode createMaterialForcesDM(SmartPetscObj<DM> &dm,
605 SmartPetscObj<DM> dm_material,
606 const std::string prb_name,
607 const int verb = QUIET);
608
609 /** \brief create DM to calculate projection matrices (sub DM of DM
610 * material)
611 * @param dm DM to project material forces on surfaces,
612 * sub-dm of dm_crack_propagation
613 * @param dm_material material DM
614 * @param prb_name problem name
615 * @param surface_ids IDs of surface meshsets
616 * @param fe_list name of elements on surface elements
617 * @param verb compilation parameter determining the amount
618 * of information printed
619 * @return error code
620 */
621 MoFEMErrorCode createSurfaceProjectionDM(
622 SmartPetscObj<DM> &dm, SmartPetscObj<DM> dm_material,
623 const std::string prb_name, const std::vector<int> surface_ids,
624 const std::vector<std::string> fe_list, const int verb = QUIET);
625
626 /** \brief create DM to calculate Griffith energy
627 * @param dm dm to project material forces on crack
628 * surfaces, sub-dm of dm_crack_propagation
629 * @param dm_material material DM
630 * @param prb_name problem name
631 * @param verb compilation parameter determining the amount
632 * of information printed
633 * @return error code
634 */
635 MoFEMErrorCode createCrackFrontAreaDM(SmartPetscObj<DM> &dm,
636 SmartPetscObj<DM> dm_material,
637 const std::string prb_name,
638 const bool verb = QUIET);
639
640 /**@}*/
641
642 /** \name Create finite element instances */
643
644 /**@{*/
645
646 /** \brief create elastic finite element instance for spatial assembly
647 * @param mwls_stress_tag_name Name of internal stress tag
648 * @param close_crack if true, crack surface is closed
649 * @param verb compilation parameter determining the amount
650 * of information printed (not in use here)
651 * @param debug flag for debugging
652 * @return error code
653 */
654 MoFEMErrorCode assembleElasticDM(const std::string mwls_stress_tag_name,
655 const int verb = QUIET,
656 const bool debug = false);
657
658 /** \brief create elastic finite element instance for spatial assembly
659 * @param verb compilation parameter determining the amount
660 * of information printed (not in use here)
661 * @param debug flag for debugging
662 * @return error code
663 */
664 MoFEMErrorCode assembleElasticDM(const int verb = QUIET,
665 const bool debug = false);
666
667 // \brief add elastic element instances to SNES
668 MoFEMErrorCode addElasticFEInstancesToSnes(
669 DM dm, Mat m, Vec q, Vec f,
670 boost::shared_ptr<FEMethod> arc_method = boost::shared_ptr<FEMethod>(),
671 boost::shared_ptr<ArcLengthCtx> arc_ctx = nullptr, const int verb = QUIET,
672 const bool debug = false);
673
674 /// \brief create material element instance
675 MoFEMErrorCode assembleMaterialForcesDM(DM dm, const int verb = QUIET,
676 const bool debug = false);
677
678 /// \brief create smoothing element instance
679 MoFEMErrorCode assembleSmootherForcesDM(DM dm, const std::vector<int> ids,
680 const int verb = QUIET,
681 const bool debug = false);
682
683 /// \brief assemble coupling element instances
684 MoFEMErrorCode assembleCouplingForcesDM(DM dm, const int verb = QUIET,
685 const bool debug = false);
686
687 /// \brief add material elements instances to SNES
688 MoFEMErrorCode addMaterialFEInstancesToSnes(DM dm, const bool fix_crack_front,
689 const int verb = QUIET,
690 const bool debug = false);
691
692 /// \brief add softening elements instances to SNES
693 MoFEMErrorCode addSmoothingFEInstancesToSnes(DM dm,
694 const bool fix_crack_front,
695 const int verb = QUIET,
696 const bool debug = false);
697
698 /// \brief add finite element to SNES for crack propagation problem
699 MoFEMErrorCode
700 addPropagationFEInstancesToSnes(DM dm, boost::shared_ptr<FEMethod> arc_method,
701 boost::shared_ptr<ArcLengthCtx> arc_ctx,
702 const std::vector<int> &surface_ids,
703 const int verb = QUIET,
704 const bool debug = false);
705
706 /** \brief test LHS Jacobians
707 * @param bit bit level
708 * @param mask bitref mask
709 * @param test_nb name of the test
710 */
711 MoFEMErrorCode testJacobians(const BitRefLevel bit, const BitRefLevel mask,
712 tangent_tests test);
713
714 MoFEMErrorCode
715 addMWLSStressOperators(boost::shared_ptr<CrackFrontElement> &fe_rhs,
716 boost::shared_ptr<CrackFrontElement> &fe_lhs);
717
718 MoFEMErrorCode
719 addMWLSDensityOperators(boost::shared_ptr<CrackFrontElement> &fe_rhs,
720 boost::shared_ptr<CrackFrontElement> &fe_lhs);
721
722 /// \brief Update fixed nodes
723 MoFEMErrorCode updateMaterialFixedNode(const bool fix_front,
724 const bool fix_small_g,
725 const bool debug = false);
726
727 /**@}*/
728
729 /** \name Assemble and calculate */
730
731 /**@{*/
732
733 /// \brief assemble elastic part of matrix, by running elastic finite element
734 /// instance
735 MoFEMErrorCode calculateElasticEnergy(DM dm, const std::string msg = "");
736
737 /** \brief assemble material forces, by running material finite element
738 * instance
739 * @param dm dm used to calculate material forces, sub-dm
740 * of dm_crack_propagation
741 * @param q DoF sub-vector for crack propagation DM
742 * @param f Rhs sub-vector for crack propagation DM
743 * @param verb compilation parameter determining the amount
744 * of information printed
745 * @param debug flag for debugging
746 * @return error code
747 */
748 MoFEMErrorCode calculateMaterialForcesDM(DM dm, Vec q, Vec f,
749 const int verb = QUIET,
750 const bool debug = false);
751
752 /** \brief assemble smoothing forces, by running material finite element
753 * instance
754 * @param dm crack propagation dm that is composition of
755 * elastic DM and material DM
756 * @param q DoF sub-vector for crack propagation DM
757 * @param f Rhs sub-vector for crack propagation DM
758 * @param verb compilation parameter determining the amount
759 * of information printed
760 * @param debug flag for debugging
761 * @return error code
762 */
763 MoFEMErrorCode calculateSmoothingForcesDM(DM dm, Vec q, Vec f,
764 const int verb = QUIET,
765 const bool debug = false);
766
767 MoFEMErrorCode calculateSmoothingForceFactor(const int verb = QUIET,
768 const bool debug = true);
769
770 /// \brief assemble projection matrices
771 MoFEMErrorCode calculateSurfaceProjectionMatrix(DM dm_front, DM dm_project,
772 const std::vector<int> &ids,
773 const int verb = QUIET,
774 const bool debug = false);
775
776 /// \brief assemble crack front projection matrix (that constrains crack area
777 /// growth)
778 MoFEMErrorCode calculateFrontProjectionMatrix(DM dm_surface, DM dm_project,
779 const int verb = QUIET,
780 const bool debug = false);
781
782 /** \brief project material forces along the crack elongation direction
783 * @param dm_project material forces DM, sub-dm of
784 * dm_crack_propagation
785 * @param f_proj Rhs sub-vector of projected material forces
786 * for material forces DM
787 * @param f Rhs sub-vector for material forces DM
788 * @param verb compilation parameter determining the amount
789 * of information printed
790 * @param debug flag for debugging
791 * @return error code
792 */
793 MoFEMErrorCode projectMaterialForcesDM(DM dm_project, Vec f, Vec f_proj,
794 const int verb = QUIET,
795 const bool debug = false);
796
797 /// \brief calculate Griffith (driving) force
798 MoFEMErrorCode calculateGriffithForce(DM dm, const double gc, Vec f_griffith,
799 const int verb = QUIET,
800 const bool debug = false);
801
802 /// \brief project Griffith forces
803 MoFEMErrorCode projectGriffithForce(DM dm, Vec f_griffith,
804 Vec f_griffith_proj,
805 const int verb = QUIET,
806 const bool debug = false);
807
808 /// \brief calculate release energy
809 MoFEMErrorCode calculateReleaseEnergy(DM dm, Vec f_material_proj,
810 Vec f_griffith_proj, Vec f_lambda,
811 const double gc, const int verb = QUIET,
812 const bool debug = true);
813
814 /** \name Solvers */
815
816 /**@{*/
817
818 /** \brief solve elastic problem
819 * @param dm elastic DM
820 * @param snes snes solver for elastic DM
821 * @param m stiffness sub-matrix for elastic DM
822 * @param q DoF sub-vector for elastic DM
823 * @param f Rhs sub-vector for material DM
824 * @param snes_set_up boolean flag wheather to determin snes solver
825 * @param shell_m stiffness sub-matrix for elastic DM
826 * @return error code
827 */
828 MoFEMErrorCode solveElasticDM(DM dm, SNES snes, Mat m, Vec q, Vec f,
829 bool snes_set_up, Mat *shell_m);
830
831 /** \brief solve crack propagation problem
832 * @param dm crack propagation DM that is composed of
833 * elastic DM and material DM
834 * @param dm_elastic elastic DM
835 * @param snes snes solver for crack propagation DM
836 * @param m stiffness sub-matrix for crack propagation DM
837 * @param q DoF sub-vector for crack propagation DM
838 * @param f Rhs sub-vector for crack propagation DM
839 * @param snes_set_up boolean flag wheather to determin snes solver
840 * @return error code
841 */
842 MoFEMErrorCode solvePropagationDM(DM dm, DM dm_elastic, SNES snes, Mat m,
843 Vec q, Vec f);
844
845 /**@}*/
846
847 /**@}*/
848
849 /** \name Postprocess results */
850
851 /**@{*/
852
853 MoFEMErrorCode savePositionsOnCrackFrontDM(DM dm, Vec q,
854 const int verb = QUIET,
855 const bool debug = false);
856
857 MoFEMErrorCode writeCrackFont(const BitRefLevel &bit, const int step);
858
859 /// \brief post-process results for elastic solution
860 MoFEMErrorCode postProcessDM(DM dm, const int step, const std::string fe_name,
861 const bool approx_internal_stress);
862
863 /**@}*/
864
865 /** \name Set data from mesh to vectors and back */
866
867 /**@{*/
868
869 /// \brief set field from node positions
870 MoFEMErrorCode setFieldFromCoords(const std::string field_name);
871
872 /// \brief set material field from nodes
873 MoFEMErrorCode setMaterialPositionFromCoords();
874
875 /// \brief set spatial field from nodes
876 MoFEMErrorCode setSpatialPositionFromCoords();
877
878 /// \brief set coordinates from field
879 MoFEMErrorCode
880 setCoordsFromField(const std::string field_name = "MESH_NODE_POSITIONS");
881
882 /// \brief set singular dofs (on edges adjacent to crack front) from geometry
883 MoFEMErrorCode setSingularDofs(const string field_name,
884 const int verb = QUIET);
885
886 /// \brief set singular dofs of material field (on edges adjacent to crack
887 /// front) from geometry
888 MoFEMErrorCode setSingularElementMatrialPositions(const int verb = QUIET);
889
890 /// \brief remove singularity
892
893 /// \brief set maetrial field order to one
895
896 /**@}*/
897
898 /** \name Auxiliary method */
899
900 /**@{*/
901
902 boost::shared_ptr<DofEntity> arcLengthDof;
903
904 /// \brief set pointer to arc-length DOF
905 MoFEMErrorCode getArcLengthDof();
906
907 inline boost::shared_ptr<ArcLengthCtx> &getArcCtx() { return arcCtx; }
908 inline boost::shared_ptr<ArcLengthCtx> &getEigenArcCtx() {
909 return arcEigenCtx;
910 }
911
912 inline boost::shared_ptr<FEMethod>
913 getFrontArcLengthControl(boost::shared_ptr<ArcLengthCtx> arc_ctx) {
914 boost::shared_ptr<FEMethod> arc_method(
917 griffithForceElement->commonData, "LAMBDA_CRACKFRONT_AREA",
918 arc_ctx));
919 return arc_method;
920 }
921
922 inline boost::shared_ptr<AnalyticalDirichletBC::DirichletBC>
925 }
926
927 boost::shared_ptr<MWLSApprox> &getMWLSApprox() { return mwlsApprox; }
928
929 /**@}*/
930
931 /**
932 * \brief operator to post-process results on crack front nodes
933 */
935
937 Vec F;
938 std::string tagName;
939 PostProcVertexMethod(MoFEM::Interface &m_field, Vec f, std::string tag_name)
940 : mField(m_field), F(f), tagName(tag_name) {}
941
942 Tag tH;
943 double *aRray;
944 MoFEMErrorCode preProcess();
945
946 MoFEMErrorCode postProcess();
947
948 MoFEMErrorCode operator()();
949 };
950
951 /**
952 * \brief Constrains material displacement on both sides
953 *
954 * Both sides of crack surface have the same material displacements.
955 *
956 */
958
959 /**
960 * @brief Construct a new Both Surface Constrains object
961 *
962 * @param m_field
963 * @param lambda_field_name lagrange multiplayers
964 * @param spatial_field_name constrained spatial field
965 */
967 MoFEM::Interface &m_field,
968 const std::string lambda_field_name = "LAMBDA_BOTH_SIDES",
969 const std::string spatial_field_name = "MESH_NODE_POSITIONS")
970 : mField(m_field), lambdaFieldName(lambda_field_name),
971 spatialFieldName(spatial_field_name), aLpha(1) {
972 ierr = getOptions();
973 CHKERRABORT(PETSC_COMM_WORLD, ierr);
974 }
975
976 MoFEMErrorCode getOptions() {
978 ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
979 "Get Booth surface constrains element scaling",
980 "none");
981 CHKERRQ(ierr);
982 CHKERR PetscOptionsScalar("-booth_side_alpha", "scaling parameter", "",
983 aLpha, &aLpha, PETSC_NULL);
984 ierr = PetscOptionsEnd();
985 CHKERRQ(ierr);
987 }
988
989 MoFEMErrorCode preProcess() { return 0; }
990 MoFEMErrorCode operator()();
991 MoFEMErrorCode postProcess() { return 0; }
992
993 double aLpha; ///! Scaling parameter
994 Range masterNodes; ///! Node on which lagrange multiplier is set
995
996 private:
998 const std::string lambdaFieldName;
999 const std::string spatialFieldName;
1000 };
1001
1002 /// \brief Determine face orientation
1007 BitRefLevel defaultBitLevel;
1008 FaceOrientation(bool crack_front, Range contact_faces,
1009 BitRefLevel bit_level)
1010 : useProjectionFromCrackFront(crack_front), contactFaces(contact_faces),
1011 defaultBitLevel(bit_level) {}
1012 virtual ~FaceOrientation() {}
1013 MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field,
1014 const EntityHandle face,
1015 const BitRefLevel bit);
1017 const FEMethod *fe_method_ptr) {
1019 EntityHandle face = fe_method_ptr->numeredEntFiniteElementPtr->getEnt();
1020 BitRefLevel bit = fe_method_ptr->problemPtr->getBitRefLevel();
1021 CHKERR getElementOrientation(m_field, face, bit);
1023 }
1025 const EntityHandle face) {
1029 }
1030 };
1031
1032 struct ArcLengthSnesCtx : public SnesCtx {
1033
1034 ArcLengthSnesCtx(MoFEM::Interface &m_field, const std::string problem_name,
1035 boost::shared_ptr<ArcLengthCtx> &arc_ptr)
1036 : SnesCtx(m_field, problem_name), arcPtr(arc_ptr.get()) {
1037 MOFEM_LOG("CPWorld", Sev::noisy) << "ArcLengthSnesCtx create";
1038 }
1040 MOFEM_LOG("CPWorld", Sev::noisy) << "ArcLengthSnesCtx destroyed";
1041 };
1042
1043 inline auto getArcPtr() { return arcPtr; }
1044 inline void setVecDiagM(Vec diag_m) { diagM = diag_m; }
1045 inline auto getVecDiagM() { return diagM; }
1046
1047 private:
1050 };
1051
1052 // private:
1053
1054 /** \name Pointers to finite element instances */
1055
1056 /**@{*/
1057
1058 bool onlyHooke; ///< True if only Hooke material is applied
1059 PetscBool onlyHookeFromOptions; ///< True if only Hooke material is applied
1060
1061 /// Integrate spatial stresses, matrix and vector (not only element but all
1062 /// data structure)
1063 boost::shared_ptr<NonlinearElasticElement> elasticFe;
1064 boost::shared_ptr<NonlinearElasticElement> materialFe;
1065 boost::shared_ptr<CrackFrontElement> feLhs; ///< Integrate elastic FE
1066 boost::shared_ptr<CrackFrontElement> feRhs; ///< Integrate elastic FE
1067 boost::shared_ptr<CrackFrontElement>
1068 feMaterialRhs; ///< Integrate material stresses, assemble vector
1069 boost::shared_ptr<CrackFrontElement>
1070 feMaterialLhs; ///< Integrate material stresses, assemble matrix
1071 boost::shared_ptr<CrackFrontElement> feEnergy; ///< Integrate energy
1072 boost::shared_ptr<ConstrainMatrixCtx>
1073 projSurfaceCtx; ///< Data structure to project on the body surface
1074 boost::shared_ptr<ConstrainMatrixCtx>
1075 projFrontCtx; ///< Data structure to project on crack front
1076 boost::shared_ptr<MWLSApprox>
1077 mwlsApprox; ///< Data struture for MWLS approximation
1078 /// Surface elment to calculate tractions in material space
1079 boost::shared_ptr<NeumannForcesSurface::MyTriangleFE>
1081 boost::shared_ptr<DirichletSpatialPositionsBc>
1082 spatialDirichletBc; ///< apply Dirichlet BC to sparial positions
1083 boost::shared_ptr<AnalyticalDirichletBC::DirichletBC>
1084 analyticalDirichletBc; ///< apply analytical Dirichlet BC to sparial
1085 ///< positions
1086
1087 boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
1088 surfaceForces; ///< assemble surface forces
1089 boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
1090 surfaceForceAle; ///< assemble surface pressure (ALE)
1091 boost::shared_ptr<NeumannForcesSurface::DataAtIntegrationPts>
1092 commonDataSurfaceForceAle; ///< common data at integration points (ALE)
1093 boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
1094 surfacePressure; ///< assemble surface pressure
1095 boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
1096 surfacePressureAle; ///< assemble surface pressure (ALE)
1097 boost::shared_ptr<NeumannForcesSurface::DataAtIntegrationPts>
1098 commonDataSurfacePressureAle; ///< common data at integration points (ALE)
1099 boost::shared_ptr<boost::ptr_map<string, EdgeForce>>
1100 edgeForces; ///< assemble edge forces
1101 boost::shared_ptr<boost::ptr_map<string, NodalForce>>
1102 nodalForces; ///< assemble nodal forces
1103 boost::shared_ptr<FEMethod> assembleFlambda; ///< assemble F_lambda vector
1104 boost::shared_ptr<FEMethod> zeroFlambda; ///< assemble F_lambda vector
1105
1106 boost::shared_ptr<CrackFrontElement>
1107 feCouplingElasticLhs; ///< FE instance to assemble coupling terms
1108 boost::shared_ptr<CrackFrontElement>
1109 feCouplingMaterialLhs; ///< FE instance to assemble coupling terms
1110
1111 boost::shared_ptr<Smoother> smootherFe;
1112 boost::shared_ptr<Smoother::MyVolumeFE>
1113 feSmootherRhs; ///< Integrate smoothing operators
1114 boost::shared_ptr<Smoother::MyVolumeFE>
1115 feSmootherLhs; ///< Integrate smoothing operators
1116 boost::shared_ptr<VolumeLengthQuality<double>> volumeLengthDouble;
1117 boost::shared_ptr<VolumeLengthQuality<adouble>> volumeLengthAdouble;
1118
1119 boost::shared_ptr<
1121 tangentConstrains; ///< Constrains crack front in tangent directiona
1122 boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
1124 boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
1126 boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
1128 map<int, boost::shared_ptr<SurfaceSlidingConstrains>> surfaceConstrain;
1129 map<int, boost::shared_ptr<EdgeSlidingConstrains>> edgeConstrains;
1130
1131 boost::shared_ptr<BothSurfaceConstrains> bothSidesConstrains;
1132 boost::shared_ptr<BothSurfaceConstrains> closeCrackConstrains;
1133 boost::shared_ptr<BothSurfaceConstrains> bothSidesContactConstrains;
1134 boost::shared_ptr<DirichletFixFieldAtEntitiesBc> fixMaterialEnts;
1135
1136 boost::shared_ptr<GriffithForceElement> griffithForceElement;
1137 boost::shared_ptr<GriffithForceElement::MyTriangleFE> feGriffithForceRhs;
1138 boost::shared_ptr<GriffithForceElement::MyTriangleFE> feGriffithForceLhs;
1139 boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrainsDelta>
1141 boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrains>
1143 boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrains>
1145
1146 boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
1147 analiticalSurfaceElement; ///< Finite element to integrate tractions given
1148 ///< by analytical formula
1149 boost::shared_ptr<FaceElementForcesAndSourcesCore> feSpringLhsPtr;
1150 boost::shared_ptr<FaceElementForcesAndSourcesCore> feSpringRhsPtr;
1151 /**@}*/
1152
1153 /** \name Entity ranges */
1154
1155 /**@{*/
1156
1168
1170 arcLengthVertex; ///< Vertex associated with dof for arc-length control
1171
1172 /**@}*/
1173
1174 /** \name Auxiliary data */
1175
1176 /**@{*/
1177
1179
1183
1184 double maxG1;
1185 double maxJ;
1187 map<EntityHandle, double> mapG1; ///< hashmap of g1 - release energy at nodes
1188 map<EntityHandle, double> mapG3; ///< hashmap of g3 - release energy at nodes
1189 map<EntityHandle, double> mapJ; ///< hashmap of J - release energy at nodes
1190 map<EntityHandle, VectorDouble3>
1191 mapMatForce; ///< hashmap of material force at nodes
1192 map<EntityHandle, VectorDouble3>
1193 mapGriffith; ///< hashmap of Griffith energy at nodes
1194
1195 map<EntityHandle, double> mapSmoothingForceFactor;
1196
1197 bool doSurfaceProjection; ///< If crack surface is immersed in the body do not
1198 ///< projection
1199
1200 double smootherAlpha; ///< Controls mesh smoothing
1201 double initialSmootherAlpha; ///< Controls mesh smoothing, set at start of
1202 ///< analysis
1203 double smootherGamma; ///< Controls mesh smoothing
1204
1205 /**@}*/
1206
1207 /** \name arc-length options */
1208
1209 /**@{*/
1210
1211 double arcAlpha;
1212 double arcBeta;
1213 double arcS;
1214
1215 boost::shared_ptr<ArcLengthCtx> arcCtx;
1216 boost::shared_ptr<ArcLengthCtx> arcEigenCtx;
1217
1218 /**@}*/
1219
1220 /** \name Contact parameters, entity ranges and data structures */
1221
1222 /**@{*/
1223
1224 int contactOrder; ///< Approximation order for spatial positions
1225 ///< used in contact elements
1226 int contactLambdaOrder; ///< Approximation order for field of Lagrange
1227 ///< multipliers used to enforce contact
1228 double rValue; ///< Parameter for regularizing absolute value
1229 ///< function in the complementarity function
1230 double cnValue; ///< Augmentation parameter in the complementarity function
1231 ///< used to enforce KKT constraints
1232
1233 PetscBool ignoreContact; ///< If set true, contact interfaces are ignored
1234 ///< and prism interfaces are not created
1235 PetscBool fixContactNodes; ///< If set true, contact nodes are fixed in
1236 ///< the material configuration
1238
1240
1244
1248
1249 boost::shared_ptr<ContactSearchKdTree::ContactCommonData_multiIndex>
1255
1257
1258 boost::shared_ptr<FaceElementForcesAndSourcesCore> feLhsSpringALE;
1259 boost::shared_ptr<FaceElementForcesAndSourcesCore> feLhsSpringALEMaterial;
1260 boost::shared_ptr<FaceElementForcesAndSourcesCore> feRhsSpringALEMaterial;
1261 boost::shared_ptr<FaceElementForcesAndSourcesCore>
1263 boost::shared_ptr<FaceElementForcesAndSourcesCore>
1265 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
1267
1268 boost::shared_ptr<SimpleContactProblem> contactProblem;
1269 boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1271 boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1273 boost::shared_ptr<SimpleContactProblem::CommonDataSimpleContact>
1275
1276 boost::shared_ptr<MortarContactProblem> mortarContactProblemPtr;
1277
1278 boost::shared_ptr<MortarContactProblem::MortarContactElement>
1280 boost::shared_ptr<MortarContactProblem::MortarContactElement>
1282 boost::shared_ptr<MortarContactProblem::CommonDataMortarContact>
1284 boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1286 boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1288 boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1290 boost::shared_ptr<SimpleContactProblem::CommonDataSimpleContact>
1292
1294 boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1296 boost::shared_ptr<MortarContactProblem::MortarContactElement>
1299 moab::Interface &contactPostProcMoab;
1300
1301 Range constrainedInterface; ///< Range of faces on the constrained interface
1302
1303 /**@}*/
1304
1305 /** \name Data for bone */
1306
1307 /**@{*/
1308
1309 double rHo0; ///< Reference density if bone is analyzed
1310 double nBone; ///< Exponent parameter in bone density
1311
1312 /**@}*/
1313
1314 /** \name Interfaces */
1315
1316 /**@{*/
1317
1318 const boost::scoped_ptr<UnknownInterface> cpSolversPtr;
1319 const boost::scoped_ptr<UnknownInterface> cpMeshCutPtr;
1320
1321 /**@}*/
1322
1323 /** \name Cut/Refine/Split options */
1324
1325 /**@{*/
1326
1327 PetscBool doCutMesh;
1330
1331 /**@}*/
1332};
1333
1334} // namespace FractureMechanics
1335
1336#endif // __CRACK_PROPAGATION_HPP__
static PetscErrorCode ierr
@ QUIET
Definition: definitions.h:208
#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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
FTensor::Index< 'm', SPACE_DIM > m
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
auto bit
set bit
constexpr int nb_levels
Definition: level_set.cpp:41
MoFEMErrorCode broadcast_entities(moab::Interface &moab, moab::Interface &moab_tmp, const int from_proc, Range &entities, const bool adjacencies, const bool tags)
MoFEMErrorCode clean_pcomms(moab::Interface &moab, boost::shared_ptr< WrapMPIComm > moab_comm_wrap)
AdolcTags
Tapes numbers used by ADOL-C.
tangent_tests
Names of tangent matrices tests.
constexpr auto field_name
static constexpr int approx_order
Store variables for ArcLength analysis.
ArcLengthSnesCtx(MoFEM::Interface &m_field, const std::string problem_name, boost::shared_ptr< ArcLengthCtx > &arc_ptr)
Constrains material displacement on both sides.
BothSurfaceConstrains(MoFEM::Interface &m_field, const std::string lambda_field_name="LAMBDA_BOTH_SIDES", const std::string spatial_field_name="MESH_NODE_POSITIONS")
Construct a new Both Surface Constrains object.
MoFEM::Interface & mField
! Node on which lagrange multiplier is set
FaceOrientation(bool crack_front, Range contact_faces, BitRefLevel bit_level)
MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const EntityHandle face)
MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const FEMethod *fe_method_ptr)
MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const EntityHandle face, const BitRefLevel bit)
operator to post-process results on crack front nodes
MoFEM::Interface & mField
Tag tH
Vec F
double * aRray
MoFEMErrorCode operator()()
function is run for every finite element
std::string tagName
MoFEMErrorCode postProcess()
function is run at the end of loop
MoFEMErrorCode preProcess()
function is run at the beginning of loop
PostProcVertexMethod(MoFEM::Interface &m_field, Vec f, std::string tag_name)
boost::shared_ptr< FaceElementForcesAndSourcesCore > feSpringRhsPtr
MoFEMErrorCode zeroLambdaFields()
Zero fields with lagrange multipliers.
boost::shared_ptr< BothSurfaceConstrains > bothSidesConstrains
MoFEMErrorCode solveElasticDM(DM dm, SNES snes, Mat m, Vec q, Vec f, bool snes_set_up, Mat *shell_m)
solve elastic problem
boost::shared_ptr< AnalyticalDirichletBC::DirichletBC > analyticalDirichletBc
MoFEMErrorCode declareSpringsAleFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
Declare FE for spring BC in ALE formulation (in material domain)
boost::shared_ptr< CrackFrontElement > feMaterialRhs
Integrate material stresses, assemble vector.
boost::shared_ptr< FaceElementForcesAndSourcesCore > feSpringLhsPtr
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfacePressure
assemble surface pressure
boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > contactOrientation
MoFEMErrorCode getArcLengthDof()
set pointer to arc-length DOF
MoFEMErrorCode declareSimpleContactAleFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
Declare FE for pressure BC in ALE formulation (in material domain)
MoFEMErrorCode createMaterialDM(SmartPetscObj< DM > &dm, const std::string prb_name, const BitRefLevel bit, const BitRefLevel mask, const std::vector< std::string > fe_list, const bool debug=false)
Create DM fto calculate material problem.
MoFEMErrorCode deleteEntities(const int verb=QUIET, const bool debug=false)
MoFEMErrorCode declareSmoothingFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
declare mesh smoothing finite elements
MoFEMErrorCode declareBothSidesFE(const BitRefLevel bit_spatial, const BitRefLevel bit_material, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
boost::shared_ptr< DirichletSpatialPositionsBc > spatialDirichletBc
apply Dirichlet BC to sparial positions
MoFEMErrorCode declareExternalForcesFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
MoFEMErrorCode solvePropagationDM(DM dm, DM dm_elastic, SNES snes, Mat m, Vec q, Vec f)
solve crack propagation problem
std::string mwlsApproxFile
Name of file with internal stresses.
map< int, boost::shared_ptr< SurfaceSlidingConstrains > > surfaceConstrain
MoFEMErrorCode readMedFile()
read mesh file
boost::shared_ptr< FaceElementForcesAndSourcesCore > feRhsSurfaceForceALEMaterial
boost::shared_ptr< CrackFrontElement > feRhs
Integrate elastic FE.
MoFEMErrorCode buildBothSidesFieldId(const BitRefLevel bit_spatial, const BitRefLevel bit_material, const bool proc_only=false, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
Lagrange multipliers field which constrains material displacements.
PetscBool propagateCrack
If true crack propagation is calculated.
boost::shared_ptr< MortarContactProblem::MortarContactElement > fePostProcMortarContact
std::map< std::string, BitRefLevel > mapBitLevel
MoFEMErrorCode declareArcLengthFE(const BitRefLevel bits, const int verb=QUIET)
create arc-length element entity and declare elemets
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrains > feGriffithConstrainsLhs
boost::shared_ptr< Smoother::MyVolumeFE > feSmootherLhs
Integrate smoothing operators.
double rHo0
Reference density if bone is analyzed.
double smootherGamma
Controls mesh smoothing.
MoFEMErrorCode setCrackFrontBitLevel(BitRefLevel from_bit, BitRefLevel bit, const int nb_levels=2, const bool debug=false)
Set bit ref level for entities adjacent to crack front.
MoFEMErrorCode testJacobians(const BitRefLevel bit, const BitRefLevel mask, tangent_tests test)
test LHS Jacobians
MoFEMErrorCode declareEdgeFE(std::string fe_name, const BitRefLevel bit, const BitRefLevel mask, const bool proc_only=true, const int verb=QUIET, const bool debug=false)
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feRhsSimpleContact
MoFEMErrorCode createSurfaceProjectionDM(SmartPetscObj< DM > &dm, SmartPetscObj< DM > dm_material, const std::string prb_name, const std::vector< int > surface_ids, const std::vector< std::string > fe_list, const int verb=QUIET)
create DM to calculate projection matrices (sub DM of DM material)
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > analiticalSurfaceElement
boost::shared_ptr< ConstrainMatrixCtx > projFrontCtx
Data structure to project on crack front.
PetscBool onlyHookeFromOptions
True if only Hooke material is applied.
MoFEMErrorCode getInterfaceVersion(Version &version) const
MoFEMErrorCode buildProblemFiniteElements(BitRefLevel bit1, BitRefLevel bit2, std::vector< int > surface_ids, const int verb=QUIET, const bool debug=false)
Build problem finite elements.
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrainsDelta > feGriffithConstrainsDelta
boost::shared_ptr< NonlinearElasticElement > elasticFe
boost::shared_ptr< MWLSApprox > & getMWLSApprox()
MoFEMErrorCode calculateMaterialForcesDM(DM dm, Vec q, Vec f, const int verb=QUIET, const bool debug=false)
assemble material forces, by running material finite element instance
MoFEMErrorCode declareFrontFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
MoFEMErrorCode buildProblemFields(const BitRefLevel &bit1, const BitRefLevel &mask1, const BitRefLevel &bit2, const int verb=QUIET, const bool debug=false)
Build problem fields.
MoFEMErrorCode postProcessDM(DM dm, const int step, const std::string fe_name, const bool approx_internal_stress)
post-process results for elastic solution
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feLhsSimpleContactALE
boost::shared_ptr< SimpleContactProblem::CommonDataSimpleContact > commonDataSimpleContact
double betaGc
heterogeneous Griffith energy exponent
EntityHandle arcLengthVertex
Vertex associated with dof for arc-length control.
double griffithR
Griffith regularisation parameter.
std::string mwlsEigenStressTagName
Name of tag with eigen stresses.
MoFEMErrorCode partitionMesh(BitRefLevel bit1, BitRefLevel bit2, int verb=QUIET, const bool debug=false)
partotion mesh
PetscBool isGcHeterogeneous
flag for heterogeneous gc
MoFEMErrorCode addMWLSStressOperators(boost::shared_ptr< CrackFrontElement > &fe_rhs, boost::shared_ptr< CrackFrontElement > &fe_lhs)
boost::shared_ptr< Smoother::MyVolumeFE > feSmootherRhs
Integrate smoothing operators.
std::string mwlsStressTagName
Name of tag with internal stresses.
map< EntityHandle, double > mapJ
hashmap of J - release energy at nodes
boost::shared_ptr< DirichletFixFieldAtEntitiesBc > fixMaterialEnts
boost::shared_ptr< FEMethod > zeroFlambda
assemble F_lambda vector
boost::shared_ptr< FaceElementForcesAndSourcesCore > feLhsSurfaceForceALEMaterial
boost::shared_ptr< NonlinearElasticElement > materialFe
boost::shared_ptr< boost::ptr_map< string, EdgeForce > > edgeForces
assemble edge forces
MoFEMErrorCode buildSurfaceFields(const BitRefLevel bit, const bool proc_only=true, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
build fields with Lagrange multipliers to constrain surfaces
MoFEMErrorCode setSingularDofs(const string field_name, const int verb=QUIET)
set singular dofs (on edges adjacent to crack front) from geometry
MoFEMErrorCode savePositionsOnCrackFrontDM(DM dm, Vec q, const int verb=QUIET, const bool debug=false)
boost::shared_ptr< NeumannForcesSurface::MyTriangleFE > feMaterialAnaliticalTraction
Surface elment to calculate tractions in material space.
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrains > feGriffithConstrainsRhs
MoFEMErrorCode getOptions()
Get options form command line.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Getting interface of core database.
MoFEMErrorCode assembleMaterialForcesDM(DM dm, const int verb=QUIET, const bool debug=false)
create material element instance
map< EntityHandle, double > mapSmoothingForceFactor
MoFEMErrorCode calculateGriffithForce(DM dm, const double gc, Vec f_griffith, const int verb=QUIET, const bool debug=false)
calculate Griffith (driving) force
boost::shared_ptr< FEMethod > getFrontArcLengthControl(boost::shared_ptr< ArcLengthCtx > arc_ctx)
boost::shared_ptr< MortarContactProblem > mortarContactProblemPtr
boost::shared_ptr< FaceElementForcesAndSourcesCore > feLhsSpringALEMaterial
boost::shared_ptr< Smoother > smootherFe
MoFEMErrorCode setSpatialPositionFromCoords()
set spatial field from nodes
boost::shared_ptr< MortarContactProblem::CommonDataMortarContact > commonDataMortarContact
boost::shared_ptr< FaceElementForcesAndSourcesCore > feRhsSpringALEMaterial
MoFEMErrorCode updateMaterialFixedNode(const bool fix_front, const bool fix_small_g, const bool debug=false)
Update fixed nodes.
const boost::scoped_ptr< UnknownInterface > cpSolversPtr
MoFEMErrorCode createBcDM(SmartPetscObj< DM > &dm, const std::string prb_name, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set())
Create problem to calculate boundary conditions.
double smootherAlpha
Controls mesh smoothing.
MoFEMErrorCode createEigenElasticDM(SmartPetscObj< DM > &dm, const std::string prb_name, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set())
Create elastic problem DM.
boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > skinOrientation
MoFEMErrorCode buildElasticFields(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
Declate fields for elastic analysis.
MoFEMErrorCode setSingularElementMatrialPositions(const int verb=QUIET)
set singular dofs of material field (on edges adjacent to crack front) from geometry
boost::shared_ptr< DofEntity > arcLengthDof
MoFEMErrorCode declarePressureAleFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
Declare FE for pressure BC in ALE formulation (in material domain)
MoFEMErrorCode setMaterialPositionFromCoords()
set material field from nodes
MoFEMErrorCode addElasticFEInstancesToSnes(DM dm, Mat m, Vec q, Vec f, boost::shared_ptr< FEMethod > arc_method=boost::shared_ptr< FEMethod >(), boost::shared_ptr< ArcLengthCtx > arc_ctx=nullptr, const int verb=QUIET, const bool debug=false)
map< EntityHandle, VectorDouble3 > mapMatForce
hashmap of material force at nodes
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfacePressureAle
assemble surface pressure (ALE)
MoFEMErrorCode calculateElasticEnergy(DM dm, const std::string msg="")
assemble elastic part of matrix, by running elastic finite element instance
MoFEMErrorCode createMaterialForcesDM(SmartPetscObj< DM > &dm, SmartPetscObj< DM > dm_material, const std::string prb_name, const int verb=QUIET)
Create DM for calculation of material forces (sub DM of DM material)
boost::shared_ptr< VolumeLengthQuality< adouble > > volumeLengthAdouble
MoFEMErrorCode projectGriffithForce(DM dm, Vec f_griffith, Vec f_griffith_proj, const int verb=QUIET, const bool debug=false)
project Griffith forces
double griffithE
Griffith stability parameter.
MoFEMErrorCode declareSurfaceFE(std::string fe_name, const BitRefLevel bit, const BitRefLevel mask, const std::vector< int > &ids, const bool proc_only=true, const int verb=QUIET, const bool debug=false)
declare surface sliding elements
boost::shared_ptr< ArcLengthCtx > arcEigenCtx
Range constrainedInterface
Range of faces on the constrained interface.
MoFEMErrorCode resolveSharedBitRefLevel(const BitRefLevel bit, const int verb=QUIET, const bool debug=false)
resole shared entities by bit level
MoFEMErrorCode assembleSmootherForcesDM(DM dm, const std::vector< int > ids, const int verb=QUIET, const bool debug=false)
create smoothing element instance
boost::shared_ptr< boost::ptr_map< string, NodalForce > > nodalForces
assemble nodal forces
boost::shared_ptr< VolumeLengthQuality< double > > volumeLengthDouble
MoFEMErrorCode setCoordsFromField(const std::string field_name="MESH_NODE_POSITIONS")
set coordinates from field
MoFEMErrorCode projectMaterialForcesDM(DM dm_project, Vec f, Vec f_proj, const int verb=QUIET, const bool debug=false)
project material forces along the crack elongation direction
MoFEMErrorCode buildEdgeFields(const BitRefLevel bit, const bool proc_only=true, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
build fields with Lagrange multipliers to constrain edges
map< EntityHandle, double > mapG3
hashmap of g3 - release energy at nodes
PetscBool isPressureAle
If true surface pressure is considered in ALE.
MoFEMErrorCode createElasticDM(SmartPetscObj< DM > &dm, const std::string prb_name, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set())
Create elastic problem DM.
map< EntityHandle, double > mapG1
hashmap of g1 - release energy at nodes
MoFEMErrorCode declareCrackSurfaceFE(std::string fe_name, const BitRefLevel bit, const BitRefLevel mask, const bool proc_only=true, const int verb=QUIET, const bool debug=false)
MoFEMErrorCode createCrackPropagationDM(SmartPetscObj< DM > &dm, const std::string prb_name, SmartPetscObj< DM > dm_elastic, SmartPetscObj< DM > dm_material, const BitRefLevel bit, const BitRefLevel mask, const std::vector< std::string > fe_list)
Create DM by composition of elastic DM and material DM.
boost::shared_ptr< MortarContactProblem::MortarContactElement > feLhsMortarContact
MoFEMErrorCode tetsingReleaseEnergyCalculation()
This is run with ctest.
MoFEMErrorCode saveEachPart(const std::string prefix, const Range &ents)
Save entities on ech processor.
MoFEMErrorCode calculateSmoothingForceFactor(const int verb=QUIET, const bool debug=true)
boost::shared_ptr< ConstrainMatrixCtx > projSurfaceCtx
Data structure to project on the body surface.
static FTensor::Tensor2_symmetric< double, 3 > analyticalStrainFunction(FTensor::Tensor1< FTensor::PackPtr< double *, 1 >, 3 > &t_coords)
boost::shared_ptr< ArcLengthCtx > arcCtx
MoFEMErrorCode cleanSingularElementMatrialPositions()
set maetrial field order to one
boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > crackOrientation
map< int, boost::shared_ptr< EdgeSlidingConstrains > > edgeConstrains
PetscBool solveEigenStressProblem
Solve eigen problem.
MoFEMErrorCode writeCrackFont(const BitRefLevel &bit, const int step)
MoFEMErrorCode setFieldFromCoords(const std::string field_name)
set field from node positions
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > commonDataSpringsALE
MoFEMErrorCode addMWLSDensityOperators(boost::shared_ptr< CrackFrontElement > &fe_rhs, boost::shared_ptr< CrackFrontElement > &fe_lhs)
MoFEMErrorCode createCrackFrontAreaDM(SmartPetscObj< DM > &dm, SmartPetscObj< DM > dm_material, const std::string prb_name, const bool verb=QUIET)
create DM to calculate Griffith energy
MoFEMErrorCode declareElasticFE(const BitRefLevel bit1, const BitRefLevel mask1, const BitRefLevel bit2, const BitRefLevel mask2, const bool add_forces=true, const bool proc_only=true, const int verb=QUIET)
declare elastic finite elements
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feLhsSimpleContact
boost::shared_ptr< CrackFrontElement > feCouplingMaterialLhs
FE instance to assemble coupling terms.
double nBone
Exponent parameter in bone density.
boost::shared_ptr< CrackFrontElement > feCouplingElasticLhs
FE instance to assemble coupling terms.
boost::shared_ptr< NeumannForcesSurface::DataAtIntegrationPts > commonDataSurfaceForceAle
common data at integration points (ALE)
MoFEMErrorCode assembleElasticDM(const std::string mwls_stress_tag_name, const int verb=QUIET, const bool debug=false)
create elastic finite element instance for spatial assembly
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfaceForces
assemble surface forces
boost::shared_ptr< GriffithForceElement::MyTriangleFE > feGriffithForceLhs
const boost::scoped_ptr< UnknownInterface > cpMeshCutPtr
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feLhsSimpleContactALEMaterial
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfaceForceAle
assemble surface pressure (ALE)
boost::shared_ptr< GriffithForceElement > griffithForceElement
MoFEMErrorCode assembleCouplingForcesDM(DM dm, const int verb=QUIET, const bool debug=false)
assemble coupling element instances
MoFEMErrorCode declareSurfaceForceAleFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
Declare FE for pressure BC in ALE formulation (in material domain)
boost::shared_ptr< MWLSApprox > mwlsApprox
MoFEMErrorCode addMaterialFEInstancesToSnes(DM dm, const bool fix_crack_front, const int verb=QUIET, const bool debug=false)
add material elements instances to SNES
MoFEMErrorCode calculateReleaseEnergy(DM dm, Vec f_material_proj, Vec f_griffith_proj, Vec f_lambda, const double gc, const int verb=QUIET, const bool debug=true)
calculate release energy
bool onlyHooke
True if only Hooke material is applied.
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > fePostProcSimpleContact
MoFEMErrorCode addPropagationFEInstancesToSnes(DM dm, boost::shared_ptr< FEMethod > arc_method, boost::shared_ptr< ArcLengthCtx > arc_ctx, const std::vector< int > &surface_ids, const int verb=QUIET, const bool debug=false)
add finite element to SNES for crack propagation problem
boost::shared_ptr< CrackFrontElement > feLhs
Integrate elastic FE.
MoFEMErrorCode calculateSurfaceProjectionMatrix(DM dm_front, DM dm_project, const std::vector< int > &ids, const int verb=QUIET, const bool debug=false)
assemble projection matrices
boost::shared_ptr< BothSurfaceConstrains > bothSidesContactConstrains
MoFEMErrorCode buildArcLengthField(const BitRefLevel bit, const bool build_fields=true, const int verb=QUIET)
Declate field for arc-length.
int postProcLevel
level of postprocessing (amount of output files)
MoFEMErrorCode resolveShared(const Range &tets, Range &proc_ents, const int verb=QUIET, const bool debug=false)
resolve shared entities
boost::shared_ptr< AnalyticalDirichletBC::DirichletBC > getAnalyticalDirichletBc()
boost::shared_ptr< FaceElementForcesAndSourcesCore > feLhsSpringALE
boost::shared_ptr< FEMethod > assembleFlambda
assemble F_lambda vector
boost::shared_ptr< BothSurfaceConstrains > closeCrackConstrains
PetscBool areSpringsAle
If true surface spring is considered in ALE.
boost::shared_ptr< ObosleteUsersModules::TangentWithMeshSmoothingFrontConstrain > tangentConstrains
Constrains crack front in tangent directiona.
MoFEMErrorCode calculateSmoothingForcesDM(DM dm, Vec q, Vec f, const int verb=QUIET, const bool debug=false)
assemble smoothing forces, by running material finite element instance
boost::shared_ptr< WrapMPIComm > moabCommWorld
MoFEMErrorCode createProblemDataStructures(const std::vector< int > surface_ids, const int verb=QUIET, const bool debug=false)
Construct problem data structures.
boost::shared_ptr< ContactSearchKdTree::ContactCommonData_multiIndex > contactSearchMultiIndexPtr
boost::shared_ptr< MortarContactProblem::MortarContactElement > feRhsMortarContact
map< EntityHandle, VectorDouble3 > mapGriffith
hashmap of Griffith energy at nodes
MoFEMErrorCode addSmoothingFEInstancesToSnes(DM dm, const bool fix_crack_front, const int verb=QUIET, const bool debug=false)
add softening elements instances to SNES
MoFEMErrorCode declareMaterialFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
declare material finite elements
boost::shared_ptr< ArcLengthCtx > & getEigenArcCtx()
MoFEMErrorCode createDMs(SmartPetscObj< DM > &dm_elastic, SmartPetscObj< DM > &dm_eigen_elastic, SmartPetscObj< DM > &dm_material, SmartPetscObj< DM > &dm_crack_propagation, SmartPetscObj< DM > &dm_material_forces, SmartPetscObj< DM > &dm_surface_projection, SmartPetscObj< DM > &dm_crack_srf_area, std::vector< int > surface_ids, std::vector< std::string > fe_surf_proj_list)
Crate DMs for all problems and sub problems.
boost::shared_ptr< SimpleContactProblem::CommonDataSimpleContact > commonDataSimpleContactALE
PetscBool ignoreMaterialForce
If true surface pressure is considered in ALE.
MoFEMErrorCode buildCrackSurfaceFieldId(const BitRefLevel bit, const bool proc_only=true, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
declare crack surface files
boost::shared_ptr< GriffithForceElement::MyTriangleFE > feGriffithForceRhs
MoFEMErrorCode buildCrackFrontFieldId(const BitRefLevel bit, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
declare crack surface files
boost::shared_ptr< CrackFrontElement > feEnergy
Integrate energy.
boost::shared_ptr< SimpleContactProblem > contactProblem
boost::shared_ptr< CrackFrontElement > feMaterialLhs
Integrate material stresses, assemble matrix.
int residualStressBlock
Block on which residual stress is applied.
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feRhsSimpleContactALEMaterial
MoFEMErrorCode unsetSingularElementMatrialPositions()
remove singularity
MoFEMErrorCode calculateFrontProjectionMatrix(DM dm_surface, DM dm_project, const int verb=QUIET, const bool debug=false)
assemble crack front projection matrix (that constrains crack area growth)
PetscBool isSurfaceForceAle
If true surface pressure is considered in ALE.
std::string mwlsRhoTagName
Name of tag with density.
boost::shared_ptr< NeumannForcesSurface::DataAtIntegrationPts > commonDataSurfacePressureAle
common data at integration points (ALE)
boost::shared_ptr< ArcLengthCtx > & getArcCtx()
Deprecated interface functions.
Data structure to exchange data between mofem and User Loop Methods on entities.
base class for all interface classes
Class implemented by user to detect face orientation.