v0.10.0
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  * @brief Names of tangent matrices tests
60  *
61  */
67 };
68 /**
69  * Crack propagation data and methods
70  * \nosubgrouping
71  *
72  * - Functions strarting with build, for example buildElastic,
73  * add and build fields.
74  *
75  * - Functions starting with declare, for example declareElasticFE,
76  * add and build finite elements.
77  *
78  * - Functions starting with assemble, for example assembleElasticDM,
79  * create finite element instances and operators.
80  *
81  * - Functions starting with create, for example createElasticDM, create
82  * problems.
83  *
84  */
86 
87  static bool debug;
88 
89  MoFEMErrorCode getInterfaceVersion(Version &version) const {
91  version = Version(FM_VERSION_MAJOR, FM_VERSION_MINOR, FM_VERSION_BUILD);
93  }
94 
95  /**
96  * \brief Getting interface of core database
97  * @param uuid unique ID of interface that can be either for
98  * CrackPropagation, CPSolvers or CPMeshCut interface
99  * @param iface returned pointer to interface
100  * @return error code
101  */
102  MoFEMErrorCode query_interface(const MOFEMuuid &uuid,
103  UnknownInterface **iface) const;
104 
105  MPI_Comm moabCommWorld;
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 
127  std::string mwlsApproxFile; ///< Name of file with internal stresses
128  std::string mwlsStressTagName; ///< Name of tag with internal stresses
129  std::string mwlsRhoTagName; ///< Name of tag with density
130  int residualStressBlock; ///< Block on which residual stress is applied
131  int densityMapBlock; ///< Block on which density is applied (Note: if used
132  ///< with internal stress, it has to be the same ID!)
133  std::string medFileName;
134  std::string configFile;
136 
138 
140 
141  );
142 
143  PetscBool
144  addAnalyticalInternalStressOperators; ///< Add analytical internal stress
145  ///< operators for MWLS test
146 
147  int postProcLevel; ///< level of postprocessing (amount of output files)
149 
152 
153  /**
154  * Map associating string literals to bit refinement levels
155  * The following convention is currently used:
156  * "mesh_cut" : level 0
157  * "spatial_domain" : level 1
158  * "material_domain" : level 2
159  * */
160  std::map<std::string, BitRefLevel> mapBitLevel;
161 
162  /**
163  * \brief Constructor of Crack Propagation, prints fracture module version and
164  registers the three interfaces
165  * used to solve fracture problem, i.e. CrackPropagation, CPSolvers and
166  * CPMeshCut
167  * @param m_field interface to MoAB database
168  * @param approx_order order of approximation space, default is 2
169  * @param geometry_order order of geometry, default is 1 (can be up to 2)
170  */
171  CrackPropagation(MoFEM::Interface &m_field, const int approx_order = 2,
172  const int geometry_order = 1);
173 
174  virtual ~CrackPropagation();
175 
176  /**
177  * \brief Get options form command line
178  * @return error code
179  */
181 
182  /**
183  * \brief This is run with ctest
184  * @return error code
185  */
187 
188  inline int &getNbLoadSteps() { return nbLoadSteps; }
189  inline int getNbLoadSteps() const { return nbLoadSteps; }
190  inline double &getLoadScale() { return loadScale; }
191  inline double getLoadScale() const { return loadScale; }
192  inline int &getNbCutSteps() { return nbCutSteps; }
193  inline int getNbCutSteps() const { return nbCutSteps; }
194 
195  /// \brief read mesh file
197 
198  /// \brief read configuration file
200 
201  /**
202  * @brief Save entities on ech processor.
203  *
204  * Usually used for debugging, to check if meshes are consistent on all
205  * processors.
206  *
207  * @param prefix file name prefix
208  * @param ents entities to save
209  * @return MoFEMErrorCode error code
210  */
211  MoFEMErrorCode saveEachPart(const std::string prefix, const Range &ents);
212 
213  /** \name Resolve shared entities and partition mesh */
214 
215  /**@{*/
216 
217  /// \brief partotion mesh
219  int verb = QUIET, const bool debug = false);
220 
221  /// \brief resolve shared entities
222  MoFEMErrorCode resolveShared(const Range &tets, Range &proc_ents,
223  const int verb = QUIET,
224  const bool debug = false);
225 
226  /// \brief resole shared entities by bit level
228  const int verb = QUIET,
229  const bool debug = false);
230 
231  /**@}*/
232 
233  /** \name Set bit level for material problem */
234 
235  /**@{*/
236 
237  /**
238  * @brief Set bit ref level for entities adjacent to crack front
239  *
240  * @param from_bit
241  * @param bit
242  * @param nb_levels adjacency bridge levels
243  * @return MoFEMErrorCode
244  */
246  const int nb_levels = 2,
247  const bool debug = false);
248 
249  /**@}*/
250 
251  /** \name Build problem */
252 
253  /**@{*/
254 
255  /**
256  * \brief Build problem fields
257  * @param bit1 bit level for computational domain
258  * @param mask1 bit mask
259  * @param bit2 bit level for elements near crack front
260  * @param mask2 bit mask 2
261  * @param verb verbosity level
262  * @param debug
263  * @return error code
264  *
265  */
267  const BitRefLevel &mask1,
268  const BitRefLevel &bit2,
269  const int verb = QUIET,
270  const bool debug = false);
271 
272  /**
273  * \brief Build problem finite elements
274  * @param bit1 bit level for computational domain
275  * @param bit2 bit level for elements near crack front
276  * @param surface_ids meshsets Ids on body skin
277  * @param verb verbosity level
278  * @param debug
279  * @return error code
280  *
281  */
283  std::vector<int> surface_ids,
284  const int verb = QUIET,
285  const bool debug = false);
286 
287  /**
288  * \brief Construct problem data structures
289  * @param surface_ids surface_ids meshsets Ids on body skin
290  * @param verb verbosity level
291  * @param debug
292  * @return error code
293  */
294  MoFEMErrorCode createProblemDataStructures(const std::vector<int> surface_ids,
295  const int verb = QUIET,
296  const bool debug = false);
297 
298  /**
299  * \brief Crate DMs for all problems and sub problems
300  * @param dm_elastic elastic DM
301  * @param dm_material material DM
302  * @param dm_crack_propagation composition of elastic DM and material DM
303  * @param dm_material_forces used to calculate material forces, sub-dm of
304  * dm_crack_propagation
305  * @param dm_surface_projection dm to project material forces on surfaces,
306  * sub-dm of dm_crack_propagation
307  * @param dm_crack_srf_area dm to project material forces on crack
308  * surfaces, sub-dm of dm_crack_propagation
309  * @param surface_ids IDs of surface meshsets
310  * @param fe_surf_proj_list name of elements on surface elements
311  * @return error code
312  */
313  MoFEMErrorCode createDMs(DM *dm_elastic, DM *dm_material,
314  DM *dm_crack_propagation, DM *dm_material_forces,
315  DM *dm_surface_projection, DM *dm_crack_srf_area,
316  std::vector<int> surface_ids,
317  std::vector<std::string> fe_surf_proj_list);
318 
319  /// \brief destroy DMs
320  MoFEMErrorCode destroyDMs(DM *dm_elastic, DM *dm_material,
321  DM *dm_crack_propagation, DM *dm_material_forces,
322  DM *dm_surface_projection, DM *dm_crack_srf_area);
323 
324  MoFEMErrorCode deleteEntities(const int verb = QUIET,
325  const bool debug = false);
326 
327  /**@}*/
328 
329  /** \name Build fields */
330 
331  /**@{*/
332 
333  /**
334  * \brief Declate fields for elastic analysis
335  * @param bit bit refinement level
336  * @return error code
337  */
339  const BitRefLevel bit, const BitRefLevel mask = BitRefLevel().set(),
340  const bool proc_only = true, const bool build_fields = true,
341  const int verb = QUIET, const bool debug = false);
342 
343  /**
344  * \brief Declate field for arc-length
345  * @param bit bit refinement level
346  * @return error code
347  */
349  const bool build_fields = true,
350  const int verb = QUIET);
351 
352  /// \brief build fields with Lagrange multipliers to constrain surfaces
354  const bool proc_only = true,
355  const bool build_fields = true,
356  const int verb = QUIET,
357  const bool debug = false);
358 
359  /// \brief build fields with Lagrange multipliers to constrain edges
361  const bool proc_only = true,
362  const bool build_fields = true,
363  const int verb = QUIET,
364  const bool debug = false);
365 
366  /// \brief declare crack surface files
368  const bool proc_only = true,
369  const bool build_fields = true,
370  const int verb = QUIET,
371  const bool debug = false);
372 
373  /// \brief declare crack surface files
375  const bool build_fields = true,
376  const int verb = QUIET,
377  const bool debug = false);
378 
379  /**
380  * \brief Lagrange multipliers field which constrains material displacements
381  *
382  * Nodes on both sides are constrainded to have the same material
383  * displacements
384  *
385  */
387  const bool proc_only = false,
388  const bool build_fields = true,
389  const int verb = QUIET,
390  const bool debug = false);
391 
392  /** \brief Zero fields with lagrange multipliers
393  */
395 
396  /**@}*/
397 
398  /** \name Build finite elements */
399 
400  /**@{*/
401 
402  /// \brief declare elastic finite elements
404  declareElasticFE(const BitRefLevel bit1, const BitRefLevel mask1,
405  const BitRefLevel bit2, const BitRefLevel mask2,
406  const bool add_forces = true, const bool proc_only = true,
407  const int verb = QUIET);
408 
409  /// \brief create arc-length element entity and declare elemets
411  const int verb = QUIET);
412 
415  const BitRefLevel mask = BitRefLevel().set(),
416  const bool proc_only = true);
417 
418  /** \brief Declare FE for pressure BC in ALE formulation (in material domain)
419  * @param bit material domain bit ref level
420  * @param mask bit ref level mask for the problem
421  * @return error code
422  */
425  const BitRefLevel mask = BitRefLevel().set(),
426  const bool proc_only = true);
427 
428  /** \brief Declare FE for pressure BC in ALE formulation (in material domain)
429  * @param bit material domain bit ref level
430  * @param mask bit ref level mask for the problem
431  * @return error code
432  */
435  const BitRefLevel mask = BitRefLevel().set(),
436  const bool proc_only = true);
437 
438  /// \brief declare material finite elements
440  const BitRefLevel mask = BitRefLevel().set(),
441  const bool proc_only = true,
442  const bool verb = QUIET);
443 
445  const BitRefLevel mask = BitRefLevel().set(),
446  const bool proc_only = true,
447  const bool verb = QUIET);
448 
451  const BitRefLevel mask = BitRefLevel().set(),
452  const bool proc_only = true, const bool verb = QUIET);
453 
454  /// \brief declare mesh smoothing finite elements
457  const BitRefLevel mask = BitRefLevel().set(),
458  const bool proc_only = true, const bool verb = QUIET);
459 
460  /// \brief declare surface sliding elements
461  MoFEMErrorCode declareSurfaceFE(std::string fe_name, const BitRefLevel bit,
462  const BitRefLevel mask,
463  const std::vector<int> &ids,
464  const bool proc_only = true,
465  const int verb = QUIET,
466  const bool debug = false);
467 
468  // \brief declare edge sliding elements
469  MoFEMErrorCode declareEdgeFE(std::string fe_name, const BitRefLevel bit,
470  const BitRefLevel mask,
471  const bool proc_only = true,
472  const int verb = QUIET,
473  const bool debug = false);
474 
476  declareCrackSurfaceFE(std::string fe_name, const BitRefLevel bit,
477  const BitRefLevel mask, const bool proc_only = true,
478  const int verb = QUIET, const bool debug = false);
479 
480  /**@}*/
481 
482  /** \name Create DMs & problems */
483 
484  /**@{*/
485 
486  /** \brief Create DM by composition of elastic DM and material DM
487  * @param dm composition of elastic DM and material DM
488  * @param prb_name problem name
489  * @param dm_elastic elastic DM
490  * @param dm_material material DM
491  * @param bit domain bit ref level
492  * @param mask mask ref level for problem
493  * @return error code
494  */
496  createCrackPropagationDM(DM *dm, const std::string prb_name, DM dm_elastic,
497  DM dm_material, const BitRefLevel bit,
498  const BitRefLevel mask,
499  const std::vector<std::string> fe_list);
500 
501  /** \brief Create elastic problem DM
502  * @param dm elastic DM
503  * @param prb_name problem name
504  * @param bit domain bit ref level
505  * @param mask mask ref level for problem
506  * @return error code
507  */
508  MoFEMErrorCode createElasticDM(DM *dm, const std::string prb_name,
509  const BitRefLevel bit,
510  const BitRefLevel mask = BitRefLevel().set());
511 
512  /** \brief Create problem to calculate boundary conditions
513  * @param dm newly created boundary condition DM
514  * @param prb_name problem name
515  * @param bit auxiliary bit ref level
516  * @param mask mask ref level for problem
517  * @return error code
518  * \note BcDM is only useed for testing, to enforce boundary conditions
519 staisifing anlylitical solution to which numerical solution is compared.
520  */
521  MoFEMErrorCode createBcDM(DM *dm, const std::string prb_name,
522  const BitRefLevel bit,
523  const BitRefLevel mask = BitRefLevel().set());
524 
525  /** \brief Create DM fto calculate material problem
526  * @param dm used to calculate material forces, sub-dm of
527  * dm_crack_propagation
528  * @param prb_name problem name
529  * @param bit material domain
530  * @param mask dof mask ref level for problem
531  * @param fe_list name of elements on surface elements
532  * @param debug flag for debugging
533  * @return error code
534  */
535  MoFEMErrorCode createMaterialDM(DM *dm, const std::string prb_name,
536  const BitRefLevel bit, const BitRefLevel mask,
537  const std::vector<std::string> fe_list,
538  const bool debug = false);
539 
540  /** \brief Create DM for calculation of material forces (sub DM of DM
541  * material)
542  * @param dm used to calculate material forces, sub-dm of
543  * dm_crack_propagation
544  * @param prb_name problem name
545  * @param verb compilation parameter determining the amount
546  * of information printed
547  * @return error code
548  */
549  MoFEMErrorCode createMaterialForcesDM(DM *dm, DM dm_material,
550  const std::string prb_name,
551  const int verb = QUIET);
552 
553  /** \brief create DM to calculate projection matrices (sub DM of DM
554  * material)
555  * @param dm DM to project material forces on surfaces,
556  * sub-dm of dm_crack_propagation
557  * @param dm_material material DM
558  * @param prb_name problem name
559  * @param surface_ids IDs of surface meshsets
560  * @param fe_list name of elements on surface elements
561  * @param verb compilation parameter determining the amount
562  * of information printed
563  * @return error code
564  */
566  createSurfaceProjectionDM(DM *dm, DM dm_material, const std::string prb_name,
567  const std::vector<int> surface_ids,
568  const std::vector<std::string> fe_list,
569  const int verb = QUIET);
570 
571  /** \brief create DM to calculate Griffith energy
572  * @param dm dm to project material forces on crack
573  * surfaces, sub-dm of dm_crack_propagation
574  * @param dm_material material DM
575  * @param prb_name problem name
576  * @param verb compilation parameter determining the amount
577  * of information printed
578  * @return error code
579  */
580  MoFEMErrorCode createCrackFrontAreaDM(DM *dm, DM dm_material,
581  const std::string prb_name,
582  const bool verb = QUIET);
583 
584  /**@}*/
585 
586  /** \name Create finite element instances */
587 
588  /**@{*/
589 
590  /** \brief create elastic finite element instance for spatial assembly
591  * @param dm elastic DM
592  * @param m stiffness sub-matrix for elastic DM
593  * @param q DoF sub-vector for elastic DM
594  * @param f Rhs sub-vector for material DM
595  * @param verb compilation parameter determining the amount
596  * of information printed (not in use here)
597  * @param debug flag for debugging
598  * @return error code
599  */
600  MoFEMErrorCode assembleElasticDM(DM dm, Mat m, Vec q, Vec f,
601  const int verb = QUIET,
602  const bool debug = false);
603 
604  // \brief add elastic element instances to SNES
606  DM dm, Mat m, Vec q, Vec f,
607  boost::shared_ptr<FEMethod> arc_method = boost::shared_ptr<FEMethod>(),
608  const int verb = QUIET, const bool debug = false);
609 
610  /// \brief create material element instance
611  MoFEMErrorCode assembleMaterialForcesDM(DM dm, const int verb = QUIET,
612  const bool debug = false);
613 
614  /// \brief create smoothing element instance
615  MoFEMErrorCode assembleSmootherForcesDM(DM dm, const std::vector<int> ids,
616  const int verb = QUIET,
617  const bool debug = false);
618 
619  /// \brief assemble coupling element instances
620  MoFEMErrorCode assembleCouplingForcesDM(DM dm, const int verb = QUIET,
621  const bool debug = false);
622 
623  /// \brief add material elements instances to SNES
624  MoFEMErrorCode addMaterialFEInstancesToSnes(DM dm, const bool fix_crack_front,
625  const int verb = QUIET,
626  const bool debug = false);
627 
628  /// \brief add softening elements instances to SNES
630  const bool fix_crack_front,
631  const int verb = QUIET,
632  const bool debug = false);
633 
634  /// \brief add finite element to SNES for crack propagation problem
636  DM dm, boost::shared_ptr<FEMethod> &arc_method,
637  const std::vector<int> &surface_ids, const int verb = QUIET,
638  const bool debug = false);
639 
640  /** \brief test LHS Jacobians
641  * @param bit bit level
642  * @param mask bitref mask
643  * @param test_nb name of the test
644  */
645  MoFEMErrorCode testJacobians(const BitRefLevel bit, const BitRefLevel mask,
646  tangent_tests test);
647 
649  addMWLSStressOperators(boost::shared_ptr<CrackFrontElement> &fe_rhs,
650  boost::shared_ptr<CrackFrontElement> &fe_lhs);
651 
653  addMWLSDensityOperators(boost::shared_ptr<CrackFrontElement> &fe_rhs,
654  boost::shared_ptr<CrackFrontElement> &fe_lhs);
655 
656  /// \brief Update fixed nodes
657  MoFEMErrorCode updateMaterialFixedNode(const bool fix_front,
658  const bool fix_small_g,
659  const bool debug = false);
660 
661  /**@}*/
662 
663  /** \name Assemble and calculate */
664 
665  /**@{*/
666 
667  /// \brief assemble elastic part of matrix, by running elastic finite element
668  /// instance
669  MoFEMErrorCode calculateElasticEnergy(DM dm, const std::string msg = "");
670 
671  /** \brief assemble material forces, by running material finite element
672  * instance
673  * @param dm dm used to calculate material forces, sub-dm
674  * of dm_crack_propagation
675  * @param q DoF sub-vector for crack propagation DM
676  * @param f Rhs sub-vector for crack propagation DM
677  * @param verb compilation parameter determining the amount
678  * of information printed
679  * @param debug flag for debugging
680  * @return error code
681  */
682  MoFEMErrorCode calculateMaterialForcesDM(DM dm, Vec q, Vec f,
683  const int verb = QUIET,
684  const bool debug = false);
685 
686  /** \brief assemble smoothing forces, by running material finite element
687  * instance
688  * @param dm crack propagation dm that is composition of
689  * elastic DM and material DM
690  * @param q DoF sub-vector for crack propagation DM
691  * @param f Rhs sub-vector for crack propagation DM
692  * @param verb compilation parameter determining the amount
693  * of information printed
694  * @param debug flag for debugging
695  * @return error code
696  */
697  MoFEMErrorCode calculateSmoothingForcesDM(DM dm, Vec q, Vec f,
698  const int verb = QUIET,
699  const bool debug = false);
700 
702  const bool debug = true);
703 
704  /// \brief assemble projection matrices
705  MoFEMErrorCode calculateSurfaceProjectionMatrix(DM dm_front, DM dm_project,
706  const std::vector<int> &ids,
707  const int verb = QUIET,
708  const bool debug = false);
709 
710  /// \brief assemble crack front projection matrix (that constrains crack area
711  /// growth)
712  MoFEMErrorCode calculateFrontProjectionMatrix(DM dm_surface, DM dm_project,
713  const int verb = QUIET,
714  const bool debug = false);
715 
716  /** \brief project material forces along the crack elongation direction
717  * @param dm_project material forces DM, sub-dm of
718  * dm_crack_propagation
719  * @param f_proj Rhs sub-vector of projected material forces
720  * for material forces DM
721  * @param f Rhs sub-vector for material forces DM
722  * @param verb compilation parameter determining the amount
723  * of information printed
724  * @param debug flag for debugging
725  * @return error code
726  */
727  MoFEMErrorCode projectMaterialForcesDM(DM dm_project, Vec f, Vec f_proj,
728  const int verb = QUIET,
729  const bool debug = false);
730 
731  /// \brief calculate Griffith (driving) force
732  MoFEMErrorCode calculateGriffithForce(DM dm, const double gc, Vec f_griffith,
733  const int verb = QUIET,
734  const bool debug = false);
735 
736  /// \brief project Griffith forces
737  MoFEMErrorCode projectGriffithForce(DM dm, Vec f_griffith,
738  Vec f_griffith_proj,
739  const int verb = QUIET,
740  const bool debug = false);
741 
742  /// \brief calculate release energy
743  MoFEMErrorCode calculateReleaseEnergy(DM dm, Vec f_material_proj,
744  Vec f_griffith_proj, Vec f_lambda,
745  const double gc, const int verb = QUIET,
746  const bool debug = true);
747 
748  /** \name Solvers */
749 
750  /**@{*/
751 
752  /** \brief solve elastic problem
753  * @param dm elastic DM
754  * @param snes snes solver for elastic DM
755  * @param m stiffness sub-matrix for elastic DM
756  * @param q DoF sub-vector for elastic DM
757  * @param f Rhs sub-vector for material DM
758  * @param snes_set_up boolean flag wheather to determin snes solver
759  * @param shell_m stiffness sub-matrix for elastic DM
760  * @return error code
761  */
762  MoFEMErrorCode solveElasticDM(DM dm, SNES snes, Mat m, Vec q, Vec f,
763  bool snes_set_up, Mat *shell_m);
764 
765  /** \brief solve crack propagation problem
766  * @param dm crack propagation DM that is composed of
767  * elastic DM and material DM
768  * @param dm_elastic elastic DM
769  * @param snes snes solver for crack propagation DM
770  * @param m stiffness sub-matrix for crack propagation DM
771  * @param q DoF sub-vector for crack propagation DM
772  * @param f Rhs sub-vector for crack propagation DM
773  * @param snes_set_up boolean flag wheather to determin snes solver
774  * @return error code
775  */
776  MoFEMErrorCode solvePropagationDM(DM dm, DM dm_elastic, SNES snes, Mat m,
777  Vec q, Vec f, bool snes_set_up,
778  Mat *shell_m, Mat *m_elastic);
779 
780  /**@}*/
781 
782  /**@}*/
783 
784  /** \name Postprocess results */
785 
786  /**@{*/
787 
789  const int verb = QUIET,
790  const bool debug = false);
791 
792  MoFEMErrorCode writeCrackFont(const BitRefLevel &bit, const int step);
793 
794  /// \brief post-process results for elastic solution
795  MoFEMErrorCode postProcessDM(DM dm, const int step, const std::string fe_name,
796  const bool approx_internal_stress);
797 
798  /**@}*/
799 
800  /** \name Set data from mesh to vectors and back */
801 
802  /**@{*/
803 
804  /// \brief set field from node positions
805  MoFEMErrorCode setFieldFromCoords(const std::string field_name);
806 
807  /// \brief set material field from nodes
809 
810  /// \brief set spatial field from nodes
812 
813  /// \brief set coordinates from field
815  setCoordsFromField(const std::string field_name = "MESH_NODE_POSITIONS");
816 
817  /// \brief set singular dofs (on edges adjacent to crack front) from geometry
818  MoFEMErrorCode setSingularDofs(const string field_name,
819  const int verb = QUIET);
820 
821  /// \brief set singular dofs of material field (on edges adjacent to crack
822  /// front) from geometry
824 
825  /// \brief remove singularity
827 
828  /// \brief set maetrial field order to one
830 
831  /**@}*/
832 
833  /** \name Auxiliary method */
834 
835  /**@{*/
836 
837  boost::shared_ptr<DofEntity> arcLengthDof;
838 
839  /// \brief set pointer to arc-length DOF
841 
842  inline boost::shared_ptr<ArcLengthCtx> &getArcCtx() { return arcCtx; }
843 
844  inline boost::shared_ptr<FEMethod> getFrontArcLengthControl() {
845  boost::shared_ptr<FEMethod> arc_method(
847  ARC_LENGTH_TAG, griffithForceElement->blockData[0],
848  griffithForceElement->commonData, "LAMBDA_CRACKFRONT_AREA",
849  arcCtx));
850  return arc_method;
851  }
852 
853  inline boost::shared_ptr<AnalyticalDirichletBC::DirichletBC>
855  return analyticalDirichletBc;
856  }
857 
858  boost::shared_ptr<MWLSApprox> &getMWLSApprox() { return mwlsApprox; }
859 
860  /**@}*/
861 
862  /**
863  * \brief operator to post-process results on crack front nodes
864  */
866 
868  Vec F;
869  std::string tagName;
870  PostProcVertexMethod(MoFEM::Interface &m_field, Vec f, std::string tag_name)
871  : mField(m_field), F(f), tagName(tag_name) {}
872 
873  Tag tH;
874  double *aRray;
876 
878 
880  };
881 
882  /**
883  * \brief Constrains material displacement on both sides
884  *
885  * Both sides of crack surface have the same material displacements.
886  *
887  */
888  struct BothSurfaceConstrains : public FEMethod {
889 
891  double aLpha;
892  const std::string lambdaFieldName;
893  Range masterNodes;
895  MoFEM::Interface &m_field,
896  const std::string lambda_field_name = "LAMBDA_BOTH_SIDES")
897  : mField(m_field), lambdaFieldName(lambda_field_name), aLpha(1) {
898  ierr = getOptions();
899  CHKERRABORT(PETSC_COMM_WORLD, ierr);
900  }
901 
904  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "",
905  "Get Booth surface constrains element scaling",
906  "none");
907  CHKERRQ(ierr);
908  CHKERR PetscOptionsScalar("-booth_side_alpha", "scaling parameter", "",
909  aLpha, &aLpha, PETSC_NULL);
910  ierr = PetscOptionsEnd();
911  CHKERRQ(ierr);
913  }
914 
916 
920  };
921 
922  /// \brief Determine face orientation
926  FaceOrientation(bool crack_front)
927  : useProjectionFromCrackFront(crack_front) {}
928  virtual ~FaceOrientation() {}
930  const EntityHandle face,
931  const BitRefLevel bit);
933  const FEMethod *fe_method_ptr) {
935  EntityHandle face = fe_method_ptr->numeredEntFiniteElementPtr->getEnt();
936  BitRefLevel bit = fe_method_ptr->problemPtr->getBitRefLevel();
937  CHKERR getElementOrientation(m_field, face, bit);
939  }
940  };
941 
942  struct ArcLengthSnesCtx : public SnesCtx {
943  boost::shared_ptr<ArcLengthCtx> arcPtr;
946  SmartPetscObj<Vec> diagM;
947  ArcLengthSnesCtx(MoFEM::Interface &m_field, const std::string &problem_name,
948  boost::shared_ptr<ArcLengthCtx> &arc_ptr,
949  DMMGViaApproxOrdersCtx *dm_ctx,
950  DMMGViaApproxOrdersCtx *dm_ctx_elastic = NULL)
951  : SnesCtx(m_field, problem_name), arcPtr(arc_ptr), dmCtx(dm_ctx),
952  dmCtxElastic(dm_ctx_elastic) {}
953  };
954 
955  // private:
956 
957  /** \name Pointers to finite element instances */
958 
959  /**@{*/
960 
961  bool onlyHooke; ///< True if only Hooke material is applied
962  PetscBool onlyHookeFromOptions; ///< True if only Hooke material is applied
963 
964  /// Integrate spatial stresses, matrix and vector (not only element but all
965  /// data structure)
966  boost::shared_ptr<NonlinearElasticElement> elasticFe;
967  boost::shared_ptr<NonlinearElasticElement> materialFe;
968  boost::shared_ptr<CrackFrontElement> feLhs; ///< Integrate elastic FE
969  boost::shared_ptr<CrackFrontElement> feRhs; ///< Integrate elastic FE
970  boost::shared_ptr<CrackFrontElement>
971  feMaterialRhs; ///< Integrate material stresses, assemble vector
972  boost::shared_ptr<CrackFrontElement>
973  feMaterialLhs; ///< Integrate material stresses, assemble matrix
974  boost::shared_ptr<CrackFrontElement> feEnergy; ///< Integrate energy
975  boost::shared_ptr<ConstrainMatrixCtx>
976  projSurfaceCtx; ///< Data structure to project on the body surface
977  boost::shared_ptr<ConstrainMatrixCtx>
978  projFrontCtx; ///< Data structure to project on crack front
979  boost::shared_ptr<MWLSApprox>
980  mwlsApprox; ///< Data struture for MWLS approximation
981  /// Surface elment to calculate tractions in material space
982  boost::shared_ptr<NeumannForcesSurface::MyTriangleFE>
984  boost::shared_ptr<DirichletSpatialPositionsBc>
985  spatialDirichletBc; ///< apply Dirichlet BC to sparial positions
986  boost::shared_ptr<AnalyticalDirichletBC::DirichletBC>
987  analyticalDirichletBc; ///< apply analytical Dirichlet BC to sparial
988  ///< positions
989  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
990  surfaceForces; ///< assemble surface forces
991  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
992  surfacePressure; ///< assemble surface pressure
993  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
994  surfacePressureAle; ///< assemble surface pressure (ALE)
995  boost::shared_ptr<NeumannForcesSurface::DataAtIntegrationPts>
996  commonDataSurfacePressureAle; ///< common data at integration points (ALE)
997  boost::shared_ptr<boost::ptr_map<string, EdgeForce>>
998  edgeForces; ///< assemble edge forces
999  boost::shared_ptr<boost::ptr_map<string, NodalForce>>
1000  nodalForces; ///< assemble nodal forces
1001  boost::shared_ptr<FEMethod> assembleFlambda; ///< assemble F_lambda vector
1002  boost::shared_ptr<FEMethod> zeroFlambda; ///< assemble F_lambda vector
1003 
1004  boost::shared_ptr<CrackFrontElement>
1005  feCouplingElasticLhs; ///< FE instance to assemble coupling terms
1006  boost::shared_ptr<CrackFrontElement>
1007  feCouplingMaterialLhs; ///< FE instance to assemble coupling terms
1008 
1009  boost::shared_ptr<Smoother> smootherFe;
1010  boost::shared_ptr<Smoother::MyVolumeFE>
1011  feSmootherRhs; ///< Integrate smoothing operators
1012  boost::shared_ptr<Smoother::MyVolumeFE>
1013  feSmootherLhs; ///< Integrate smoothing operators
1014  boost::shared_ptr<VolumeLengthQuality<double>> volumeLengthDouble;
1015  boost::shared_ptr<VolumeLengthQuality<adouble>> volumeLengthAdouble;
1016 
1017  boost::shared_ptr<
1019  tangentConstrains; ///< Constrains crack front in tangent directiona
1020  boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
1022  boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
1024  map<int, boost::shared_ptr<SurfaceSlidingConstrains>> surfaceConstrain;
1025  map<int, boost::shared_ptr<EdgeSlidingConstrains>> edgeConstrains;
1026 
1027  boost::shared_ptr<BothSurfaceConstrains> bothSidesConstrains;
1028  boost::shared_ptr<BothSurfaceConstrains> bothSidesContactConstrains;
1029  boost::shared_ptr<DirichletFixFieldAtEntitiesBc> fixMaterialEnts;
1030 
1031  boost::shared_ptr<GriffithForceElement> griffithForceElement;
1032  boost::shared_ptr<GriffithForceElement::MyTriangleFE> feGriffithForceRhs;
1033  boost::shared_ptr<GriffithForceElement::MyTriangleFE> feGriffithForceLhs;
1034  boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrainsDelta>
1036  boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrains>
1038  boost::shared_ptr<GriffithForceElement::MyTriangleFEConstrains>
1040 
1041  boost::shared_ptr<boost::ptr_map<string, NeumannForcesSurface>>
1042  analiticalSurfaceElement; ///< Finite element to integrate tractions given
1043  ///< by analytical formula
1044  boost::shared_ptr<FaceElementForcesAndSourcesCore> feSpringLhsPtr;
1045  boost::shared_ptr<FaceElementForcesAndSourcesCore> feSpringRhsPtr;
1046  /**@}*/
1047 
1048  /** \name Entity ranges */
1049 
1050  /**@{*/
1051 
1052  Range bitEnts;
1054  Range bodySkin;
1055  Range crackFaces;
1058  Range crackFront;
1063 
1064  EntityHandle
1065  arcLengthVertex; ///< Vertex associated with dof for arc-length control
1066 
1067  /**@}*/
1068 
1069  /** \name Auxiliary data */
1070 
1071  /**@{*/
1072 
1074 
1077  double loadScale;
1078 
1079  double maxG1;
1080  double maxJ;
1082  map<EntityHandle, double> mapG1; ///< hashmap of g1 - release energy at nodes
1083  map<EntityHandle, double> mapG3; ///< hashmap of g3 - release energy at nodes
1084  map<EntityHandle, double> mapJ; ///< hashmap of J - release energy at nodes
1085  map<EntityHandle, VectorDouble3>
1086  mapMatForce; ///< hashmap of material force at nodes
1087  map<EntityHandle, VectorDouble3>
1088  mapGriffith; ///< hashmap of Griffith energy at nodes
1089 
1090  map<EntityHandle, double> mapSmoothingForceFactor;
1091 
1092  bool doSurfaceProjection; ///< If crack surface is immersed in the body do not
1093  ///< projection
1094 
1095  double smootherAlpha; ///< Controls mesh smoothing
1096  double initialSmootherAlpha; ///< Controls mesh smoothing, set at start of
1097  ///< analysis
1098  double smootherGamma; ///< Controls mesh smoothing
1099 
1100  /**@}*/
1101 
1102  /** \name arc-length options */
1103 
1104  /**@{*/
1105 
1106  double arcAlpha;
1107  double arcBeta;
1108  double arcS;
1109 
1110  boost::shared_ptr<ArcLengthMatShell> arcMatCtx;
1111  boost::shared_ptr<ArcLengthCtx> arcCtx;
1112  boost::shared_ptr<PCArcLengthCtx> pcArcLengthCtx;
1113 
1114  /**@}*/
1115 
1116  /** \name Contact parameters, entity ranges and data structures */
1117 
1118  /**@{*/
1119 
1120  int contactOrder; ///< Approximation order for spatial positions
1121  ///< used in contact elements
1122  int contactLambdaOrder; ///< Approximation order for field of Lagrange
1123  ///< multipliers used to enforce contact
1124  double rValue; ///< Parameter for regularizing absolute value
1125  ///< function in the complementarity function
1126  double cnValue; ///< Augmentation parameter in the complementarity function
1127  ///< used to enforce KKT constraints
1128 
1129  PetscBool ignoreContact; ///< If set true, contact interfaces are ignored
1130  ///< and prism interfaces are not created
1131  PetscBool fixContactNodes; ///< If set true, contact nodes are fixed in
1132  ///< the material configuration
1133  PetscBool printContactState;
1134 
1136 
1144 
1145  boost::shared_ptr<SimpleContactProblem> contactProblem;
1146 
1147  boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1149  boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1151  boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1153  boost::shared_ptr<PostProcFaceOnRefinedMesh> fePostProcFaceSimpleContact;
1154  boost::shared_ptr<SimpleContactProblem::CommonDataSimpleContact>
1156 
1157  boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1159  boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1161  boost::shared_ptr<SimpleContactProblem::SimpleContactElement>
1163  boost::shared_ptr<SimpleContactProblem::CommonDataSimpleContact>
1165 
1166  /**@}*/
1167 
1168  /** \name Data for bone */
1169 
1170  /**@{*/
1171 
1172  double rHo0; ///< Reference density if bone is analyzed
1173  double nBone; ///< Exponent parameter in bone density
1174 
1175  /**@}*/
1176 
1177  /** \name Interfaces */
1178 
1179  /**@{*/
1180 
1181  const boost::scoped_ptr<UnknownInterface> cpSolversPtr;
1182  const boost::scoped_ptr<UnknownInterface> cpMeshCutPtr;
1183 
1184  /**@}*/
1185 
1186  /** \name Cut/Refine/Split options */
1187 
1188  /**@{*/
1189 
1190  PetscBool doCutMesh;
1192 
1193  /**@}*/
1194 };
1195 
1196 } // namespace FractureMechanics
1197 
1198 #endif // __CRACK_PROPAGATION_HPP__
boost::shared_ptr< SimpleContactProblem::CommonDataSimpleContact > commonDataSimpleContact
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 addMWLSDensityOperators(boost::shared_ptr< CrackFrontElement > &fe_rhs, boost::shared_ptr< CrackFrontElement > &fe_lhs)
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.
map< EntityHandle, VectorDouble3 > mapGriffith
hashmap of Griffith energy at nodes
boost::shared_ptr< FEMethod > getFrontArcLengthControl()
double smootherGamma
Controls mesh smoothing.
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.
static FTensor::Tensor2_symmetric< double, 3 > analyticalStrainFunction(FTensor::Tensor1< FTensor::PackPtr< double *, 1 >, 3 > &t_coords)
boost::shared_ptr< ConstrainMatrixCtx > projFrontCtx
Data structure to project on crack front.
MoFEMErrorCode deleteEntities(const int verb=QUIET, const bool debug=false)
Vec F
PetscBool isGcHeterogeneous
flag for heterogeneous gc
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
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feRhsSimpleContact
MoFEMErrorCode setFieldFromCoords(const std::string field_name)
set field from node positions
boost::shared_ptr< NeumannForcesSurface::MyTriangleFE > feMaterialAnaliticalTraction
Surface elment to calculate tractions in material space.
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 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 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.
PetscBool isPressureAle
If true surface pressure is considered in ALE.
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.
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfacePressureAle
assemble surface pressure (ALE)
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:509
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)
boost::shared_ptr< NeumannForcesSurface::DataAtIntegrationPts > commonDataSurfacePressureAle
common data at integration points (ALE)
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
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:485
MoFEMErrorCode calculateSmoothingForceFactor(const int verb=QUIET, const bool debug=true)
double betaGc
heterogeneous Griffith energy exponent
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 declareBothSidesFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feRhsSimpleContactALEMaterial
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:516
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 postProcess()
function is run at the end of loop
boost::shared_ptr< Smoother::MyVolumeFE > feSmootherLhs
Integrate smoothing operators.
boost::shared_ptr< PostProcFaceOnRefinedMesh > fePostProcFaceSimpleContact
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
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feLhsSimpleContact
map< EntityHandle, double > mapG3
hashmap of g3 - release energy at nodes
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feLhsSimpleContactALE
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > fePostProcSimpleContact
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
boost::shared_ptr< SimpleContactProblem::SimpleContactElement > feLhsSimpleContactALEMaterial
MoFEMErrorCode savePositionsOnCrackFrontDM(DM dm, Vec q, const int verb=QUIET, const bool debug=false)
static Index< 'm', 3 > m
PetscBool fixContactNodes
the material configuration
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
boost::shared_ptr< ObosleteUsersModules::TangentWithMeshSmoothingFrontConstrain > tangentConstrains
Constrains crack front in tangent directiona.
boost::shared_ptr< GriffithForceElement > griffithForceElement
Structure for DM for multi-grid via approximation orders.
MoFEMErrorCode createDMs(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 buildBothSidesFieldId(const BitRefLevel bit, 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.
MoFEMErrorCode readMedFile()
read mesh file
boost::shared_ptr< GriffithForceElement::MyTriangleFEConstrains > feGriffithConstrainsRhs
MoFEMErrorCode declareSmoothingFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
declare mesh smoothing finite elements
BothSurfaceConstrains(MoFEM::Interface &m_field, const std::string lambda_field_name="LAMBDA_BOTH_SIDES")
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.
std::map< std::string, BitRefLevel > mapBitLevel
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
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.
boost::shared_ptr< BothSurfaceConstrains > bothSidesContactConstrains
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)
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.
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 createProblemDataStructures(const std::vector< int > surface_ids, const int verb=QUIET, const bool debug=false)
Construct problem data structures.
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 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)
boost::shared_ptr< SimpleContactProblem::CommonDataSimpleContact > commonDataSimpleContactALE
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:604
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 addMWLSStressOperators(boost::shared_ptr< CrackFrontElement > &fe_rhs, boost::shared_ptr< CrackFrontElement > &fe_lhs)
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.
MoFEMErrorCode declareFrontFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true, const bool verb=QUIET)
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
tangent_tests
Names of tangent matrices tests.
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
double griffithR
Griffith regularisation parameter.
MoFEMErrorCode getArcLengthDof()
set pointer to arc-length DOF
boost::shared_ptr< SimpleContactProblem > contactProblem
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 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 resolveSharedBitRefLevel(const BitRefLevel bit, const int verb=QUIET, const bool debug=false)
resole shared entities by bit level
boost::shared_ptr< MWLSApprox > mwlsApprox
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:415
MoFEMErrorCode declareArcLengthFE(const BitRefLevel bits, const int verb=QUIET)
create arc-length element entity and declare elemets
MoFEMErrorCode testJacobians(const BitRefLevel bit, const BitRefLevel mask, tangent_tests test)
test LHS Jacobians
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.
MoFEMErrorCode createBcDM(DM *dm, const std::string prb_name, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set())
Create problem to calculate boundary conditions.
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
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< DirichletFixFieldAtEntitiesBc > fixMaterialEnts
MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const FEMethod *fe_method_ptr)
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 declareExternalForcesFE(const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), const bool proc_only=true)
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
boost::shared_ptr< boost::ptr_map< string, NeumannForcesSurface > > surfaceForces
assemble surface forces
MoFEMErrorCode createElasticDM(DM *dm, const std::string prb_name, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set())
Create elastic problem DM.
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 buildCrackFrontFieldId(const BitRefLevel bit, const bool build_fields=true, const int verb=QUIET, const bool debug=false)
declare crack surface files
Tag tH