v0.9.1
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 
22 namespace FractureMechanics {
23 
24 /**
25  * Three interfaces used to solve fracture problem
26  */
31 };
32 
33 /**
34  * Unknown interface ID variable shared by all structs and used to register
35  * interface for CrackPropagation at CrackPropagation constructor
36  */
37 static const MOFEMuuid IDD_MOFEMCrackPropagationInterface =
39 
40 /**
41  * @brief Tapes numbers used by ADOL-C.
42  *
43  */
44 enum AdolcTags {
56 };
57 
58 /**
59  * Crack propagation data and methods
60  * \nosubgrouping
61  *
62  * - Functions strarting with build, for example buildElastic,
63  * add and build fields.
64  *
65  * - Functions starting with declare, for example declareElasticFE,
66  * add and build finite elements.
67  *
68  * - Functions starting with assemble, for example assembleElasticDM,
69  * create finite element instances and operators.
70  *
71  * - Functions starting with create, for example createElasticDM, create
72  * problems.
73  *
74  */
76 
77  static bool debug;
78 
79  MoFEMErrorCode getInterfaceVersion(Version &version) const {
81  version = Version(FM_VERSION_MAJOR, FM_VERSION_MINOR, FM_VERSION_BUILD);
83  }
84 
85  /**
86  * \brief Getting interface of core database
87  * @param uuid unique ID of interface that can be either for CrackPropagation,
88  * CPSolvers or CPMeshCut interface
89  * @param iface returned pointer to interface
90  * @return error code
91  */
92  MoFEMErrorCode query_interface(const MOFEMuuid &uuid,
93  UnknownInterface **iface) const;
94 
95  MPI_Comm moabCommWorld;
100 
101  double gC; ///< Griffith energy
102  double griffithE; ///< Griffith stability parameter
103  double griffithR; ///< Griffith regularisation parameter
104 
105  PetscBool isPartitioned;
106  PetscBool propagateCrack; ///< If true crack propagation is calculated
109  PetscBool addSingularity;
111 
112  std::string mwlsApproxFile; ///< Name of file with internal stresses
113  std::string mwlsStressTagName; ///< Name of tag with internal stresses
114  std::string mwlsRhoTagName; ///< Name of tag with density
115  int residualStressBlock; ///< Block on which residual stress is applied
116  int densityMapBlock; ///< Block on which density is applied (Note: if used
117  ///< with internal stress, it has to be the same ID!)
118  std::string medFileName;
119  std::string configFile;
121 
122  int postProcLevel; ///< level of postprocessing (amount of output files)
124 
127 
128  /**
129  * \brief Constructor of Crack Propagation, prints fracture module version and
130  registers the three interfaces
131  * used to solve fracture problem, i.e. CrackPropagation, CPSolvers and
132  * CPMeshCut
133  * @param m_field interface to MoAB database
134  * @param approx_order order of approximation space, default is 2
135  * @param geometry_order order of geometry, default is 1 (can be up to 2)
136  */
137  CrackPropagation(MoFEM::Interface &m_field, const int approx_order = 2,
138  const int geometry_order = 1);
139 
140  virtual ~CrackPropagation();
141 
142  /**
143  * \brief Get options form command line
144  * @return error code
145  */
147 
148  /**
149  * \brief This is run with ctest
150  * @return error code
151  */
153 
154  inline int &getNbLoadSteps() { return nbLoadSteps; }
155  inline int getNbLoadSteps() const { return nbLoadSteps; }
156  inline double &getLoadScale() { return loadScale; }
157  inline double getLoadScale() const { return loadScale; }
158  inline int &getNbCutSteps() { return nbCutSteps; }
159  inline int getNbCutSteps() const { return nbCutSteps; }
160 
161  /// \brief read mesh file
163 
164  /// \brief read configuration file
166 
167  /**
168  * @brief Save entities on ech processor.
169  *
170  * Usually used for debugging, to check if meshes are consistent on all
171  * processors.
172  *
173  * @param prefix file name prefix
174  * @param ents entities to save
175  * @return MoFEMErrorCode error code
176  */
177  MoFEMErrorCode saveEachPart(const std::string prefix, const Range &ents);
178 
179  /** \name Resolve shared entities and partition mesh */
180 
181  /**@{*/
182 
183  /// \brief partotion mesh
185  int verb = QUIET, const bool debug = false);
186 
187  /// \brief resolve shared entities
188  MoFEMErrorCode resolveShared(const Range &tets, Range &proc_ents,
189  const int verb = QUIET,
190  const bool debug = false);
191 
192  /// \brief resole shared entities by bit level
194  const int verb = QUIET,
195  const bool debug = false);
196 
197  /**@}*/
198 
199  /** \name Set bit level for material problem */
200 
201  /**@{*/
202 
203  /**
204  * @brief Set bit ref level for entities adjacent to crack front
205  *
206  * @param from_bit
207  * @param bit
208  * @param nb_levels adjacency bridge levels
209  * @return MoFEMErrorCode
210  */
212  const int nb_levels = 2,
213  const bool debug = false);
214 
215  /**@}*/
216 
217  /** \name Build problem */
218 
219  /**@{*/
220 
221  /**
222  * \brief Build problem fields
223  * @param bit1 bit level for computational domain
224  * @param mask1 bit mask
225  * @param bit2 bit level for elements near crack front
226  * @param mask2 bit mask 2
227  * @param verb verbosity level
228  * @param debug
229  * @return error code
230  *
231  */
233  const BitRefLevel &mask1,
234  const BitRefLevel &bit2,
235  const int verb = QUIET,
236  const bool debug = false);
237 
238  /**
239  * \brief Build problem finite elements
240  * @param bit1 bit level for computational domain
241  * @param bit2 bit level for elements near crack front
242  * @param surface_ids meshsets Ids on body skin
243  * @param verb verbosity level
244  * @param debug
245  * @return error code
246  *
247  */
249  std::vector<int> surface_ids,
250  const int verb = QUIET,
251  const bool debug = false);
252 
253  /**
254  * \brief Construct problem data structures
255  * @param bit0 Auxiliary bit ref level
256  * @param bit1 bit level for computational domain
257  * @param bit2 bit level for elements near crack front
258  * @param surface_ids surface_ids meshsets Ids on body skin
259  * @param verb verbosity level
260  * @param debug
261  * @return error code
262  */
264  const BitRefLevel bit1,
265  const BitRefLevel bit2,
266  const std::vector<int> surface_ids,
267  const int verb = QUIET,
268  const bool debug = false);
269 
270  /**
271  * \brief Crate DMs for all problems and sub problems
272  * @param bit1 domain bit ref level
273  * @param bit2 material domain
274  * @param dm_elastic elastic DM
275  * @param dm_material material DM
276  * @param dm_crack_propagation composition of elastic DM and material DM
277  * @param dm_material_forces used to calculate material forces, sub-dm of
278  * dm_crack_propagation
279  * @param dm_surface_projection dm to project material forces on surfaces,
280  * sub-dm of dm_crack_propagation
281  * @param dm_crack_srf_area dm to project material forces on crack
282  * surfaces, sub-dm of dm_crack_propagation
283  * @param surface_ids IDs of surface meshsets
284  * @param fe_surf_proj_list name of elements on surface elements
285  * @return error code
286  */
287  MoFEMErrorCode createDMs(BitRefLevel bit1, BitRefLevel bit2, DM *dm_elastic,
288  DM *dm_material, DM *dm_crack_propagation,
289  DM *dm_material_forces, DM *dm_surface_projection,
290  DM *dm_crack_srf_area, std::vector<int> surface_ids,
291  std::vector<std::string> fe_surf_proj_list);
292 
293  /// \brief destroy DMs
294  MoFEMErrorCode destroyDMs(DM *dm_elastic, DM *dm_material,
295  DM *dm_crack_propagation, DM *dm_material_forces,
296  DM *dm_surface_projection, DM *dm_crack_srf_area);
297 
298  MoFEMErrorCode deleteEntities(const int verb = QUIET,
299  const bool debug = false);
300 
301  /**@}*/
302 
303  /** \name Build fields */
304 
305  /**@{*/
306 
307  /**
308  * \brief Declate fields for elastic analysis
309  * @param bit bit refinement level
310  * @return error code
311  */
313  buildElasticFields(const BitRefLevel bit = BitRefLevel().set(0),
314  const BitRefLevel mask = BitRefLevel().set(),
315  const bool proc_only = true,
316  const bool build_fields = true,const int verb = QUIET);
317 
318  /**
319  * \brief Declate field for arc-length
320  * @param bit bit refinement level
321  * @return error code
322  */
324  const bool build_fields = true,
325  const int verb = QUIET);
326 
327  /// \brief build fields with Lagrange multipliers to constrain surfaces
329  buildSurfaceFields(const BitRefLevel bit = BitRefLevel().set(0),
330  const bool proc_only = true,
331  const bool build_fields = true, const int verb = QUIET,
332  const bool debug = false);
333 
334  /// \brief build fields with Lagrange multipliers to constrain edges
336  const bool proc_only = true,
337  const bool build_fields = true,
338  const int verb = QUIET,
339  const bool debug = false);
340 
341  /// \brief declare crack surface files
343  const bool proc_only = true,
344  const bool build_fields = true,
345  const int verb = QUIET,
346  const bool debug = false);
347 
348  /// \brief declare crack surface files
350  buildCrackFrontFieldId(const BitRefLevel bit = BitRefLevel().set(0),
351  const bool build_fields = true, const int verb = QUIET,
352  const bool debug = false);
353 
354  /**
355  * \brief Lagrange multipliers field which constrains material displacements
356  *
357  * Nodes on both sides are constrainded to have the same material
358  * displacements
359  *
360  */
362  buildBothSidesFieldId(const BitRefLevel bit = BitRefLevel().set(0),
363  const bool proc_only = false,
364  const bool build_fields = true,
365  const int verb = QUIET, const bool debug = false);
366 
367  /** \brief Zero fields with lagrange multipliers
368  */
370 
371  /**@}*/
372 
373  /** \name Build finite elements */
374 
375  /**@{*/
376 
377  /// \brief declare elastic finite elements
379  declareElasticFE(const BitRefLevel bit1, const BitRefLevel mask1,
380  const BitRefLevel bit2, const BitRefLevel mask2,
381  const bool add_forces = true, const bool proc_only = true,
382  const int verb = QUIET);
383 
384  /// \brief create arc-length element entity and declare elemets
386  const int verb = QUIET);
387 
389  declareExternalForcesFE(const BitRefLevel bit = BitRefLevel().set(0),
390  const BitRefLevel mask = BitRefLevel().set(),
391  const bool proc_only = true);
392 
393  /// \brief declare material finite elements
395  const BitRefLevel mask = BitRefLevel().set(),
396  const bool proc_only = true,
397  const bool verb = QUIET);
398 
400  const BitRefLevel mask = BitRefLevel().set(),
401  const bool proc_only = true,
402  const bool verb = QUIET);
403 
405  declareBothSidesFE(const BitRefLevel bit = BitRefLevel().set(0),
406  const BitRefLevel mask = BitRefLevel().set(),
407  const bool proc_only = true, const bool verb = QUIET);
408 
409  /// \brief declare mesh smoothing finite elements
411  declareSmoothingFE(const BitRefLevel bit = BitRefLevel().set(0),
412  const BitRefLevel mask = BitRefLevel().set(),
413  const bool proc_only = true, const bool verb = QUIET);
414 
415  /// \brief declare surface sliding elements
416  MoFEMErrorCode declareSurfaceFE(std::string fe_name, const BitRefLevel bit,
417  const BitRefLevel mask,
418  const std::vector<int> &ids,
419  const bool proc_only = true,
420  const int verb = QUIET,
421  const bool debug = false);
422 
423  // \brief declare edge sliding elements
424  MoFEMErrorCode declareEdgeFE(std::string fe_name, const BitRefLevel bit,
425  const BitRefLevel mask,
426  const bool proc_only = true,
427  const int verb = QUIET,
428  const bool debug = false);
429 
431  declareCrackSurfaceFE(std::string fe_name, const BitRefLevel bit,
432  const BitRefLevel mask, const bool proc_only = true,
433  const int verb = QUIET, const bool debug = false);
434 
435  /**@}*/
436 
437  /** \name Create DMs & problems */
438 
439  /**@{*/
440 
441  /** \brief Create DM by composition of elastic DM and material DM
442  * @param dm composition of elastic DM and material DM
443  * @param prb_name problem name
444  * @param dm_elastic elastic DM
445  * @param dm_material material DM
446  * @param bit domain bit ref level
447  * @param mask mask ref level for problem
448  * @return error code
449  */
451  createCrackPropagationDM(DM *dm, const std::string prb_name, DM dm_elastic,
452  DM dm_material, const BitRefLevel bit,
453  const BitRefLevel mask,
454  const std::vector<std::string> fe_list);
455 
456  /** \brief Create elastic problem DM
457  * @param dm elastic DM
458  * @param prb_name problem name
459  * @param bit domain bit ref level
460  * @param mask mask ref level for problem
461  * @return error code
462  */
463  MoFEMErrorCode createElasticDM(DM *dm, const std::string prb_name,
464  const BitRefLevel bit = BitRefLevel().set(0),
465  const BitRefLevel mask = BitRefLevel().set());
466 
467  /** \brief Create problem to calculate boundary conditions
468  * @param dm newly created boundary condition DM
469  * @param prb_name problem name
470  * @param bit auxiliary bit ref level
471  * @param mask mask ref level for problem
472  * @return error code
473  * \note BcDM is only useed for testing, to enforce boundary conditions
474 staisifing anlylitical solution to which numerical solution is compared.
475  */
476  MoFEMErrorCode createBcDM(DM *dm, const std::string prb_name,
477  const BitRefLevel bit = BitRefLevel().set(0),
478  const BitRefLevel mask = BitRefLevel().set());
479 
480 
481  /** \brief Create DM fto calculate material problem
482  * @param dm used to calculate material forces, sub-dm of
483  * dm_crack_propagation
484  * @param prb_name problem name
485  * @param bit material domain
486  * @param mask dof mask ref level for problem
487  * @param fe_list name of elements on surface elements
488  * @param debug flag for debugging
489  * @return error code
490  */
491  MoFEMErrorCode createMaterialDM(DM *dm, const std::string prb_name,
492  const BitRefLevel bit, const BitRefLevel mask,
493  const std::vector<std::string> fe_list,
494  const bool debug = false);
495 
496  /** \brief Create DM for calculation of material forces (sub DM of DM
497  * material)
498  * @param dm used to calculate material forces, sub-dm of
499  * dm_crack_propagation
500  * @param prb_name problem name
501  * @param verb compilation parameter determining the amount of
502  * information printed
503  * @return error code
504  */
505  MoFEMErrorCode createMaterialForcesDM(DM *dm, DM dm_material,
506  const std::string prb_name,
507  const int verb = QUIET);
508 
509  /** \brief create DM to calculate projection matrices (sub DM of DM
510  * material)
511  * @param dm DM to project material forces on surfaces,
512  * sub-dm of dm_crack_propagation
513  * @param dm_material material DM
514  * @param prb_name problem name
515  * @param surface_ids IDs of surface meshsets
516  * @param fe_list name of elements on surface elements
517  * @param verb compilation parameter determining the amount
518  * of information printed
519  * @return error code
520  */
522  createSurfaceProjectionDM(DM *dm, DM dm_material, const std::string prb_name,
523  const std::vector<int> surface_ids,
524  const std::vector<std::string> fe_list,
525  const int verb = QUIET);
526 
527  /** \brief create DM to calculate Griffith energy
528  * @param dm dm to project material forces on crack
529  * surfaces, sub-dm of dm_crack_propagation
530  * @param dm_material material DM
531  * @param prb_name problem name
532  * @param verb compilation parameter determining the amount
533  * of information printed
534  * @return error code
535  */
536  MoFEMErrorCode createCrackFrontAreaDM(DM *dm, DM dm_material,
537  const std::string prb_name,
538  const bool verb = QUIET);
539 
540  /**@}*/
541 
542  /** \name Create finite element instances */
543 
544  /**@{*/
545 
546  /** \brief create elastic finite element instance for spatial assembly
547  * @param dm elastic DM
548  * @param m stiffness sub-matrix for elastic DM
549  * @param q DoF sub-vector for elastic DM
550  * @param f Rhs sub-vector for material DM
551  * @param verb compilation parameter determining the amount
552  * of information printed (not in use here)
553  * @param debug flag for debugging
554  * @return error code
555  */
556  MoFEMErrorCode assembleElasticDM(DM dm, Mat m, Vec q, Vec f,
557  const int verb = QUIET,
558  const bool debug = false);
559 
560  // \brief add elastic element instances to SNES
562  DM dm, Mat m, Vec q, Vec f,
563  boost::shared_ptr<FEMethod> arc_method = boost::shared_ptr<FEMethod>(),
564  const int verb = QUIET, const bool debug = false);
565 
566  /// \brief create material element instance
567  MoFEMErrorCode assembleMaterialForcesDM(DM dm, const int verb = QUIET,
568  const bool debug = false);
569 
570  /// \brief create smoothing element instance
571  MoFEMErrorCode assembleSmootherForcesDM(DM dm, const std::vector<int> ids,
572  const int verb = QUIET,
573  const bool debug = false);
574 
575  /// \brief assemble coupling element instances
576  MoFEMErrorCode assembleCouplingForcesDM(DM dm, const int verb = QUIET,
577  const bool debug = false);
578 
579  /// \brief add material elements instances to SNES
580  MoFEMErrorCode addMaterialFEInstancesToSnes(DM dm, const bool fix_crack_front,
581  const int verb = QUIET,
582  const bool debug = false);
583 
584  /// \brief add softening elements instances to SNES
585  MoFEMErrorCode addSmoothingFEInstancesToSnes(DM dm, const bool fix_crack_front,
586  const int verb = QUIET,
587  const bool debug = false);
588 
589  /// \brief add finite element to SNES for crack propagation problem
591  addPropagationFEInstancesToSnes(DM dm, boost::shared_ptr<FEMethod> &arc_method,
592  const std::vector<int> &surface_ids,
593  const int verb = QUIET,
594  const bool debug = false);
595 
596  /// \brief test MWLS
598  const BitRefLevel mask);
599  /// \brief test MWLS
601  const BitRefLevel mask);
602 
603  /// \brief Update fixed nodes
604  MoFEMErrorCode updateMaterialFixedNode(const bool fix_front,
605  const bool fix_small_g,
606  const bool debug = false);
607 
608  /**@}*/
609 
610  /** \name Assemble and calculate */
611 
612  /**@{*/
613 
614  /// \brief assemble elastic part of matrix, by running elastic finite element
615  /// instance
616  MoFEMErrorCode calculateElasticEnergy(DM dm, const std::string msg = "");
617 
618  /** \brief assemble material forces, by running material finite element
619  * instance
620  * @param dm dm used to calculate material forces, sub-dm
621  * of dm_crack_propagation
622  * @param q DoF sub-vector for crack propagation DM
623  * @param f Rhs sub-vector for crack propagation DM
624  * @param verb compilation parameter determining the amount
625  * of information printed
626  * @param debug flag for debugging
627  * @return error code
628  */
629  MoFEMErrorCode calculateMaterialForcesDM(DM dm, Vec q, Vec f,
630  const int verb = QUIET,
631  const bool debug = false);
632 
633  /** \brief assemble smoothing forces, by running material finite element
634  * instance
635  * @param dm crack propagation dm that is composition of
636  * elastic DM and material DM
637  * @param q DoF sub-vector for crack propagation DM
638  * @param f Rhs sub-vector for crack propagation DM
639  * @param verb compilation parameter determining the amount
640  * of information printed
641  * @param debug flag for debugging
642  * @return error code
643  */
644  MoFEMErrorCode calculateSmoothingForcesDM(DM dm, Vec q, Vec f,
645  const int verb = QUIET,
646  const bool debug = false);
647 
649  const bool debug = true);
650 
651  /// \brief assemble projection matrices
652  MoFEMErrorCode calculateSurfaceProjectionMatrix(DM dm_front, DM dm_project,
653  const std::vector<int> &ids,
654  const int verb = QUIET,
655  const bool debug = false);
656 
657  /// \brief assemble crack front projection matrix (that constrains crack area
658  /// growth)
659  MoFEMErrorCode calculateFrontProjectionMatrix(DM dm_surface, DM dm_project,
660  const int verb = QUIET,
661  const bool debug = false);
662 
663  /** \brief project material forces along the crack elongation direction
664  * @param dm_project material forces DM, sub-dm of
665  * dm_crack_propagation
666  * @param f_proj Rhs sub-vector of projected material forces
667  * for material forces DM
668  * @param f Rhs sub-vector for material forces DM
669  * @param verb compilation parameter determining the amount
670  * of information printed
671  * @param debug flag for debugging
672  * @return error code
673  */
674  MoFEMErrorCode projectMaterialForcesDM(DM dm_project, Vec f, Vec f_proj,
675  const int verb = QUIET,
676  const bool debug = false);
677 
678  /// \brief calculate Griffith (driving) force
679  MoFEMErrorCode calculateGriffithForce(DM dm, const double gc, Vec f_griffith,
680  const int verb = QUIET,
681  const bool debug = false);
682 
683  /// \brief project Griffith forces
684  MoFEMErrorCode projectGriffithForce(DM dm, Vec f_griffith, Vec f_griffith_proj,
685  const int verb = QUIET,
686  const bool debug = false);
687 
688  /// \brief calculate release energy
689  MoFEMErrorCode calculateReleaseEnergy(DM dm, Vec f_material_proj,
690  Vec f_griffith_proj, Vec f_lambda,
691  const double gc, const int verb = QUIET,
692  const bool debug = true);
693 
694  /** \name Solvers */
695 
696  /**@{*/
697 
698  /** \brief solve elastic problem
699  * @param dm elastic DM
700  * @param snes snes solver for elastic DM
701  * @param m stiffness sub-matrix for elastic DM
702  * @param q DoF sub-vector for elastic DM
703  * @param f Rhs sub-vector for material DM
704  * @param snes_set_up boolean flag wheather to determin snes solver
705  * @param shell_m stiffness sub-matrix for elastic DM
706  * @return error code
707  */
708  MoFEMErrorCode solveElasticDM(DM dm, SNES snes, Mat m, Vec q, Vec f,
709  bool snes_set_up, Mat *shell_m);
710 
711  /** \brief solve crack propagation problem
712  * @param dm crack propagation DM that is composed of
713  * elastic DM and material DM
714  * @param dm_elastic elastic DM
715  * @param snes snes solver for crack propagation DM
716  * @param m stiffness sub-matrix for crack propagation DM
717  * @param q DoF sub-vector for crack propagation DM
718  * @param f Rhs sub-vector for crack propagation DM
719  * @param snes_set_up boolean flag wheather to determin snes solver
720  * @return error code
721  */
722  MoFEMErrorCode solvePropagationDM(DM dm, DM dm_elastic, SNES snes, Mat m,
723  Vec q, Vec f, bool snes_set_up, Mat *shell_m,
724  Mat *m_elastic);
725 
726  /**@}*/
727 
728  /**@}*/
729 
730  /** \name Postprocess results */
731 
732  /**@{*/
733 
735  const int verb = QUIET,
736  const bool debug = false);
737 
738  MoFEMErrorCode writeCrackFont(const BitRefLevel &bit, const int step);
739 
740  /// \brief post-process results for elastic solution
741  MoFEMErrorCode postProcessDM(DM dm, const int step, const std::string fe_name,
742  const bool approx_internal_stress);
743 
744  /**@}*/
745 
746  /** \name Set data from mesh to vectors and back */
747 
748  /**@{*/
749 
750  /// \brief set field from node positions
751  MoFEMErrorCode setFieldFromCoords(const std::string field_name);
752 
753  /// \brief set material field from nodes
755 
756  /// \brief set spatial field from nodes
758 
759  /// \brief set coordinates from field
761  setCoordsFromField(const std::string field_name = "MESH_NODE_POSITIONS");
762 
763  /// \brief set singular dofs (on edges adjacent to crack front) from geometry
764  MoFEMErrorCode setSingularDofs(const string field_name,
765  const int verb = QUIET);
766 
767  /// \brief set singular dofs of material field (on edges adjacent to crack
768  /// front) from geometry
770 
771  /// \brief remove singularity
773 
774  /// \brief set maetrial field order to one
776 
777  /**@}*/
778 
779  /** \name Auxiliary method */
780 
781  /**@{*/
782 
783  boost::shared_ptr<DofEntity> arcLengthDof;
784 
785  /// \brief set pointer to arc-length DOF
787 
788  inline boost::shared_ptr<ArcLengthCtx> &getArcCtx() { return arcCtx; }
789 
790  inline boost::shared_ptr<FEMethod> getFrontArcLengthControl() {
791  boost::shared_ptr<FEMethod> arc_method(
793  ARC_LENGTH_TAG, griffithForceElement->blockData[0],
794  griffithForceElement->commonData, "LAMBDA_CRACKFRONT_AREA", arcCtx));
795  return arc_method;
796  }
797 
798  inline boost::shared_ptr<AnalyticalDirichletBC::DirichletBC>
800  return analyticalDirichletBc;
801  }
802 
803  boost::shared_ptr<MWLSApprox> &getMWLSApprox() { return mwlsApprox; }
804 
805  /**@}*/
806 
807  /**
808  * \brief operator to post-process results on crack front nodes
809  */
811 
813  Vec F;
814  std::string tagName;
815  PostProcVertexMethod(MoFEM::Interface &m_field, Vec f, std::string tag_name)
816  : mField(m_field), F(f), tagName(tag_name) {}
817 
818  Tag tH;
819  double *aRray;
821 
823 
825  };
826 
827  /**
828  * \brief Constrains material displacement on both sides
829  *
830  * Both sides of crack surface have the same material displacements.
831  *
832  */
833  struct BothSurfaceConstrains : public FEMethod {
834 
836  double aLpha;
838  : mField(m_field), aLpha(1) {
839  ierr = getOptions();
840  CHKERRABORT(PETSC_COMM_WORLD, ierr);
841  }
842 
845  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
846  "Get Booth surface constrains element scaling",
847  "none");
848  CHKERRQ(ierr);
849  CHKERR PetscOptionsScalar("-booth_side_alpha", "scaling parameter", "",
850  aLpha, &aLpha, PETSC_NULL);
851  ierr = PetscOptionsEnd();
852  CHKERRQ(ierr);
854  }
855 
857 
861 
862  };
863 
864  /// \brief Determine face orientation
868  FaceOrientation(bool crack_front)
869  : useProjectionFromCrackFront(crack_front) {}
870  virtual ~FaceOrientation() {}
872  const EntityHandle face,
873  const BitRefLevel bit);
875  const FEMethod *fe_method_ptr) {
877  EntityHandle face = fe_method_ptr->numeredEntFiniteElementPtr->getEnt();
878  BitRefLevel bit = fe_method_ptr->problemPtr->getBitRefLevel();
879  CHKERR getElementOrientation(m_field,face,bit);
881  }
882  };
883 
884  struct ArcLengthSnesCtx : public SnesCtx {
885  boost::shared_ptr<ArcLengthCtx> arcPtr;
888  SmartPetscObj<Vec> diagM;
889  ArcLengthSnesCtx(MoFEM::Interface &m_field, const std::string &problem_name,
890  boost::shared_ptr<ArcLengthCtx> &arc_ptr,
891  DMMGViaApproxOrdersCtx *dm_ctx,
892  DMMGViaApproxOrdersCtx *dm_ctx_elastic = NULL)
893  : SnesCtx(m_field, problem_name), arcPtr(arc_ptr), dmCtx(dm_ctx),
894  dmCtxElastic(dm_ctx_elastic) {}
895  };
896 
897  // private:
898 
899  /** \name Pointers to finite element instances */
900 
901  /**@{*/
902 
903  bool onlyHooke; ///< True if only Hooke material is applied
904  PetscBool onlyHookeFromOptions; ///< True if only Hooke material is applied
905 
906  /// Integrate spatial stresses, matrix and vector (not only element but all
907  /// data structure)
908  boost::shared_ptr<NonlinearElasticElement> elasticFe;
909  boost::shared_ptr<NonlinearElasticElement> materialFe;
910  boost::shared_ptr<CrackFrontElement> feLhs; ///< Integrate elastic FE
911  boost::shared_ptr<CrackFrontElement> feRhs; ///< Integrate elastic FE
912  boost::shared_ptr<CrackFrontElement>
913  feMaterialRhs; ///< Integrate material stresses, assemble vector
914  boost::shared_ptr<CrackFrontElement>
915  feMaterialLhs; ///< Integrate material stresses, assemble matrix
916  boost::shared_ptr<CrackFrontElement> feEnergy; ///< Integrate energy
917  boost::shared_ptr<ConstrainMatrixCtx>
918  projSurfaceCtx; ///< Data structure to project on the body surface
919  boost::shared_ptr<ConstrainMatrixCtx>
920  projFrontCtx; ///< Data structure to project on crack front
921  boost::shared_ptr<MWLSApprox>
922  mwlsApprox; ///< Data struture for MWLS approximation
923  /// Surface elment to calculate tractions in material space
924  boost::shared_ptr<NeumannForcesSurface::MyTriangleFE> feMaterialTraction;
925  boost::shared_ptr<DirichletSpatialPositionsBc>
926  spatialDirichletBc; ///< apply Dirichlet BC to sparial positions
927  boost::shared_ptr<AnalyticalDirichletBC::DirichletBC>
928  analyticalDirichletBc; ///< apply analytical Dirichlet BC to sparial
929  ///< positions
930  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
931  surfaceForces; ///< assemble surface forces
932  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
933  surfacePressure; ///< assemble surface pressure
934  boost::shared_ptr<boost::ptr_map<string, EdgeForce>>
935  edgeForces; ///< assemble edge forces
936  boost::shared_ptr<boost::ptr_map<string, NodalForce>>
937  nodalForces; ///< assemble nodal forces
938  boost::shared_ptr<FEMethod> assembleFlambda; ///< assemble F_lambda vector
939  boost::shared_ptr<FEMethod> zeroFlambda; ///< assemble F_lambda vector
940 
941  boost::shared_ptr<CrackFrontElement>
942  feCouplingElasticLhs; ///< FE instance to assemble coupling terms
943  boost::shared_ptr<CrackFrontElement>
944  feCouplingMaterialLhs; ///< FE instance to assemble coupling terms
945 
946  boost::shared_ptr<Smoother> smootherFe;
947  boost::shared_ptr<Smoother::MyVolumeFE>
948  feSmootherRhs; ///< Integrate smoothing operators
949  boost::shared_ptr<Smoother::MyVolumeFE>
950  feSmootherLhs; ///< Integrate smoothing operators
951  boost::shared_ptr<VolumeLengthQuality<double>> volumeLengthDouble;
952  boost::shared_ptr<VolumeLengthQuality<adouble>> volumeLengthAdouble;
953 
954  boost::shared_ptr<
956  tangentConstrains; ///< Constrains crack front in tangent directiona
957  boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
959  boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
961  map<int, boost::shared_ptr<SurfaceSlidingConstrains>> surfaceConstrain;
962  map<int, boost::shared_ptr<EdgeSlidingConstrains>> edgeConstrains;
963 
964  boost::shared_ptr<BothSurfaceConstrains> bothSidesConstrains;
965  boost::shared_ptr<DirichletFixFieldAtEntitiesBc> fixMaterialEnts;
966 
967  boost::shared_ptr<GriffithForceElement> griffithForceElement;
968  boost::shared_ptr<GriffithForceElement::MyTriangleFE> feGriffithForceRhs;
969  boost::shared_ptr<GriffithForceElement::MyTriangleFE> feGriffithForceLhs;
970  boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrainsDelta>
972  boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrains>
974  boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrains>
976 
977  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
978  analiticalSurfaceElement; ///< Finite element to integrate tractions given
979  ///< by analytical formula
980  boost::shared_ptr<FaceElementForcesAndSourcesCore> feSpringLhsPtr;
981  boost::shared_ptr<FaceElementForcesAndSourcesCore> feSpringRhsPtr;
982  /**@}*/
983 
984  /** \name Entity ranges */
985 
986  /**@{*/
987 
988  Range bitEnts;
989  Range bitProcEnts;
990  Range bodySkin;
991  Range crackFaces;
994  Range crackFront;
998 
999  EntityHandle arcLengthVertex; ///< Vertex associated with dof for arc-length control
1000 
1001  /**@}*/
1002 
1003  /** \name Auxiliary data */
1004 
1005  /**@{*/
1006 
1008  boost::shared_ptr<VectorDouble>
1009  rhoAtGaussPt; ///< density at integration points
1010  boost::shared_ptr<MatrixDouble>
1011  rhoGradAtGaussPt; ///< gradient of density at integration point
1012 
1015  double loadScale;
1016 
1017  double maxG1;
1018  double maxJ;
1020  map<EntityHandle, double> mapG1; ///< hashmap of g1 - release energy at nodes
1021  map<EntityHandle, double> mapG3; ///< hashmap of g3 - release energy at nodes
1022  map<EntityHandle, double> mapJ; ///< hashmap of J - release energy at nodes
1023  map<EntityHandle, VectorDouble3>
1024  mapMatForce; ///< hashmap of material force at nodes
1025  map<EntityHandle, VectorDouble3>
1026  mapGriffith; ///< hashmap of Griffith energy at nodes
1027 
1028  map<EntityHandle, double> mapSmoothingForceFactor;
1029 
1030  bool doSurfaceProjection; ///< If crack surface is immersed in the body do not
1031  ///< projection
1032 
1033  double smootherAlpha; ///< Controls mesh smoothing
1034  double initialSmootherAlpha; ///< Controls mesh smoothing, set at start of
1035  ///< analysis
1036  double smootherGamma; ///< Controls mesh smoothing
1037 
1038  /**@}*/
1039 
1040  /** \name arc-length options */
1041 
1042  /**@{*/
1043 
1044  double arcAlpha;
1045  double arcBeta;
1046  double arcS;
1047 
1048  boost::shared_ptr<ArcLengthMatShell> arcMatCtx;
1049  boost::shared_ptr<ArcLengthCtx> arcCtx;
1050  boost::shared_ptr<PCArcLengthCtx> pcArcLengthCtx;
1051 
1052  /**@}*/
1053 
1054  /** \name Data for bone */
1055 
1056  /**@{*/
1057 
1058  double rHo0; ///< Reference density if bone is analyzed
1059  double nBone; ///< Exponent parameter in bone density
1060 
1061  /**@}*/
1062 
1063  /** \name Interfaces */
1064 
1065  /**@{*/
1066 
1067  const boost::scoped_ptr<UnknownInterface> cpSolversPtr;
1068  const boost::scoped_ptr<UnknownInterface> cpMeshCutPtr;
1069 
1070  /**@}*/
1071 
1072  /** \name Cut/Refine/Split options */
1073 
1074  /**@{*/
1075 
1076  PetscBool doCutMesh;
1078 
1079  /**@}*/
1080 };
1081 
1082 } // namespace FractureMechanics
1083 
1084 #endif // __CRACK_PROPAGATION_HPP__
boost::shared_ptr< MatrixDouble > rhoGradAtGaussPt
gradient of density at integration point
MoFEMErrorCode updateMaterialFixedNode(const bool fix_front, const bool fix_small_g, const bool debug=false)
Update fixed nodes.
MoFEMErrorCode createMaterialDM(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 declareEdgeFE(std::string fe_name, const BitRefLevel bit, const BitRefLevel mask, const bool proc_only=true, const int verb=QUIET, const bool debug=false)
EntityHandle arcLengthVertex
Vertex associated with dof for arc-length control.
boost::shared_ptr< NeumannForcesSurface::MyTriangleFE > feMaterialTraction
Surface elment to calculate tractions in material space.
map< EntityHandle, VectorDouble3 > mapGriffith
hashmap of Griffith energy at nodes
boost::shared_ptr< FEMethod > getFrontArcLengthControl()
double smootherGamma
Controls mesh smoothing.
MoFEMErrorCode testMWLSJacobians(const BitRefLevel bit, const BitRefLevel mask)
test MWLS
MoFEMErrorCode addSmoothingFEInstancesToSnes(DM dm, const bool fix_crack_front, const int verb=QUIET, const bool debug=false)
add softening elements instances to SNES
boost::shared_ptr< ArcLengthCtx > & getArcCtx()
Deprecated interface functions.
MoFEMErrorCode buildBothSidesFieldId(const BitRefLevel bit=BitRefLevel().set(0), 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.
boost::shared_ptr< ConstrainMatrixCtx > projFrontCtx
Data structure to project on crack front.
MoFEMErrorCode deleteEntities(const int verb=QUIET, const bool debug=false)
Vec F
MoFEMErrorCode assembleCouplingForcesDM(DM dm, const int verb=QUIET, const bool debug=false)
assemble coupling element instances
boost::shared_ptr< CrackFrontElement > feLhs
Integrate elastic FE.
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 setCoordsFromField(const std::string field_name="MESH_NODE_POSITIONS")
set coordinates from field
MoFEMErrorCode setFieldFromCoords(const std::string field_name)
set field from node positions
MoFEMErrorCode addElasticFEInstancesToSnes(DM dm, Mat m, Vec q, Vec f, boost::shared_ptr< FEMethod > arc_method=boost::shared_ptr< FEMethod >(), const int verb=QUIET, const bool debug=false)
double rHo0
Reference density if bone is analyzed.
std::string tagName
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MoFEMErrorCode writeCrackFont(const BitRefLevel &bit, const int step)
MoFEMErrorCode addPropagationFEInstancesToSnes(DM dm, boost::shared_ptr< FEMethod > &arc_method, const std::vector< int > &surface_ids, const int verb=QUIET, const bool debug=false)
add finite element to SNES for crack propagation problem
MoFEMErrorCode unsetSingularElementMatrialPositions()
remove singularity
MoFEMErrorCode getInterfaceVersion(Version &version) const
double griffithE
Griffith stability parameter.
MoFEMErrorCode buildProblemFields(const BitRefLevel &bit1, const BitRefLevel &mask1, const BitRefLevel &bit2, const int verb=QUIET, const bool debug=false)
Build problem fields.
boost::shared_ptr< BothSurfaceConstrains > bothSidesConstrains
MoFEMErrorCode zeroLambdaFields()
Zero fields with lagrange multipliers.
int residualStressBlock
Block on which residual stress is applied.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
const boost::scoped_ptr< UnknownInterface > cpMeshCutPtr
boost::shared_ptr< ArcLengthMatShell > arcMatCtx
double nBone
Exponent parameter in bone density.
map< EntityHandle, VectorDouble3 > mapMatForce
hashmap of material force at nodes
static constexpr int approx_order
boost::shared_ptr< CrackFrontElement > feEnergy
Integrate energy.
MoFEMErrorCode partitionMesh(BitRefLevel bit1, BitRefLevel bit2, int verb=QUIET, const bool debug=false)
partotion mesh
boost::shared_ptr< FaceElementForcesAndSourcesCore > feSpringRhsPtr
AdolcTags
Tapes numbers used by ADOL-C.
base class for all interface classes
boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > crackOrientation
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode calculateSmoothingForceFactor(const int verb=QUIET, const bool debug=true)
MoFEMErrorCode assembleMaterialForcesDM(DM dm, const int verb=QUIET, const bool debug=false)
create material element instance
MoFEMErrorCode buildArcLengthField(const BitRefLevel bit, const bool build_fields=true, const int verb=QUIET)
Declate field for arc-length.
MoFEMErrorCode testMWLSJacobiansWithDensity(const BitRefLevel bit, const BitRefLevel mask)
test MWLS
Class implemented by user to detect face orientation.
boost::shared_ptr< PCArcLengthCtx > pcArcLengthCtx
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrainsDelta > feGriffithConstrainsDelta
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
Data structure to exchange data between mofem and User Loop Methods on entities.It allows to exchange...
MoFEMErrorCode createSurfaceProjectionDM(DM *dm, 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)
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 buildEdgeFields(const BitRefLevel bit=BitRefLevel().set(0), 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
MoFEMErrorCode postProcess()
function is run at the end of loop
boost::shared_ptr< Smoother::MyVolumeFE > feSmootherLhs
Integrate smoothing operators.
MoFEMErrorCode cleanSingularElementMatrialPositions()
set maetrial field order to one
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
MoFEMErrorCode projectGriffithForce(DM dm, Vec f_griffith, Vec f_griffith_proj, const int verb=QUIET, const bool debug=false)
project Griffith forces
MoFEMErrorCode assembleSmootherForcesDM(DM dm, const std::vector< int > ids, const int verb=QUIET, const bool debug=false)
create smoothing element instance
MoFEMErrorCode calculateElasticEnergy(DM dm, const std::string msg="")
assemble elastic part of matrix, by running elastic finite element instance
MoFEMErrorCode setSingularElementMatrialPositions(const int verb=QUIET)
set singular dofs of material field (on edges adjacent to crack front) from geometry
MoFEMErrorCode solveElasticDM(DM dm, SNES snes, Mat m, Vec q, Vec f, bool snes_set_up, Mat *shell_m)
solve elastic problem
MoFEMErrorCode buildElasticFields(const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool build_fields=true, const int verb=QUIET)
Declate fields for elastic analysis.
map< EntityHandle, double > mapG3
hashmap of g3 - release energy at nodes
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
std::string mwlsApproxFile
Name of file with internal stresses.
const boost::scoped_ptr< UnknownInterface > cpSolversPtr
MoFEM::Interface & mField
map< EntityHandle, double > mapG1
hashmap of g1 - release energy at nodes
MoFEMErrorCode savePositionsOnCrackFrontDM(DM dm, Vec q, const int verb=QUIET, const bool debug=false)
boost::shared_ptr< VolumeLengthQuality< double > > volumeLengthDouble
PetscBool propagateCrack
If true crack propagation is calculated.
MoFEMErrorCode resolveShared(const Range &tets, Range &proc_ents, const int verb=QUIET, const bool debug=false)
resolve shared entities
boost::shared_ptr< NonlinearElasticElement > materialFe
boost::shared_ptr< boost::ptr_map< string, NodalForce > > nodalForces
assemble nodal forces
boost::shared_ptr< AnalyticalDirichletBC::DirichletBC > analyticalDirichletBc
boost::shared_ptr< boost::ptr_map< string, EdgeForce > > edgeForces
assemble edge forces
boost::shared_ptr< FaceElementForcesAndSourcesCore > feSpringLhsPtr
boost::shared_ptr< DirichletSpatialPositionsBc > spatialDirichletBc
apply Dirichlet BC to sparial positions
MoFEMErrorCode declareFrontFE(const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
boost::shared_ptr< ObosleteUsersModules::TangentWithMeshSmoothingFrontConstrain > tangentConstrains
Constrains crack front in tangent directiona.
boost::shared_ptr< GriffithForceElement > griffithForceElement
MoFEMErrorCode declareSmoothingFE(const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
declare mesh smoothing finite elements
Structure for DM for multi-grid via approximation orders.
MoFEMErrorCode readMedFile()
read mesh file
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrains > feGriffithConstrainsRhs
MoFEMErrorCode buildCrackFrontFieldId(const BitRefLevel bit=BitRefLevel().set(0), const bool build_fields=true, const int verb=QUIET, const bool debug=false)
declare crack surface files
boost::shared_ptr< ArcLengthCtx > arcCtx
CHKERRQ(ierr)
MoFEMErrorCode setSingularDofs(const string field_name, const int verb=QUIET)
set singular dofs (on edges adjacent to crack front) from geometry
MoFEMErrorCode destroyDMs(DM *dm_elastic, DM *dm_material, DM *dm_crack_propagation, DM *dm_material_forces, DM *dm_surface_projection, DM *dm_crack_srf_area)
destroy DMs
MoFEMErrorCode createCrackPropagationDM(DM *dm, const std::string prb_name, DM dm_elastic, 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.
FTensor::Index< 'm', 2 > m
Definition: PlasticOps.hpp:67
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrains > feGriffithConstrainsLhs
double * aRray
std::string mwlsStressTagName
Name of tag with internal stresses.
CrackPropagation(MoFEM::Interface &m_field, const int approx_order=2, const int geometry_order=1)
Constructor of Crack Propagation, prints fracture module version and registers the three interfaces u...
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
boost::shared_ptr< ConstrainMatrixCtx > projSurfaceCtx
Data structure to project on the body surface.
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< DofEntity > arcLengthDof
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)
boost::shared_ptr< AnalyticalDirichletBC::DirichletBC > getAnalyticalDirichletBc()
std::string mwlsRhoTagName
Name of tag with density.
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
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
double smootherAlpha
Controls mesh smoothing.
MoFEMErrorCode calculateSurfaceProjectionMatrix(DM dm_front, DM dm_project, const std::vector< int > &ids, const int verb=QUIET, const bool debug=false)
assemble projection matrices
MoFEMErrorCode buildSurfaceFields(const BitRefLevel bit=BitRefLevel().set(0), 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 tetsingReleaseEnergyCalculation()
This is run with ctest.
boost::shared_ptr< Smoother > smootherFe
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > analiticalSurfaceElement
PostProcVertexMethod(MoFEM::Interface &m_field, Vec f, std::string tag_name)
MoFEMErrorCode assembleElasticDM(DM dm, Mat m, Vec q, Vec f, const int verb=QUIET, const bool debug=false)
create elastic finite element instance for spatial assembly
boost::shared_ptr< MWLSApprox > & getMWLSApprox()
#define CHKERR
Inline error check.
Definition: definitions.h:601
MoFEMErrorCode createCrackFrontAreaDM(DM *dm, DM dm_material, const std::string prb_name, const bool verb=QUIET)
create DM to calculate Griffith energy
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfacePressure
assemble surface pressure
boost::shared_ptr< NonlinearElasticElement > elasticFe
Constrains material displacement on both sides.
MoFEMErrorCode createBcDM(DM *dm, const std::string prb_name, const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set())
Create problem to calculate boundary conditions.
MoFEMErrorCode solvePropagationDM(DM dm, DM dm_elastic, SNES snes, Mat m, Vec q, Vec f, bool snes_set_up, Mat *shell_m, Mat *m_elastic)
solve crack propagation problem
bool onlyHooke
True if only Hooke material is applied.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
static const MOFEMuuid IDD_MOFEMCrackPropagationInterface
MoFEMErrorCode operator()()
function is run for every finite element
boost::shared_ptr< Smoother::MyVolumeFE > feSmootherRhs
Integrate smoothing operators.
MoFEMErrorCode setMaterialPositionFromCoords()
set material field from nodes
MoFEMErrorCode declareBothSidesFE(const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
double griffithR
Griffith regularisation parameter.
MoFEMErrorCode getArcLengthDof()
set pointer to arc-length DOF
MoFEMErrorCode saveEachPart(const std::string prefix, const Range &ents)
Save entities on ech processor.
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.
map< EntityHandle, double > mapJ
hashmap of J - release energy at nodes
MoFEMErrorCode resolveSharedBitRefLevel(const BitRefLevel bit, const int verb=QUIET, const bool debug=false)
resole shared entities by bit level
boost::shared_ptr< MWLSApprox > mwlsApprox
MoFEMErrorCode declareMaterialFE(const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
declare material finite elements
boost::shared_ptr< CrackFrontElement > feMaterialRhs
Integrate material stresses, assemble vector.
boost::shared_ptr< FEMethod > assembleFlambda
assemble F_lambda vector
MoFEMErrorCode readConfigFileOptions()
read configuration file
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
MoFEMErrorCode declareArcLengthFE(const BitRefLevel bits, const int verb=QUIET)
create arc-length element entity and declare elemets
MoFEMErrorCode getOptions()
Get options form command line.
MoFEMErrorCode createMaterialForcesDM(DM *dm, DM dm_material, const std::string prb_name, const int verb=QUIET)
Create DM for calculation of material forces (sub DM of DM material)
MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const EntityHandle face, const BitRefLevel bit)
MoFEMErrorCode postProcessDM(DM dm, const int step, const std::string fe_name, const bool approx_internal_stress)
post-process results for elastic solution
ArcLengthSnesCtx(MoFEM::Interface &m_field, const std::string &problem_name, boost::shared_ptr< ArcLengthCtx > &arc_ptr, DMMGViaApproxOrdersCtx *dm_ctx, DMMGViaApproxOrdersCtx *dm_ctx_elastic=NULL)
PetscBool onlyHookeFromOptions
True if only Hooke material is applied.
map< int, boost::shared_ptr< EdgeSlidingConstrains > > edgeConstrains
boost::shared_ptr< FEMethod > zeroFlambda
assemble F_lambda vector
boost::shared_ptr< GriffithForceElement::MyTriangleFE > feGriffithForceRhs
map< EntityHandle, double > mapSmoothingForceFactor
boost::shared_ptr< DirichletFixFieldAtEntitiesBc > fixMaterialEnts
MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const FEMethod *fe_method_ptr)
MoFEMErrorCode declareExternalForcesFE(const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
operator to post-process results on crack front nodes
boost::shared_ptr< GriffithForceElement::MyTriangleFE > feGriffithForceLhs
boost::shared_ptr< CrackFrontElement > feRhs
Integrate elastic FE.
boost::shared_ptr< VolumeLengthQuality< adouble > > volumeLengthAdouble
boost::shared_ptr< CrackFrontElement > feCouplingMaterialLhs
FE instance to assemble coupling terms.
boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > skinOrientation
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< CrackFrontElement > feCouplingElasticLhs
FE instance to assemble coupling terms.
MoFEMErrorCode createDMs(BitRefLevel bit1, BitRefLevel bit2, DM *dm_elastic, DM *dm_material, DM *dm_crack_propagation, DM *dm_material_forces, DM *dm_surface_projection, 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.
MoFEMErrorCode addMaterialFEInstancesToSnes(DM dm, const bool fix_crack_front, const int verb=QUIET, const bool debug=false)
add material elements instances to SNES
int postProcLevel
level of postprocessing (amount of output files)
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
map< int, boost::shared_ptr< SurfaceSlidingConstrains > > surfaceConstrain
MoFEMErrorCode setSpatialPositionFromCoords()
set spatial field from nodes
MoFEMErrorCode createElasticDM(DM *dm, const std::string prb_name, const BitRefLevel bit=BitRefLevel().set(0), const BitRefLevel mask=BitRefLevel().set())
Create elastic problem DM.
boost::shared_ptr< VectorDouble > rhoAtGaussPt
density at integration points
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfaceForces
assemble surface forces
boost::shared_ptr< CrackFrontElement > feMaterialLhs
Integrate material stresses, assemble matrix.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Getting interface of core database.
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
MoFEMErrorCode createProblemDataStructures(const BitRefLevel bit0, const BitRefLevel bit1, const BitRefLevel bit2, const std::vector< int > surface_ids, const int verb=QUIET, const bool debug=false)
Construct problem data structures.
Tag tH