www.mooseframework.org
Public Member Functions | Public Attributes | Protected Attributes | Private Member Functions | Private Attributes | List of all members
XFEM Class Reference

This is the XFEM class. More...

#include <XFEM.h>

Inheritance diagram for XFEM:
[legend]

Public Member Functions

 XFEM (const InputParameters &params)
 
 ~XFEM ()
 
void addGeometricCut (GeometricCutUserObject *geometric_cut)
 
void addStateMarkedElem (unsigned int elem_id, RealVectorValue &normal)
 
void addStateMarkedElem (unsigned int elem_id, RealVectorValue &normal, unsigned int marked_side)
 
void addStateMarkedFrag (unsigned int elem_id, RealVectorValue &normal)
 
void clearStateMarkedElems ()
 
void addGeomMarkedElem2D (const unsigned int elem_id, const Xfem::GeomMarkedElemInfo2D geom_info, const unsigned int interface_id)
 Add information about a new cut to be performed on a specific 2d element. More...
 
void addGeomMarkedElem3D (const unsigned int elem_id, const Xfem::GeomMarkedElemInfo3D geom_info, const unsigned int interface_id)
 Add information about a new cut to be performed on a specific 3d element. More...
 
void clearGeomMarkedElems ()
 Clear out the list of elements to be marked for cutting. More...
 
virtual bool update (Real time, const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux) override
 
virtual void initSolution (const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux) override
 
void buildEFAMesh ()
 
bool markCuts (Real time)
 
bool markCutEdgesByGeometry ()
 
bool markCutEdgesByState (Real time)
 
bool markCutFacesByGeometry ()
 
bool markCutFacesByState ()
 
bool initCutIntersectionEdge (Point cut_origin, RealVectorValue cut_normal, Point &edge_p1, Point &edge_p2, Real &dist)
 
bool cutMeshWithEFA (const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux)
 
bool healMesh ()
 Potentially heal the mesh by merging some of the pairs of partial elements cut by XFEM back into single elements if indicated by the cutting objects. More...
 
virtual bool updateHeal () override
 
Point getEFANodeCoords (EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=nullptr) const
 
Real getPhysicalVolumeFraction (const Elem *elem) const
 Get the volume fraction of an element that is physical. More...
 
bool isPointInsidePhysicalDomain (const Elem *elem, const Point &point) const
 Return true if the point is inside the element physical domain Note: if this element is not cut, return true too. More...
 
Real getCutPlane (const Elem *elem, const Xfem::XFEM_CUTPLANE_QUANTITY quantity, unsigned int plane_id) const
 Get specified component of normal or origin for cut plane for a given element. More...
 
bool isElemAtCrackTip (const Elem *elem) const
 
bool isElemCut (const Elem *elem, XFEMCutElem *&xfce) const
 
bool isElemCut (const Elem *elem) const
 
void getFragmentFaces (const Elem *elem, std::vector< std::vector< Point >> &frag_faces, bool displaced_mesh=false) const
 
void storeCrackTipOriginAndDirection ()
 
void correctCrackExtensionDirection (const Elem *elem, EFAElement2D *CEMElem, EFAEdge *orig_edge, Point normal, Point crack_tip_origin, Point crack_tip_direction, Real &distance_keep, unsigned int &edge_id_keep, Point &normal_keep)
 
void getCrackTipOrigin (std::map< unsigned int, const Elem *> &elem_id_crack_tip, std::vector< Point > &crack_front_points)
 
Xfem::XFEM_QRULEgetXFEMQRule ()
 
void setXFEMQRule (std::string &xfem_qrule)
 
void setCrackGrowthMethod (bool use_crack_growth_increment, Real crack_growth_increment)
 
void setDebugOutputLevel (unsigned int debug_output_level)
 Controls amount of debugging information output. More...
 
virtual bool getXFEMWeights (MooseArray< Real > &weights, const Elem *elem, QBase *qrule, const MooseArray< Point > &q_points) override
 
virtual bool getXFEMFaceWeights (MooseArray< Real > &weights, const Elem *elem, QBase *qrule, const MooseArray< Point > &q_points, unsigned int side) override
 
virtual const ElementPairLocator::ElementPairListgetXFEMCutElemPairs (unsigned int interface_id)
 Get the list of cut element pairs corresponding to a given interface ID. More...
 
virtual const ElementPairLocator::ElementPairListgetXFEMDisplacedCutElemPairs (unsigned int interface_id)
 Get the list of cut element pairs on the displaced mesh corresponding to a given interface ID. More...
 
virtual unsigned int getGeometricCutID (const GeometricCutUserObject *gcu)
 Get the interface ID corresponding to a given GeometricCutUserObject. More...
 
virtual void getXFEMIntersectionInfo (const Elem *elem, unsigned int plane_id, Point &normal, std::vector< Point > &intersectionPoints, bool displaced_mesh=false) const
 
virtual void getXFEMqRuleOnLine (std::vector< Point > &intersection_points, std::vector< Point > &quad_pts, std::vector< Real > &quad_wts) const
 
virtual void getXFEMqRuleOnSurface (std::vector< Point > &intersection_points, std::vector< Point > &quad_pts, std::vector< Real > &quad_wts) const
 
bool has_secondary_cut ()
 
EFAElement2DgetEFAElem2D (const Elem *elem)
 Get the EFAElement2D object for a specified libMesh element. More...
 
EFAElement3DgetEFAElem3D (const Elem *elem)
 Get the EFAElement3D object for a specified libMesh element. More...
 
void getFragmentEdges (const Elem *elem, EFAElement2D *CEMElem, std::vector< std::vector< Point >> &frag_edges) const
 
void getFragmentFaces (const Elem *elem, EFAElement3D *CEMElem, std::vector< std::vector< Point >> &frag_faces) const
 
const std::map< const Elem *, std::vector< Point > > & getCrackTipOriginMap () const
 
CutSubdomainID getCutSubdomainID (const GeometricCutUserObject *gcuo, const Elem *cut_elem, const Elem *parent_elem=nullptr) const
 Determine which cut subdomain the element belongs to relative to the cut. More...
 
void setMesh (MooseMesh *mesh)
 
void setDisplacedMesh (MooseMesh *displaced_mesh)
 
void setMaterialData (const std::vector< MaterialData * > &data)
 
void setBoundaryMaterialData (const std::vector< MaterialData * > &data)
 

Public Attributes

const ConsoleStream _console
 

Protected Attributes

FEProblemBase_fe_problem
 
std::vector< MaterialData *> _material_data
 
std::vector< MaterialData *> _bnd_material_data
 
MooseMesh_moose_mesh
 
MooseMesh_moose_displaced_mesh
 
MeshBase * _mesh
 
MeshBase * _displaced_mesh
 

Private Member Functions

void storeSolutionForNode (const Node *node_to_store_to, const Node *node_to_store_from, SystemBase &sys, std::map< unique_id_type, std::vector< Real >> &stored_solution, const NumericVector< Number > &current_solution, const NumericVector< Number > &old_solution, const NumericVector< Number > &older_solution)
 Store the solution in stored_solution for a given node. More...
 
void storeSolutionForElement (const Elem *elem_to_store_to, const Elem *elem_to_store_from, SystemBase &sys, std::map< unique_id_type, std::vector< Real >> &stored_solution, const NumericVector< Number > &current_solution, const NumericVector< Number > &old_solution, const NumericVector< Number > &older_solution)
 Store the solution in stored_solution for a given element. More...
 
void setSolution (SystemBase &sys, const std::map< unique_id_type, std::vector< Real >> &stored_solution, NumericVector< Number > &current_solution, NumericVector< Number > &old_solution, NumericVector< Number > &older_solution)
 Set the solution for all locally-owned nodes/elements that have stored values. More...
 
void setSolutionForDOFs (const std::vector< Real > &stored_solution, const std::vector< dof_id_type > &stored_solution_dofs, NumericVector< Number > &current_solution, NumericVector< Number > &old_solution, NumericVector< Number > &older_solution)
 Set the solution for a set of DOFs. More...
 
std::vector< dof_id_typegetElementSolutionDofs (const Elem *elem, SystemBase &sys) const
 Get a vector of the dof indices for all components of all variables associated with an element. More...
 
std::vector< dof_id_typegetNodeSolutionDofs (const Node *node, SystemBase &sys) const
 Get a vector of the dof indices for all components of all variables associated with a node. More...
 
const GeometricCutUserObjectgetGeometricCutForElem (const Elem *elem) const
 Get the GeometricCutUserObject associated with an element. More...
 
void storeMaterialPropertiesForElementHelper (const Elem *elem, MaterialPropertyStorage &storage)
 
void storeMaterialPropertiesForElement (const Elem *parent_elem, const Elem *child_elem)
 Helper function to store the material properties of a healed element. More...
 
void loadMaterialPropertiesForElementHelper (const Elem *elem, const Xfem::CachedMaterialProperties &cached_props, MaterialPropertyStorage &storage) const
 Load the material properties. More...
 
void loadMaterialPropertiesForElement (const Elem *elem, const Elem *elem_from, std::unordered_map< const Elem *, Xfem::CutElemInfo > &cached_cei) const
 Helper function to store the material properties of a healed element. More...
 
const Node * pickFirstPhysicalNode (const Elem *e, const Elem *e0) const
 Return the first node in the provided element that is found to be in the physical domain. More...
 

Private Attributes

bool _has_secondary_cut
 
Xfem::XFEM_QRULE _XFEM_qrule
 
bool _use_crack_growth_increment
 
Real _crack_growth_increment
 
std::vector< const GeometricCutUserObject * > _geometric_cuts
 
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
 
std::set< const Elem * > _crack_tip_elems
 
std::set< const Elem * > _crack_tip_elems_to_be_healed
 
std::map< unsigned int, ElementPairLocator::ElementPairList_sibling_elems
 
std::map< unsigned int, ElementPairLocator::ElementPairList_sibling_displaced_elems
 
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
 
std::map< const Elem *, RealVectorValue_state_marked_elems
 
std::set< const Elem * > _state_marked_frags
 
std::map< const Elem *, unsigned int_state_marked_elem_sides
 
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo2D > > _geom_marked_elems_2d
 Data structure for storing information about all 2D elements to be cut by geometry. More...
 
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo3D > > _geom_marked_elems_3d
 Data structure for storing information about all 3D elements to be cut by geometry. More...
 
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
 Data structure for storing the elements cut by specific geometric cutters. More...
 
std::map< const GeometricCutUserObject *, unsigned int_geom_marker_id_map
 Data structure for storing the GeommetricCutUserObjects and their corresponding id. More...
 
ElementFragmentAlgorithm _efa_mesh
 
unsigned int _debug_output_level
 Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mesh 3: Full dump of element fragment algorithm mesh. More...
 
std::map< unique_id_type, std::vector< Real > > _cached_solution
 Data structure to store the nonlinear solution for nodes/elements affected by XFEM For each node/element, this is stored as a vector that contains all components of all applicable variables in an order defined by getElementSolutionDofs() or getNodeSolutionDofs(). More...
 
std::map< unique_id_type, std::vector< Real > > _cached_aux_solution
 Data structure to store the auxiliary solution for nodes/elements affected by XFEM For each node/element, this is stored as a vector that contains all components of all applicable variables in an order defined by getElementSolutionDofs() or getNodeSolutionDofs(). More...
 
std::unordered_map< const Elem *, Xfem::CutElemInfo_geom_cut_elems
 All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK. More...
 
std::unordered_map< const Elem *, Xfem::CutElemInfo_old_geom_cut_elems
 All geometrically cut elements and their CutElemInfo before the current execution of XFEM_MARK. More...
 

Detailed Description

This is the XFEM class.

This class implements algorithms for dynamic mesh modification in support of a phantom node approach for XFEM

Definition at line 107 of file XFEM.h.

Constructor & Destructor Documentation

◆ XFEM()

XFEM::XFEM ( const InputParameters params)
explicit

Definition at line 36 of file XFEM.C.

37  : XFEMInterface(params), _efa_mesh(Moose::out), _debug_output_level(1)
38 {
39 #ifndef LIBMESH_ENABLE_UNIQUE_ID
40  mooseError("MOOSE requires unique ids to be enabled in libmesh (configure with "
41  "--enable-unique-id) to use XFEM!");
42 #endif
43  _has_secondary_cut = false;
44 }
unsigned int _debug_output_level
Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mes...
Definition: XFEM.h:362
void mooseError(Args &&... args)
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
bool _has_secondary_cut
Definition: XFEM.h:320
XFEMInterface(const InputParameters &params)

◆ ~XFEM()

XFEM::~XFEM ( )

Definition at line 46 of file XFEM.C.

47 {
48  for (std::map<unique_id_type, XFEMCutElem *>::iterator cemit = _cut_elem_map.begin();
49  cemit != _cut_elem_map.end();
50  ++cemit)
51  delete cemit->second;
52 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329

Member Function Documentation

◆ addGeometricCut()

void XFEM::addGeometricCut ( GeometricCutUserObject geometric_cut)

Definition at line 55 of file XFEM.C.

56 {
57  _geometric_cuts.push_back(geometric_cut);
58 
59  geometric_cut->setInterfaceID(_geometric_cuts.size() - 1);
60 
61  _geom_marker_id_map[geometric_cut] = _geometric_cuts.size() - 1;
62 }
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:327
std::map< const GeometricCutUserObject *, unsigned int > _geom_marker_id_map
Data structure for storing the GeommetricCutUserObjects and their corresponding id.
Definition: XFEM.h:353
void setInterfaceID(unsigned int interface_id)
Set the interface ID for this cutting object.

◆ addGeomMarkedElem2D()

void XFEM::addGeomMarkedElem2D ( const unsigned int  elem_id,
const Xfem::GeomMarkedElemInfo2D  geom_info,
const unsigned int  interface_id 
)

Add information about a new cut to be performed on a specific 2d element.

Parameters
elem_idThe id of the element to be cut
geom_infoThe object containing information about the cut to be performed
interface_idThe ID of the interface

Definition at line 157 of file XFEM.C.

160 {
161  Elem * elem = _mesh->elem_ptr(elem_id);
162  _geom_marked_elems_2d[elem].push_back(geom_info);
163  _geom_marker_id_elems[interface_id].insert(elem_id);
164 }
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo2D > > _geom_marked_elems_2d
Data structure for storing information about all 2D elements to be cut by geometry.
Definition: XFEM.h:344
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
Data structure for storing the elements cut by specific geometric cutters.
Definition: XFEM.h:350
MeshBase * _mesh

◆ addGeomMarkedElem3D()

void XFEM::addGeomMarkedElem3D ( const unsigned int  elem_id,
const Xfem::GeomMarkedElemInfo3D  geom_info,
const unsigned int  interface_id 
)

Add information about a new cut to be performed on a specific 3d element.

Parameters
elem_idThe id of the element to be cut
geom_infoThe object containing information about the cut to be performed
interface_idThe ID of the interface

Definition at line 167 of file XFEM.C.

170 {
171  Elem * elem = _mesh->elem_ptr(elem_id);
172  _geom_marked_elems_3d[elem].push_back(geom_info);
173  _geom_marker_id_elems[interface_id].insert(elem_id);
174 }
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo3D > > _geom_marked_elems_3d
Data structure for storing information about all 3D elements to be cut by geometry.
Definition: XFEM.h:347
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
Data structure for storing the elements cut by specific geometric cutters.
Definition: XFEM.h:350
MeshBase * _mesh

◆ addStateMarkedElem() [1/2]

void XFEM::addStateMarkedElem ( unsigned int  elem_id,
RealVectorValue normal 
)

Definition at line 111 of file XFEM.C.

Referenced by addStateMarkedElem(), and addStateMarkedFrag().

112 {
113  Elem * elem = _mesh->elem_ptr(elem_id);
114  std::map<const Elem *, RealVectorValue>::iterator mit;
115  mit = _state_marked_elems.find(elem);
116  if (mit != _state_marked_elems.end())
117  mooseError(" ERROR: element ", elem->id(), " already marked for crack growth.");
118  _state_marked_elems[elem] = normal;
119 }
void mooseError(Args &&... args)
MeshBase * _mesh
std::map< const Elem *, RealVectorValue > _state_marked_elems
Definition: XFEM.h:339

◆ addStateMarkedElem() [2/2]

void XFEM::addStateMarkedElem ( unsigned int  elem_id,
RealVectorValue normal,
unsigned int  marked_side 
)

Definition at line 122 of file XFEM.C.

123 {
124  addStateMarkedElem(elem_id, normal);
125  Elem * elem = _mesh->elem_ptr(elem_id);
126  std::map<const Elem *, unsigned int>::iterator mit;
127  mit = _state_marked_elem_sides.find(elem);
128  if (mit != _state_marked_elem_sides.end())
129  mooseError(" ERROR: side of element ", elem->id(), " already marked for crack initiation.");
130 
131  _state_marked_elem_sides[elem] = marked_side;
132 }
void mooseError(Args &&... args)
MeshBase * _mesh
void addStateMarkedElem(unsigned int elem_id, RealVectorValue &normal)
Definition: XFEM.C:111
std::map< const Elem *, unsigned int > _state_marked_elem_sides
Definition: XFEM.h:341

◆ addStateMarkedFrag()

void XFEM::addStateMarkedFrag ( unsigned int  elem_id,
RealVectorValue normal 
)

Definition at line 135 of file XFEM.C.

136 {
137  addStateMarkedElem(elem_id, normal);
138  Elem * elem = _mesh->elem_ptr(elem_id);
139  std::set<const Elem *>::iterator mit;
140  mit = _state_marked_frags.find(elem);
141  if (mit != _state_marked_frags.end())
142  mooseError(
143  " ERROR: element ", elem->id(), " already marked for fragment-secondary crack initiation.");
144 
145  _state_marked_frags.insert(elem);
146 }
void mooseError(Args &&... args)
std::set< const Elem * > _state_marked_frags
Definition: XFEM.h:340
MeshBase * _mesh
void addStateMarkedElem(unsigned int elem_id, RealVectorValue &normal)
Definition: XFEM.C:111

◆ buildEFAMesh()

void XFEM::buildEFAMesh ( )

Definition at line 336 of file XFEM.C.

Referenced by update(), and updateHeal().

337 {
338  _efa_mesh.reset();
339 
340  // Load all existing elements in to EFA mesh
341  for (auto & elem : _mesh->element_ptr_range())
342  {
343  std::vector<unsigned int> quad;
344  for (unsigned int i = 0; i < elem->n_nodes(); ++i)
345  quad.push_back(elem->node_id(i));
346 
347  if (_mesh->mesh_dimension() == 2)
348  _efa_mesh.add2DElement(quad, elem->id());
349  else if (_mesh->mesh_dimension() == 3)
350  _efa_mesh.add3DElement(quad, elem->id());
351  else
352  mooseError("XFEM only works for 2D and 3D");
353  }
354 
355  // Restore fragment information for elements that have been previously cut
356  for (auto & elem : _mesh->element_ptr_range())
357  {
358  std::map<unique_id_type, XFEMCutElem *>::iterator cemit = _cut_elem_map.find(elem->unique_id());
359  if (cemit != _cut_elem_map.end())
360  {
361  XFEMCutElem * xfce = cemit->second;
362  EFAElement * CEMElem = _efa_mesh.getElemByID(elem->id());
363  _efa_mesh.restoreFragmentInfo(CEMElem, xfce->getEFAElement());
364  }
365  }
366 
367  // Must update edge neighbors before restore edge intersections. Otherwise, when we
368  // add edge intersections, we do not have neighbor information to use.
369  // Correction: no need to use neighbor info now
372 }
void mooseError(Args &&... args)
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
EFAElement * add3DElement(const std::vector< unsigned int > &quad, unsigned int id)
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
virtual const EFAElement * getEFAElement() const =0
void restoreFragmentInfo(EFAElement *const elem, const EFAElement *const from_elem)
MeshBase * _mesh
EFAElement * add2DElement(const std::vector< unsigned int > &quad, unsigned int id)
EFAElement * getElemByID(unsigned int id)

◆ clearGeomMarkedElems()

void XFEM::clearGeomMarkedElems ( )

Clear out the list of elements to be marked for cutting.

Called after cutting is done.

Definition at line 177 of file XFEM.C.

Referenced by update().

178 {
179  _geom_marked_elems_2d.clear();
180  _geom_marked_elems_3d.clear();
181 }
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo3D > > _geom_marked_elems_3d
Data structure for storing information about all 3D elements to be cut by geometry.
Definition: XFEM.h:347
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo2D > > _geom_marked_elems_2d
Data structure for storing information about all 2D elements to be cut by geometry.
Definition: XFEM.h:344

◆ clearStateMarkedElems()

void XFEM::clearStateMarkedElems ( )

Definition at line 149 of file XFEM.C.

Referenced by update().

150 {
151  _state_marked_elems.clear();
152  _state_marked_frags.clear();
153  _state_marked_elem_sides.clear();
154 }
std::set< const Elem * > _state_marked_frags
Definition: XFEM.h:340
std::map< const Elem *, unsigned int > _state_marked_elem_sides
Definition: XFEM.h:341
std::map< const Elem *, RealVectorValue > _state_marked_elems
Definition: XFEM.h:339

◆ correctCrackExtensionDirection()

void XFEM::correctCrackExtensionDirection ( const Elem *  elem,
EFAElement2D CEMElem,
EFAEdge orig_edge,
Point  normal,
Point  crack_tip_origin,
Point  crack_tip_direction,
Real distance_keep,
unsigned int edge_id_keep,
Point &  normal_keep 
)

Definition at line 444 of file XFEM.C.

Referenced by markCutEdgesByState().

453 {
454  std::vector<Point> edge_ends(2, Point(0.0, 0.0, 0.0));
455  Point edge1(0.0, 0.0, 0.0);
456  Point edge2(0.0, 0.0, 0.0);
457  Point left_angle(0.0, 0.0, 0.0);
458  Point right_angle(0.0, 0.0, 0.0);
459  Point left_angle_normal(0.0, 0.0, 0.0);
460  Point right_angle_normal(0.0, 0.0, 0.0);
461  Point crack_direction_normal(0.0, 0.0, 0.0);
462  Point edge1_to_tip(0.0, 0.0, 0.0);
463  Point edge2_to_tip(0.0, 0.0, 0.0);
464  Point edge1_to_tip_normal(0.0, 0.0, 0.0);
465  Point edge2_to_tip_normal(0.0, 0.0, 0.0);
466 
467  Real cos_45 = std::cos(45.0 / 180.0 * 3.14159);
468  Real sin_45 = std::sin(45.0 / 180.0 * 3.14159);
469 
470  left_angle(0) = cos_45 * crack_tip_direction(0) - sin_45 * crack_tip_direction(1);
471  left_angle(1) = sin_45 * crack_tip_direction(0) + cos_45 * crack_tip_direction(1);
472 
473  right_angle(0) = cos_45 * crack_tip_direction(0) + sin_45 * crack_tip_direction(1);
474  right_angle(1) = -sin_45 * crack_tip_direction(0) + cos_45 * crack_tip_direction(1);
475 
476  left_angle_normal(0) = -left_angle(1);
477  left_angle_normal(1) = left_angle(0);
478 
479  right_angle_normal(0) = -right_angle(1);
480  right_angle_normal(1) = right_angle(0);
481 
482  crack_direction_normal(0) = -crack_tip_direction(1);
483  crack_direction_normal(1) = crack_tip_direction(0);
484 
485  Real angle_min = 0.0;
486  Real distance = 0.0;
487  unsigned int nsides = CEMElem->numEdges();
488 
489  for (unsigned int i = 0; i < nsides; ++i)
490  {
491  if (!orig_edge->isPartialOverlap(*CEMElem->getEdge(i)))
492  {
493  edge_ends[0] = getEFANodeCoords(CEMElem->getEdge(i)->getNode(0), CEMElem, elem);
494  edge_ends[1] = getEFANodeCoords(CEMElem->getEdge(i)->getNode(1), CEMElem, elem);
495 
496  edge1_to_tip = (edge_ends[0] * 0.95 + edge_ends[1] * 0.05) - crack_tip_origin;
497  edge2_to_tip = (edge_ends[0] * 0.05 + edge_ends[1] * 0.95) - crack_tip_origin;
498 
499  edge1_to_tip /= pow(edge1_to_tip.norm_sq(), 0.5);
500  edge2_to_tip /= pow(edge2_to_tip.norm_sq(), 0.5);
501 
502  edge1_to_tip_normal(0) = -edge1_to_tip(1);
503  edge1_to_tip_normal(1) = edge1_to_tip(0);
504 
505  edge2_to_tip_normal(0) = -edge2_to_tip(1);
506  edge2_to_tip_normal(1) = edge2_to_tip(0);
507 
508  Real angle_edge1_normal = edge1_to_tip_normal * normal;
509  Real angle_edge2_normal = edge2_to_tip_normal * normal;
510 
511  if (std::abs(angle_edge1_normal) > std::abs(angle_min) &&
512  (edge1_to_tip * crack_tip_direction) > std::cos(45.0 / 180.0 * 3.14159))
513  {
514  edge_id_keep = i;
515  distance_keep = 0.05;
516  normal_keep = edge1_to_tip_normal;
517  angle_min = angle_edge1_normal;
518  }
519  else if (std::abs(angle_edge2_normal) > std::abs(angle_min) &&
520  (edge2_to_tip * crack_tip_direction) > std::cos(45.0 / 180.0 * 3.14159))
521  {
522  edge_id_keep = i;
523  distance_keep = 0.95;
524  normal_keep = edge2_to_tip_normal;
525  angle_min = angle_edge2_normal;
526  }
527 
529  crack_tip_origin, left_angle_normal, edge_ends[0], edge_ends[1], distance) &&
530  (!CEMElem->isEdgePhantom(i)))
531  {
532  if (std::abs(left_angle_normal * normal) > std::abs(angle_min) &&
533  (edge1_to_tip * crack_tip_direction) > std::cos(45.0 / 180.0 * 3.14159))
534  {
535  edge_id_keep = i;
536  distance_keep = distance;
537  normal_keep = left_angle_normal;
538  angle_min = left_angle_normal * normal;
539  }
540  }
541  else if (initCutIntersectionEdge(
542  crack_tip_origin, right_angle_normal, edge_ends[0], edge_ends[1], distance) &&
543  (!CEMElem->isEdgePhantom(i)))
544  {
545  if (std::abs(right_angle_normal * normal) > std::abs(angle_min) &&
546  (edge2_to_tip * crack_tip_direction) > std::cos(45.0 / 180.0 * 3.14159))
547  {
548  edge_id_keep = i;
549  distance_keep = distance;
550  normal_keep = right_angle_normal;
551  angle_min = right_angle_normal * normal;
552  }
553  }
554  else if (initCutIntersectionEdge(crack_tip_origin,
555  crack_direction_normal,
556  edge_ends[0],
557  edge_ends[1],
558  distance) &&
559  (!CEMElem->isEdgePhantom(i)))
560  {
561  if (std::abs(crack_direction_normal * normal) > std::abs(angle_min) &&
562  (crack_tip_direction * crack_tip_direction) > std::cos(45.0 / 180.0 * 3.14159))
563  {
564  edge_id_keep = i;
565  distance_keep = distance;
566  normal_keep = crack_direction_normal;
567  angle_min = crack_direction_normal * normal;
568  }
569  }
570  }
571  }
572 
573  // avoid small volume fraction cut
574  if ((distance_keep - 0.05) < 0.0)
575  {
576  distance_keep = 0.05;
577  }
578  else if ((distance_keep - 0.95) > 0.0)
579  {
580  distance_keep = 0.95;
581  }
582 }
bool isPartialOverlap(const EFAEdge &other) const
Definition: EFAEdge.C:59
EFAEdge * getEdge(unsigned int edge_id) const
unsigned int numEdges() const
bool isEdgePhantom(unsigned int edge_id) const
Real distance(const Point &p)
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
bool initCutIntersectionEdge(Point cut_origin, RealVectorValue cut_normal, Point &edge_p1, Point &edge_p2, Real &dist)
Definition: XFEM.C:898
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EFANode * getNode(unsigned int index) const
Definition: EFAEdge.C:181
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=nullptr) const
Definition: XFEM.C:1596

◆ cutMeshWithEFA()

bool XFEM::cutMeshWithEFA ( const std::vector< std::shared_ptr< NonlinearSystemBase >> &  nl,
AuxiliarySystem aux 
)

Definition at line 1101 of file XFEM.C.

Referenced by update().

1103 {
1104  if (nls.size() != 1)
1105  mooseError("XFEM does not currently support multiple nonlinear systems");
1106 
1107  std::map<unsigned int, Node *> efa_id_to_new_node;
1108  std::map<unsigned int, Node *> efa_id_to_new_node2;
1109  std::map<unsigned int, Elem *> efa_id_to_new_elem;
1110  _cached_solution.clear();
1111  _cached_aux_solution.clear();
1112 
1113  // Copy the current geometric cut element info (from last time) into the
1114  // _old_geom_cut_elems.
1116  _geom_cut_elems.clear();
1117 
1119 
1120  if (_debug_output_level > 2)
1121  {
1122  _console << "\nXFEM Element fragment algorithm mesh prior to cutting:\n";
1123  _console << std::flush;
1124  _efa_mesh.printMesh();
1125  }
1126 
1128 
1129  if (_debug_output_level > 2)
1130  {
1131  _console << "\nXFEM Element fragment algorithm mesh after cutting:\n";
1132  _console << std::flush;
1133  _efa_mesh.printMesh();
1134  }
1135 
1136  const std::vector<EFANode *> new_nodes = _efa_mesh.getNewNodes();
1137  const std::vector<EFAElement *> new_elements = _efa_mesh.getChildElements();
1138  const std::vector<EFAElement *> delete_elements = _efa_mesh.getParentElements();
1139 
1140  bool mesh_changed = (new_nodes.size() + new_elements.size() + delete_elements.size() > 0);
1141 
1142  // Prepare to cache solution on DOFs modified by XFEM
1143  if (mesh_changed)
1144  {
1145  nls[0]->serializeSolution();
1146  aux.serializeSolution();
1147  if (_debug_output_level > 1)
1148  _console << "\n";
1149  }
1150  NumericVector<Number> & current_solution = *nls[0]->system().current_local_solution;
1151  NumericVector<Number> & old_solution = nls[0]->solutionOld();
1152  NumericVector<Number> & older_solution = nls[0]->solutionOlder();
1153  NumericVector<Number> & current_aux_solution = *aux.system().current_local_solution;
1154  NumericVector<Number> & old_aux_solution = aux.solutionOld();
1155  NumericVector<Number> & older_aux_solution = aux.solutionOlder();
1156 
1157  std::map<Node *, Node *> new_nodes_to_parents;
1158 
1159  // Add new nodes
1160  for (unsigned int i = 0; i < new_nodes.size(); ++i)
1161  {
1162  unsigned int new_node_id = new_nodes[i]->id();
1163  unsigned int parent_id = new_nodes[i]->parent()->id();
1164 
1165  Node * parent_node = _mesh->node_ptr(parent_id);
1166  Node * new_node = Node::build(*parent_node, _mesh->n_nodes()).release();
1167  _mesh->add_node(new_node);
1168 
1169  new_nodes_to_parents[new_node] = parent_node;
1170 
1171  new_node->set_n_systems(parent_node->n_systems());
1172  efa_id_to_new_node.insert(std::make_pair(new_node_id, new_node));
1173  if (_debug_output_level > 1)
1174  _console << "XFEM added new node: " << new_node->id() << std::endl;
1175  if (_displaced_mesh)
1176  {
1177  const Node * parent_node2 = _displaced_mesh->node_ptr(parent_id);
1178  Node * new_node2 = Node::build(*parent_node2, _displaced_mesh->n_nodes()).release();
1179  _displaced_mesh->add_node(new_node2);
1180 
1181  new_node2->set_n_systems(parent_node2->n_systems());
1182  efa_id_to_new_node2.insert(std::make_pair(new_node_id, new_node2));
1183  }
1184  }
1185 
1186  // Add new elements
1187  std::map<unsigned int, std::vector<const Elem *>> temporary_parent_children_map;
1188 
1189  std::vector<boundary_id_type> parent_boundary_ids;
1190 
1191  for (unsigned int i = 0; i < new_elements.size(); ++i)
1192  {
1193  unsigned int parent_id = new_elements[i]->getParent()->id();
1194  unsigned int efa_child_id = new_elements[i]->id();
1195 
1196  Elem * parent_elem = _mesh->elem_ptr(parent_id);
1197  Elem * libmesh_elem = Elem::build(parent_elem->type()).release();
1198 
1199  for (unsigned int m = 0; m < _geometric_cuts.size(); ++m)
1200  {
1201  for (auto & it : _sibling_elems[_geometric_cuts[m]->getInterfaceID()])
1202  {
1203  if (parent_elem == it.first)
1204  it.first = libmesh_elem;
1205  else if (parent_elem == it.second)
1206  it.second = libmesh_elem;
1207  }
1208  }
1209 
1210  // parent has at least two children
1211  if (new_elements[i]->getParent()->numChildren() > 1)
1212  temporary_parent_children_map[parent_elem->id()].push_back(libmesh_elem);
1213 
1214  Elem * parent_elem2 = nullptr;
1215  Elem * libmesh_elem2 = nullptr;
1216  if (_displaced_mesh)
1217  {
1218  parent_elem2 = _displaced_mesh->elem_ptr(parent_id);
1219  libmesh_elem2 = Elem::build(parent_elem2->type()).release();
1220 
1221  for (unsigned int m = 0; m < _geometric_cuts.size(); ++m)
1222  {
1223  for (auto & it : _sibling_displaced_elems[_geometric_cuts[m]->getInterfaceID()])
1224  {
1225  if (parent_elem2 == it.first)
1226  it.first = libmesh_elem2;
1227  else if (parent_elem2 == it.second)
1228  it.second = libmesh_elem2;
1229  }
1230  }
1231  }
1232 
1233  for (unsigned int j = 0; j < new_elements[i]->numNodes(); ++j)
1234  {
1235  unsigned int node_id = new_elements[i]->getNode(j)->id();
1236  Node * libmesh_node;
1237 
1238  std::map<unsigned int, Node *>::iterator nit = efa_id_to_new_node.find(node_id);
1239  if (nit != efa_id_to_new_node.end())
1240  libmesh_node = nit->second;
1241  else
1242  libmesh_node = _mesh->node_ptr(node_id);
1243 
1244  if (libmesh_node->processor_id() == DofObject::invalid_processor_id)
1245  libmesh_node->processor_id() = parent_elem->processor_id();
1246 
1247  libmesh_elem->set_node(j) = libmesh_node;
1248 
1249  // Store solution for all nodes affected by XFEM (even existing nodes)
1250  if (parent_elem->is_semilocal(_mesh->processor_id()))
1251  {
1252  Node * solution_node = libmesh_node; // Node from which to store solution
1253  if (new_nodes_to_parents.find(libmesh_node) != new_nodes_to_parents.end())
1254  solution_node = new_nodes_to_parents[libmesh_node];
1255 
1256  if ((_moose_mesh->isSemiLocal(solution_node)) ||
1257  (libmesh_node->processor_id() == _mesh->processor_id()))
1258  {
1259  storeSolutionForNode(libmesh_node,
1260  solution_node,
1261  *nls[0],
1263  current_solution,
1264  old_solution,
1265  older_solution);
1266  storeSolutionForNode(libmesh_node,
1267  solution_node,
1268  aux,
1270  current_aux_solution,
1271  old_aux_solution,
1272  older_aux_solution);
1273  }
1274  }
1275 
1276  Node * parent_node = parent_elem->node_ptr(j);
1277  _mesh->get_boundary_info().boundary_ids(parent_node, parent_boundary_ids);
1278  _mesh->get_boundary_info().add_node(libmesh_node, parent_boundary_ids);
1279 
1280  if (_displaced_mesh)
1281  {
1282  std::map<unsigned int, Node *>::iterator nit2 = efa_id_to_new_node2.find(node_id);
1283  if (nit2 != efa_id_to_new_node2.end())
1284  libmesh_node = nit2->second;
1285  else
1286  libmesh_node = _displaced_mesh->node_ptr(node_id);
1287 
1288  if (libmesh_node->processor_id() == DofObject::invalid_processor_id)
1289  libmesh_node->processor_id() = parent_elem2->processor_id();
1290 
1291  libmesh_elem2->set_node(j) = libmesh_node;
1292 
1293  parent_node = parent_elem2->node_ptr(j);
1294  _displaced_mesh->get_boundary_info().boundary_ids(parent_node, parent_boundary_ids);
1295  _displaced_mesh->get_boundary_info().add_node(libmesh_node, parent_boundary_ids);
1296  }
1297  }
1298 
1299  libmesh_elem->set_p_level(parent_elem->p_level());
1300  libmesh_elem->set_p_refinement_flag(parent_elem->p_refinement_flag());
1301  _mesh->add_elem(libmesh_elem);
1302  libmesh_elem->set_n_systems(parent_elem->n_systems());
1303  libmesh_elem->subdomain_id() = parent_elem->subdomain_id();
1304  libmesh_elem->processor_id() = parent_elem->processor_id();
1305 
1306  // The crack tip origin map is stored before cut, thus the elem should be updated with new
1307  // element.
1308  std::map<const Elem *, std::vector<Point>>::iterator mit =
1309  _elem_crack_origin_direction_map.find(parent_elem);
1310  if (mit != _elem_crack_origin_direction_map.end())
1311  {
1312  std::vector<Point> crack_data = _elem_crack_origin_direction_map[parent_elem];
1314  _elem_crack_origin_direction_map[libmesh_elem] = crack_data;
1315  }
1316 
1317  if (_debug_output_level > 1)
1318  _console << "XFEM added new element: " << libmesh_elem->id() << std::endl;
1319 
1320  XFEMCutElem * xfce = nullptr;
1321  if (_mesh->mesh_dimension() == 2)
1322  {
1323  EFAElement2D * new_efa_elem2d = dynamic_cast<EFAElement2D *>(new_elements[i]);
1324  if (!new_efa_elem2d)
1325  mooseError("EFAelem is not of EFAelement2D type");
1326  xfce = new XFEMCutElem2D(libmesh_elem,
1327  new_efa_elem2d,
1328  _fe_problem->assembly(0, /*nl_sys_num=*/0).qRule()->n_points(),
1329  libmesh_elem->n_sides());
1330  }
1331  else if (_mesh->mesh_dimension() == 3)
1332  {
1333  EFAElement3D * new_efa_elem3d = dynamic_cast<EFAElement3D *>(new_elements[i]);
1334  if (!new_efa_elem3d)
1335  mooseError("EFAelem is not of EFAelement3D type");
1336  xfce = new XFEMCutElem3D(libmesh_elem,
1337  new_efa_elem3d,
1338  _fe_problem->assembly(0, /*nl_sys_num=*/0).qRule()->n_points(),
1339  libmesh_elem->n_sides());
1340  }
1341  _cut_elem_map.insert(std::pair<unique_id_type, XFEMCutElem *>(libmesh_elem->unique_id(), xfce));
1342  efa_id_to_new_elem.insert(std::make_pair(efa_child_id, libmesh_elem));
1343 
1344  if (_displaced_mesh)
1345  {
1346  libmesh_elem2->set_p_level(parent_elem2->p_level());
1347  libmesh_elem2->set_p_refinement_flag(parent_elem2->p_refinement_flag());
1348  _displaced_mesh->add_elem(libmesh_elem2);
1349  libmesh_elem2->set_n_systems(parent_elem2->n_systems());
1350  libmesh_elem2->subdomain_id() = parent_elem2->subdomain_id();
1351  libmesh_elem2->processor_id() = parent_elem2->processor_id();
1352  }
1353 
1354  unsigned int n_sides = parent_elem->n_sides();
1355  for (unsigned int side = 0; side < n_sides; ++side)
1356  {
1357  _mesh->get_boundary_info().boundary_ids(parent_elem, side, parent_boundary_ids);
1358  _mesh->get_boundary_info().add_side(libmesh_elem, side, parent_boundary_ids);
1359  }
1360  if (_displaced_mesh)
1361  {
1362  n_sides = parent_elem2->n_sides();
1363  for (unsigned int side = 0; side < n_sides; ++side)
1364  {
1365  _displaced_mesh->get_boundary_info().boundary_ids(parent_elem2, side, parent_boundary_ids);
1366  _displaced_mesh->get_boundary_info().add_side(libmesh_elem2, side, parent_boundary_ids);
1367  }
1368  }
1369 
1370  unsigned int n_edges = parent_elem->n_edges();
1371  for (unsigned int edge = 0; edge < n_edges; ++edge)
1372  {
1373  _mesh->get_boundary_info().edge_boundary_ids(parent_elem, edge, parent_boundary_ids);
1374  _mesh->get_boundary_info().add_edge(libmesh_elem, edge, parent_boundary_ids);
1375  }
1376  if (_displaced_mesh)
1377  {
1378  n_edges = parent_elem2->n_edges();
1379  for (unsigned int edge = 0; edge < n_edges; ++edge)
1380  {
1381  _displaced_mesh->get_boundary_info().edge_boundary_ids(
1382  parent_elem2, edge, parent_boundary_ids);
1383  _displaced_mesh->get_boundary_info().add_edge(libmesh_elem2, edge, parent_boundary_ids);
1384  }
1385  }
1386 
1387  // TODO: Also need to copy neighbor material data
1388  if (parent_elem->processor_id() == _mesh->processor_id())
1389  {
1390  if (_material_data[0]->getMaterialPropertyStorage().hasStatefulProperties())
1391  _material_data[0]->copy(*libmesh_elem, *parent_elem, 0);
1392 
1393  if (_bnd_material_data[0]->getMaterialPropertyStorage().hasStatefulProperties())
1394  for (unsigned int side = 0; side < parent_elem->n_sides(); ++side)
1395  {
1396  _mesh->get_boundary_info().boundary_ids(parent_elem, side, parent_boundary_ids);
1397  std::vector<boundary_id_type>::iterator it_bd = parent_boundary_ids.begin();
1398  for (; it_bd != parent_boundary_ids.end(); ++it_bd)
1399  {
1400  if (_fe_problem->needBoundaryMaterialOnSide(*it_bd, 0))
1401  _bnd_material_data[0]->copy(*libmesh_elem, *parent_elem, side);
1402  }
1403  }
1404 
1405  // Store the current information about the geometrically cut element, and load cached material
1406  // properties into the new child element, if any.
1407  const GeometricCutUserObject * gcuo = getGeometricCutForElem(parent_elem);
1408  if (gcuo && gcuo->shouldHealMesh())
1409  {
1410  CutSubdomainID gcsid = getCutSubdomainID(gcuo, libmesh_elem, parent_elem);
1411  Xfem::CutElemInfo cei(parent_elem, gcuo, gcsid);
1412  _geom_cut_elems.emplace(libmesh_elem, cei);
1413  // Find the element to copy data from.
1414  // Iterate through the old geometrically cut elements, if its parent element AND the
1415  // geometric cut user object AND the cut subdomain ID are the same as the
1416  // current element, then that must be it.
1417  for (auto old_cei : _old_geom_cut_elems)
1418  if (cei.match(old_cei.second))
1419  {
1420  loadMaterialPropertiesForElement(libmesh_elem, old_cei.first, _old_geom_cut_elems);
1421  if (_debug_output_level > 1)
1422  _console << "XFEM set material properties for element: " << libmesh_elem->id()
1423  << "\n";
1424  break;
1425  }
1426  }
1427 
1428  // Store solution for all elements affected by XFEM
1429  storeSolutionForElement(libmesh_elem,
1430  parent_elem,
1431  *nls[0],
1433  current_solution,
1434  old_solution,
1435  older_solution);
1436  storeSolutionForElement(libmesh_elem,
1437  parent_elem,
1438  aux,
1440  current_aux_solution,
1441  old_aux_solution,
1442  older_aux_solution);
1443  }
1444  }
1445 
1446  // delete elements
1447  for (std::size_t i = 0; i < delete_elements.size(); ++i)
1448  {
1449  Elem * elem_to_delete = _mesh->elem_ptr(delete_elements[i]->id());
1450 
1451  // delete the XFEMCutElem object for any elements that are to be deleted
1452  std::map<unique_id_type, XFEMCutElem *>::iterator cemit =
1453  _cut_elem_map.find(elem_to_delete->unique_id());
1454  if (cemit != _cut_elem_map.end())
1455  {
1456  delete cemit->second;
1457  _cut_elem_map.erase(cemit);
1458  }
1459 
1460  // remove the property storage of deleted element/side
1461  _material_data[0]->eraseProperty(elem_to_delete);
1462  _bnd_material_data[0]->eraseProperty(elem_to_delete);
1463 
1464  elem_to_delete->nullify_neighbors();
1465  _mesh->get_boundary_info().remove(elem_to_delete);
1466  unsigned int deleted_elem_id = elem_to_delete->id();
1467  _mesh->delete_elem(elem_to_delete);
1468  if (_debug_output_level > 1)
1469  _console << "XFEM deleted element: " << deleted_elem_id << std::endl;
1470 
1471  if (_displaced_mesh)
1472  {
1473  Elem * elem_to_delete2 = _displaced_mesh->elem_ptr(delete_elements[i]->id());
1474  elem_to_delete2->nullify_neighbors();
1475  _displaced_mesh->get_boundary_info().remove(elem_to_delete2);
1476  _displaced_mesh->delete_elem(elem_to_delete2);
1477  }
1478  }
1479 
1480  for (std::map<unsigned int, std::vector<const Elem *>>::iterator it =
1481  temporary_parent_children_map.begin();
1482  it != temporary_parent_children_map.end();
1483  ++it)
1484  {
1485  std::vector<const Elem *> & sibling_elem_vec = it->second;
1486  // TODO: for cut-node case, how to find the sibling elements?
1487  // if (sibling_elem_vec.size() != 2)
1488  // mooseError("Must have exactly 2 sibling elements");
1489 
1490  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1491  for (auto const & elem_id : _geom_marker_id_elems[_geometric_cuts[i]->getInterfaceID()])
1492  if (it->first == elem_id)
1493  _sibling_elems[_geometric_cuts[i]->getInterfaceID()].push_back(
1494  std::make_pair(sibling_elem_vec[0], sibling_elem_vec[1]));
1495  }
1496 
1497  // add sibling elems on displaced mesh
1498  if (_displaced_mesh)
1499  {
1500  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1501  {
1502  for (auto & se : _sibling_elems[_geometric_cuts[i]->getInterfaceID()])
1503  {
1504  Elem * elem = _displaced_mesh->elem_ptr(se.first->id());
1505  Elem * elem_pair = _displaced_mesh->elem_ptr(se.second->id());
1506  _sibling_displaced_elems[_geometric_cuts[i]->getInterfaceID()].push_back(
1507  std::make_pair(elem, elem_pair));
1508  }
1509  }
1510  }
1511 
1512  // clear the temporary map
1513  temporary_parent_children_map.clear();
1514 
1515  // Store information about crack tip elements
1516  if (mesh_changed)
1517  {
1518  _crack_tip_elems.clear();
1520  const std::set<EFAElement *> CrackTipElements = _efa_mesh.getCrackTipElements();
1521  std::set<EFAElement *>::const_iterator sit;
1522  for (sit = CrackTipElements.begin(); sit != CrackTipElements.end(); ++sit)
1523  {
1524  unsigned int eid = (*sit)->id();
1525  Elem * crack_tip_elem;
1526  std::map<unsigned int, Elem *>::iterator eit = efa_id_to_new_elem.find(eid);
1527  if (eit != efa_id_to_new_elem.end())
1528  crack_tip_elem = eit->second;
1529  else
1530  crack_tip_elem = _mesh->elem_ptr(eid);
1531  _crack_tip_elems.insert(crack_tip_elem);
1532 
1533  // Store the crack tip elements which are going to be healed
1534  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1535  {
1536  if (_geometric_cuts[i]->shouldHealMesh())
1537  {
1538  for (auto const & mie : _geom_marker_id_elems[_geometric_cuts[i]->getInterfaceID()])
1539  if ((*sit)->getParent() != nullptr)
1540  {
1541  if (_mesh->mesh_dimension() == 2)
1542  {
1543  EFAElement2D * efa_elem2d = dynamic_cast<EFAElement2D *>((*sit)->getParent());
1544  if (!efa_elem2d)
1545  mooseError("EFAelem is not of EFAelement2D type");
1546 
1547  for (unsigned int edge_id = 0; edge_id < efa_elem2d->numEdges(); ++edge_id)
1548  {
1549  for (unsigned int en_iter = 0; en_iter < efa_elem2d->numEdgeNeighbors(edge_id);
1550  ++en_iter)
1551  {
1552  EFAElement2D * edge_neighbor = efa_elem2d->getEdgeNeighbor(edge_id, en_iter);
1553  if (edge_neighbor != nullptr && edge_neighbor->id() == mie)
1554  _crack_tip_elems_to_be_healed.insert(crack_tip_elem);
1555  }
1556  }
1557  }
1558  else if (_mesh->mesh_dimension() == 3)
1559  {
1560  EFAElement3D * efa_elem3d = dynamic_cast<EFAElement3D *>((*sit)->getParent());
1561  if (!efa_elem3d)
1562  mooseError("EFAelem is not of EFAelement3D type");
1563 
1564  for (unsigned int face_id = 0; face_id < efa_elem3d->numFaces(); ++face_id)
1565  {
1566  for (unsigned int fn_iter = 0; fn_iter < efa_elem3d->numFaceNeighbors(face_id);
1567  ++fn_iter)
1568  {
1569  EFAElement3D * face_neighbor = efa_elem3d->getFaceNeighbor(face_id, fn_iter);
1570  if (face_neighbor != nullptr && face_neighbor->id() == mie)
1571  _crack_tip_elems_to_be_healed.insert(crack_tip_elem);
1572  }
1573  }
1574  }
1575  }
1576  }
1577  }
1578  }
1579  }
1580 
1581  if (_debug_output_level > 0)
1582  {
1583  _console << "\nXFEM mesh cutting with element fragment algorithm complete\n";
1584  _console << "# new nodes: " << new_nodes.size() << "\n";
1585  _console << "# new elements: " << new_elements.size() << "\n";
1586  _console << "# deleted elements: " << delete_elements.size() << "\n";
1587  _console << std::flush;
1588  }
1589 
1590  // store virtual nodes
1591  // store cut edge info
1592  return mesh_changed;
1593 }
unsigned int _debug_output_level
Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mes...
Definition: XFEM.h:362
const GeometricCutUserObject * getGeometricCutForElem(const Elem *elem) const
Get the GeometricCutUserObject associated with an element.
Definition: XFEM.C:2155
unsigned int numEdgeNeighbors(unsigned int edge_id) const
EFAElement3D * getFaceNeighbor(unsigned int face_id, unsigned int neighbor_id) const
CutSubdomainID getCutSubdomainID(const GeometricCutUserObject *gcuo, const Elem *cut_elem, const Elem *parent_elem=nullptr) const
Determine which cut subdomain the element belongs to relative to the cut.
Definition: XFEM.C:2271
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:327
unsigned int id() const
Definition: EFAElement.C:28
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:332
bool shouldHealMesh() const
Should the elements cut by this cutting object be healed in the current time step?
void mooseError(Args &&... args)
bool isSemiLocal(Node *const node) const
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:333
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
bool needBoundaryMaterialOnSide(BoundaryID bnd_id, const THREAD_ID tid)
const std::vector< EFAElement * > & getChildElements()
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:330
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
NumericVector< Number > & solutionOlder()
unsigned int numEdges() const
std::set< const Elem * > _crack_tip_elems_to_be_healed
Definition: XFEM.h:331
std::map< unique_id_type, std::vector< Real > > _cached_aux_solution
Data structure to store the auxiliary solution for nodes/elements affected by XFEM For each node/elem...
Definition: XFEM.h:380
unsigned int numFaces() const
unsigned int CutSubdomainID
Definition: XFEMAppTypes.h:18
FEProblemBase * _fe_problem
std::vector< MaterialData *> _bnd_material_data
void updateTopology(bool mergeUncutVirtualEdges=true)
const std::vector< EFAElement * > & getParentElements()
void loadMaterialPropertiesForElement(const Elem *elem, const Elem *elem_from, std::unordered_map< const Elem *, Xfem::CutElemInfo > &cached_cei) const
Helper function to store the material properties of a healed element.
Definition: XFEM.C:2237
const std::vector< EFANode * > & getNewNodes()
std::vector< MaterialData *> _material_data
void storeSolutionForElement(const Elem *elem_to_store_to, const Elem *elem_to_store_from, SystemBase &sys, std::map< unique_id_type, std::vector< Real >> &stored_solution, const NumericVector< Number > &current_solution, const NumericVector< Number > &old_solution, const NumericVector< Number > &older_solution)
Store the solution in stored_solution for a given element.
Definition: XFEM.C:2002
const QBase *const & qRule() const
void storeSolutionForNode(const Node *node_to_store_to, const Node *node_to_store_from, SystemBase &sys, std::map< unique_id_type, std::vector< Real >> &stored_solution, const NumericVector< Number > &current_solution, const NumericVector< Number > &old_solution, const NumericVector< Number > &older_solution)
Store the solution in stored_solution for a given node.
Definition: XFEM.C:1968
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
Data structure for storing the elements cut by specific geometric cutters.
Definition: XFEM.h:350
MeshBase * _mesh
EFAElement2D * getEdgeNeighbor(unsigned int edge_id, unsigned int neighbor_id) const
MooseMesh * _moose_mesh
virtual System & system() override
MeshBase * _displaced_mesh
std::map< unique_id_type, std::vector< Real > > _cached_solution
Data structure to store the nonlinear solution for nodes/elements affected by XFEM For each node/elem...
Definition: XFEM.h:371
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
const ConsoleStream _console
virtual void serializeSolution()
NumericVector< Number > & solutionOld()
Assembly & assembly(const THREAD_ID tid, const unsigned int nl_sys_num) override
Information about a cut element.
Definition: XFEM.h:57
const std::set< EFAElement * > & getCrackTipElements()
unsigned int numFaceNeighbors(unsigned int face_id) const
std::unordered_map< const Elem *, Xfem::CutElemInfo > _geom_cut_elems
All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK...
Definition: XFEM.h:386
std::unordered_map< const Elem *, Xfem::CutElemInfo > _old_geom_cut_elems
All geometrically cut elements and their CutElemInfo before the current execution of XFEM_MARK...
Definition: XFEM.h:392

◆ getCrackTipOrigin()

void XFEM::getCrackTipOrigin ( std::map< unsigned int, const Elem *> &  elem_id_crack_tip,
std::vector< Point > &  crack_front_points 
)

Definition at line 65 of file XFEM.C.

67 {
68  elem_id_crack_tip.clear();
69  crack_front_points.clear();
70  crack_front_points.resize(_elem_crack_origin_direction_map.size());
71 
72  unsigned int crack_tip_index = 0;
73  // This map is used to sort the order in _elem_crack_origin_direction_map such that every process
74  // has same order
75  std::map<unsigned int, const Elem *> elem_id_map;
76 
77  int m = -1;
78  for (std::map<const Elem *, std::vector<Point>>::iterator mit1 =
81  ++mit1)
82  {
83  unsigned int elem_id = mit1->first->id();
84  if (elem_id == std::numeric_limits<unsigned int>::max())
85  {
86  elem_id_map[m] = mit1->first;
87  m--;
88  }
89  else
90  elem_id_map[elem_id] = mit1->first;
91  }
92 
93  for (std::map<unsigned int, const Elem *>::iterator mit1 = elem_id_map.begin();
94  mit1 != elem_id_map.end();
95  mit1++)
96  {
97  const Elem * elem = mit1->second;
98  std::map<const Elem *, std::vector<Point>>::iterator mit2 =
100  if (mit2 != _elem_crack_origin_direction_map.end())
101  {
102  elem_id_crack_tip[crack_tip_index] = mit2->first;
103  crack_front_points[crack_tip_index] =
104  (mit2->second)[0]; // [0] stores origin coordinates and [1] stores direction
105  crack_tip_index++;
106  }
107  }
108 }
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335

◆ getCrackTipOriginMap()

const std::map<const Elem *, std::vector<Point> >& XFEM::getCrackTipOriginMap ( ) const
inline

Definition at line 304 of file XFEM.h.

305  {
307  }
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335

◆ getCutPlane()

Real XFEM::getCutPlane ( const Elem *  elem,
const Xfem::XFEM_CUTPLANE_QUANTITY  quantity,
unsigned int  plane_id 
) const

Get specified component of normal or origin for cut plane for a given element.

Definition at line 1666 of file XFEM.C.

1669 {
1670  Real comp = 0.0;
1671  Point planedata(0.0, 0.0, 0.0);
1672  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1673  it = _cut_elem_map.find(elem->unique_id());
1674  if (it != _cut_elem_map.end())
1675  {
1676  const XFEMCutElem * xfce = it->second;
1677  const EFAElement * EFAelem = xfce->getEFAElement();
1678  if (EFAelem->isPartial()) // exclude the full crack tip elements
1679  {
1680  if ((unsigned int)quantity < 3)
1681  {
1682  unsigned int index = (unsigned int)quantity;
1683  planedata = xfce->getCutPlaneOrigin(plane_id, _displaced_mesh);
1684  comp = planedata(index);
1685  }
1686  else if ((unsigned int)quantity < 6)
1687  {
1688  unsigned int index = (unsigned int)quantity - 3;
1689  planedata = xfce->getCutPlaneNormal(plane_id, _displaced_mesh);
1690  comp = planedata(index);
1691  }
1692  else
1693  mooseError("In get_cut_plane index out of range");
1694  }
1695  }
1696  return comp;
1697 }
void mooseError(Args &&... args)
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
virtual bool isPartial() const =0
virtual const EFAElement * getEFAElement() const =0
virtual Point getCutPlaneOrigin(unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const =0
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MeshBase * _displaced_mesh
virtual Point getCutPlaneNormal(unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const =0
void ErrorVector unsigned int

◆ getCutSubdomainID()

CutSubdomainID XFEM::getCutSubdomainID ( const GeometricCutUserObject gcuo,
const Elem *  cut_elem,
const Elem *  parent_elem = nullptr 
) const

Determine which cut subdomain the element belongs to relative to the cut.

Parameters
gcuoThe GeometricCutUserObject for the cut
cut_elemThe element being cut
parent_elemThe parent element

Definition at line 2271 of file XFEM.C.

Referenced by cutMeshWithEFA().

2274 {
2275  if (!parent_elem)
2276  parent_elem = cut_elem;
2277  // Pick any node from the parent element that is inside the physical domain and return its
2278  // CutSubdomainID.
2279  const Node * node = pickFirstPhysicalNode(cut_elem, parent_elem);
2280  return gcuo->getCutSubdomainID(node);
2281 }
const Node * pickFirstPhysicalNode(const Elem *e, const Elem *e0) const
Return the first node in the provided element that is found to be in the physical domain...
Definition: XFEM.C:2284
virtual CutSubdomainID getCutSubdomainID(const Node *) const
Get CutSubdomainID telling which side the node belongs to relative to the cut.

◆ getEFAElem2D()

EFAElement2D * XFEM::getEFAElem2D ( const Elem *  elem)

Get the EFAElement2D object for a specified libMesh element.

Parameters
elemPointer to the libMesh element for which the object is requested

Definition at line 1746 of file XFEM.C.

Referenced by markCutEdgesByGeometry(), and markCutEdgesByState().

1747 {
1748  EFAElement * EFAelem = _efa_mesh.getElemByID(elem->id());
1749  EFAElement2D * EFAelem2D = dynamic_cast<EFAElement2D *>(EFAelem);
1750 
1751  if (!EFAelem2D)
1752  mooseError("EFAelem is not of EFAelement2D type");
1753 
1754  return EFAelem2D;
1755 }
void mooseError(Args &&... args)
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
EFAElement * getElemByID(unsigned int id)

◆ getEFAElem3D()

EFAElement3D * XFEM::getEFAElem3D ( const Elem *  elem)

Get the EFAElement3D object for a specified libMesh element.

Parameters
elemPointer to the libMesh element for which the object is requested

Definition at line 1758 of file XFEM.C.

Referenced by markCutFacesByGeometry().

1759 {
1760  EFAElement * EFAelem = _efa_mesh.getElemByID(elem->id());
1761  EFAElement3D * EFAelem3D = dynamic_cast<EFAElement3D *>(EFAelem);
1762 
1763  if (!EFAelem3D)
1764  mooseError("EFAelem is not of EFAelement3D type");
1765 
1766  return EFAelem3D;
1767 }
void mooseError(Args &&... args)
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
EFAElement * getElemByID(unsigned int id)

◆ getEFANodeCoords()

Point XFEM::getEFANodeCoords ( EFANode CEMnode,
EFAElement CEMElem,
const Elem *  elem,
MeshBase *  displaced_mesh = nullptr 
) const

Definition at line 1596 of file XFEM.C.

Referenced by correctCrackExtensionDirection(), getFragmentEdges(), getFragmentFaces(), and markCutEdgesByState().

1600 {
1601  Point node_coor(0.0, 0.0, 0.0);
1602  std::vector<EFANode *> master_nodes;
1603  std::vector<Point> master_points;
1604  std::vector<double> master_weights;
1605 
1606  CEMElem->getMasterInfo(CEMnode, master_nodes, master_weights);
1607  for (std::size_t i = 0; i < master_nodes.size(); ++i)
1608  {
1609  if (master_nodes[i]->category() == EFANode::N_CATEGORY_PERMANENT)
1610  {
1611  unsigned int local_node_id = CEMElem->getLocalNodeIndex(master_nodes[i]);
1612  const Node * node = elem->node_ptr(local_node_id);
1613  if (displaced_mesh)
1614  node = displaced_mesh->node_ptr(node->id());
1615  Point node_p((*node)(0), (*node)(1), (*node)(2));
1616  master_points.push_back(node_p);
1617  }
1618  else
1619  mooseError("master nodes must be permanent");
1620  }
1621  for (std::size_t i = 0; i < master_nodes.size(); ++i)
1622  node_coor += master_weights[i] * master_points[i];
1623 
1624  return node_coor;
1625 }
void mooseError(Args &&... args)
virtual void getMasterInfo(EFANode *node, std::vector< EFANode *> &master_nodes, std::vector< double > &master_weights) const =0
unsigned int getLocalNodeIndex(EFANode *node) const
Definition: EFAElement.C:111

◆ getElementSolutionDofs()

std::vector< dof_id_type > XFEM::getElementSolutionDofs ( const Elem *  elem,
SystemBase sys 
) const
private

Get a vector of the dof indices for all components of all variables associated with an element.

Parameters
elemElement for which dof indices are found
sysSystem for which the dof indices are found

Definition at line 2097 of file XFEM.C.

Referenced by setSolution(), and storeSolutionForElement().

2098 {
2099  SubdomainID sid = elem->subdomain_id();
2100  const std::vector<MooseVariableFEBase *> & vars = sys.getVariables(0);
2101  std::vector<dof_id_type> solution_dofs;
2102  solution_dofs.reserve(vars.size()); // just an approximation
2103  for (auto var : vars)
2104  {
2105  if (!var->isNodal())
2106  {
2107  const std::set<SubdomainID> & var_subdomains = sys.getSubdomainsForVar(var->number());
2108  if (var_subdomains.empty() || var_subdomains.find(sid) != var_subdomains.end())
2109  {
2110  unsigned int n_comp = elem->n_comp(sys.number(), var->number());
2111  for (unsigned int icomp = 0; icomp < n_comp; ++icomp)
2112  {
2113  dof_id_type elem_dof = elem->dof_number(sys.number(), var->number(), icomp);
2114  solution_dofs.push_back(elem_dof);
2115  }
2116  }
2117  }
2118  }
2119  return solution_dofs;
2120 }
const std::vector< MooseVariableFieldBase *> & getVariables(THREAD_ID tid)
subdomain_id_type SubdomainID
unsigned int number() const
const std::set< SubdomainID > & getSubdomainsForVar(unsigned int var_number) const
uint8_t dof_id_type

◆ getFragmentEdges()

void XFEM::getFragmentEdges ( const Elem *  elem,
EFAElement2D CEMElem,
std::vector< std::vector< Point >> &  frag_edges 
) const

Definition at line 1770 of file XFEM.C.

1773 {
1774  // N.B. CEMElem here has global EFAnode
1775  frag_edges.clear();
1776  if (CEMElem->numFragments() > 0)
1777  {
1778  if (CEMElem->numFragments() > 1)
1779  mooseError("element ", elem->id(), " has more than one fragment at this point");
1780  for (unsigned int i = 0; i < CEMElem->getFragment(0)->numEdges(); ++i)
1781  {
1782  std::vector<Point> p_line(2, Point(0.0, 0.0, 0.0));
1783  p_line[0] = getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(0), CEMElem, elem);
1784  p_line[1] = getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(1), CEMElem, elem);
1785  frag_edges.push_back(p_line);
1786  }
1787  }
1788 }
void mooseError(Args &&... args)
unsigned int numEdges() const
virtual unsigned int numFragments() const
Definition: EFAElement2D.C:207
EFAFragment2D * getFragment(unsigned int frag_id) const
EFANode * getNode(unsigned int index) const
Definition: EFAEdge.C:181
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=nullptr) const
Definition: XFEM.C:1596
EFAEdge * getFragmentEdge(unsigned int frag_id, unsigned int edge_id) const

◆ getFragmentFaces() [1/2]

void XFEM::getFragmentFaces ( const Elem *  elem,
std::vector< std::vector< Point >> &  frag_faces,
bool  displaced_mesh = false 
) const

Definition at line 1729 of file XFEM.C.

1732 {
1733  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1734  it = _cut_elem_map.find(elem->unique_id());
1735  if (it != _cut_elem_map.end())
1736  {
1737  const XFEMCutElem * xfce = it->second;
1738  if (displaced_mesh)
1739  xfce->getFragmentFaces(frag_faces, _displaced_mesh);
1740  else
1741  xfce->getFragmentFaces(frag_faces);
1742  }
1743 }
virtual void getFragmentFaces(std::vector< std::vector< Point >> &frag_faces, MeshBase *displaced_mesh=nullptr) const =0
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
MeshBase * _displaced_mesh

◆ getFragmentFaces() [2/2]

void XFEM::getFragmentFaces ( const Elem *  elem,
EFAElement3D CEMElem,
std::vector< std::vector< Point >> &  frag_faces 
) const

Definition at line 1791 of file XFEM.C.

1794 {
1795  // N.B. CEMElem here has global EFAnode
1796  frag_faces.clear();
1797  if (CEMElem->numFragments() > 0)
1798  {
1799  if (CEMElem->numFragments() > 1)
1800  mooseError("element ", elem->id(), " has more than one fragment at this point");
1801  for (unsigned int i = 0; i < CEMElem->getFragment(0)->numFaces(); ++i)
1802  {
1803  unsigned int num_face_nodes = CEMElem->getFragmentFace(0, i)->numNodes();
1804  std::vector<Point> p_line(num_face_nodes, Point(0.0, 0.0, 0.0));
1805  for (unsigned int j = 0; j < num_face_nodes; ++j)
1806  p_line[j] = getEFANodeCoords(CEMElem->getFragmentFace(0, i)->getNode(j), CEMElem, elem);
1807  frag_faces.push_back(p_line);
1808  }
1809  }
1810 }
EFAFragment3D * getFragment(unsigned int frag_id) const
void mooseError(Args &&... args)
virtual unsigned int numFragments() const
Definition: EFAElement3D.C:272
EFANode * getNode(unsigned int node_id) const
Definition: EFAFace.C:99
unsigned int numNodes() const
Definition: EFAFace.C:87
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
EFAFace * getFragmentFace(unsigned int frag_id, unsigned int face_id) const
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=nullptr) const
Definition: XFEM.C:1596
unsigned int numFaces() const

◆ getGeometricCutForElem()

const GeometricCutUserObject * XFEM::getGeometricCutForElem ( const Elem *  elem) const
private

Get the GeometricCutUserObject associated with an element.

Parameters
elemThe element
Returns
A constant pointer to the GeometricCutUserObject, nullptr if nothing found

Definition at line 2155 of file XFEM.C.

Referenced by cutMeshWithEFA().

2156 {
2157  for (auto gcmit : _geom_marker_id_elems)
2158  {
2159  std::set<unsigned int> elems = gcmit.second;
2160  if (elems.find(elem->id()) != elems.end())
2161  return _geometric_cuts[gcmit.first];
2162  }
2163  return nullptr;
2164 }
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:327
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
Data structure for storing the elements cut by specific geometric cutters.
Definition: XFEM.h:350

◆ getGeometricCutID()

virtual unsigned int XFEM::getGeometricCutID ( const GeometricCutUserObject gcu)
inlinevirtual

Get the interface ID corresponding to a given GeometricCutUserObject.

Parameters
gcupointer to the GeometricCutUserObject
Returns
the interface ID

Definition at line 267 of file XFEM.h.

268  {
269  return _geom_marker_id_map[gcu];
270  };
std::map< const GeometricCutUserObject *, unsigned int > _geom_marker_id_map
Data structure for storing the GeommetricCutUserObjects and their corresponding id.
Definition: XFEM.h:353

◆ getNodeSolutionDofs()

std::vector< dof_id_type > XFEM::getNodeSolutionDofs ( const Node *  node,
SystemBase sys 
) const
private

Get a vector of the dof indices for all components of all variables associated with a node.

Parameters
nodeNode for which dof indices are found
sysSystem for which the dof indices are found

Definition at line 2123 of file XFEM.C.

Referenced by setSolution(), and storeSolutionForNode().

2124 {
2125  const std::set<SubdomainID> & sids = _moose_mesh->getNodeBlockIds(*node);
2126  const std::vector<MooseVariableFEBase *> & vars = sys.getVariables(0);
2127  std::vector<dof_id_type> solution_dofs;
2128  solution_dofs.reserve(vars.size()); // just an approximation
2129  for (auto var : vars)
2130  {
2131  if (var->isNodal())
2132  {
2133  const std::set<SubdomainID> & var_subdomains = sys.getSubdomainsForVar(var->number());
2134  std::set<SubdomainID> intersect;
2135  set_intersection(var_subdomains.begin(),
2136  var_subdomains.end(),
2137  sids.begin(),
2138  sids.end(),
2139  std::inserter(intersect, intersect.begin()));
2140  if (var_subdomains.empty() || !intersect.empty())
2141  {
2142  unsigned int n_comp = node->n_comp(sys.number(), var->number());
2143  for (unsigned int icomp = 0; icomp < n_comp; ++icomp)
2144  {
2145  dof_id_type node_dof = node->dof_number(sys.number(), var->number(), icomp);
2146  solution_dofs.push_back(node_dof);
2147  }
2148  }
2149  }
2150  }
2151  return solution_dofs;
2152 }
const std::vector< MooseVariableFieldBase *> & getVariables(THREAD_ID tid)
const std::set< SubdomainID > & getNodeBlockIds(const Node &node) const
unsigned int number() const
MooseMesh * _moose_mesh
const std::set< SubdomainID > & getSubdomainsForVar(unsigned int var_number) const
uint8_t dof_id_type

◆ getPhysicalVolumeFraction()

Real XFEM::getPhysicalVolumeFraction ( const Elem *  elem) const

Get the volume fraction of an element that is physical.

Definition at line 1628 of file XFEM.C.

Referenced by markCutEdgesByState().

1629 {
1630  Real phys_volfrac = 1.0;
1631  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1632  it = _cut_elem_map.find(elem->unique_id());
1633  if (it != _cut_elem_map.end())
1634  {
1635  XFEMCutElem * xfce = it->second;
1636  const EFAElement * EFAelem = xfce->getEFAElement();
1637  if (EFAelem->isPartial())
1638  { // exclude the full crack tip elements
1640  phys_volfrac = xfce->getPhysicalVolumeFraction();
1641  }
1642  }
1643 
1644  return phys_volfrac;
1645 }
Real getPhysicalVolumeFraction() const
Returns the volume fraction of the element fragment.
Definition: XFEMCutElem.C:51
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
virtual bool isPartial() const =0
virtual const EFAElement * getEFAElement() const =0
virtual void computePhysicalVolumeFraction()=0
Computes the volume fraction of the element fragment.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ getXFEMCutElemPairs()

virtual const ElementPairLocator::ElementPairList* XFEM::getXFEMCutElemPairs ( unsigned int  interface_id)
inlinevirtual

Get the list of cut element pairs corresponding to a given interface ID.

Parameters
interface_idThe ID of the interface
Returns
the list of elements cut by that interface

Definition at line 245 of file XFEM.h.

246  {
247  return &_sibling_elems[interface_id];
248  }
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:332

◆ getXFEMDisplacedCutElemPairs()

virtual const ElementPairLocator::ElementPairList* XFEM::getXFEMDisplacedCutElemPairs ( unsigned int  interface_id)
inlinevirtual

Get the list of cut element pairs on the displaced mesh corresponding to a given interface ID.

Parameters
interface_idThe ID of the interface
Returns
the list of elements cut by that interface

Definition at line 257 of file XFEM.h.

258  {
259  return &_sibling_displaced_elems[interface_id];
260  }
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:333

◆ getXFEMFaceWeights()

bool XFEM::getXFEMFaceWeights ( MooseArray< Real > &  weights,
const Elem *  elem,
QBase *  qrule,
const MooseArray< Point > &  q_points,
unsigned int  side 
)
overridevirtual

Implements XFEMInterface.

Definition at line 1860 of file XFEM.C.

1865 {
1866  bool have_weights = false;
1867  XFEMCutElem * xfce = nullptr;
1868  if (isElemCut(elem, xfce))
1869  {
1870  mooseAssert(xfce != nullptr, "Must have valid XFEMCutElem object here");
1871  xfce->getFaceWeightMultipliers(weights, qrule, getXFEMQRule(), q_points, side);
1872  have_weights = true;
1873  }
1874  return have_weights;
1875 }
void getFaceWeightMultipliers(MooseArray< Real > &face_weights, QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points, unsigned int side)
Definition: XFEMCutElem.C:77
Xfem::XFEM_QRULE & getXFEMQRule()
Definition: XFEM.C:1813
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1706

◆ getXFEMIntersectionInfo()

void XFEM::getXFEMIntersectionInfo ( const Elem *  elem,
unsigned int  plane_id,
Point &  normal,
std::vector< Point > &  intersectionPoints,
bool  displaced_mesh = false 
) const
virtual

Definition at line 1878 of file XFEM.C.

1883 {
1884  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1885  it = _cut_elem_map.find(elem->unique_id());
1886  if (it != _cut_elem_map.end())
1887  {
1888  const XFEMCutElem * xfce = it->second;
1889  if (displaced_mesh)
1890  xfce->getIntersectionInfo(plane_id, normal, intersectionPoints, _displaced_mesh);
1891  else
1892  xfce->getIntersectionInfo(plane_id, normal, intersectionPoints);
1893  }
1894 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
virtual void getIntersectionInfo(unsigned int plane_id, Point &normal, std::vector< Point > &intersectionPoints, MeshBase *displaced_mesh=nullptr) const =0
MeshBase * _displaced_mesh

◆ getXFEMQRule()

Xfem::XFEM_QRULE & XFEM::getXFEMQRule ( )

Definition at line 1813 of file XFEM.C.

Referenced by getXFEMFaceWeights(), and getXFEMWeights().

1814 {
1815  return _XFEM_qrule;
1816 }
Xfem::XFEM_QRULE _XFEM_qrule
Definition: XFEM.h:322

◆ getXFEMqRuleOnLine()

void XFEM::getXFEMqRuleOnLine ( std::vector< Point > &  intersection_points,
std::vector< Point > &  quad_pts,
std::vector< Real > &  quad_wts 
) const
virtual

Definition at line 1897 of file XFEM.C.

1900 {
1901  Point p1 = intersection_points[0];
1902  Point p2 = intersection_points[1];
1903 
1904  // number of quadrature points
1905  std::size_t num_qpoints = 2;
1906 
1907  // quadrature coordinates
1908  Real xi0 = -std::sqrt(1.0 / 3.0);
1909  Real xi1 = std::sqrt(1.0 / 3.0);
1910 
1911  quad_wts.resize(num_qpoints);
1912  quad_pts.resize(num_qpoints);
1913 
1914  Real integ_jacobian = pow((p1 - p2).norm_sq(), 0.5) * 0.5;
1915 
1916  quad_wts[0] = 1.0 * integ_jacobian;
1917  quad_wts[1] = 1.0 * integ_jacobian;
1918 
1919  quad_pts[0] = (1.0 - xi0) / 2.0 * p1 + (1.0 + xi0) / 2.0 * p2;
1920  quad_pts[1] = (1.0 - xi1) / 2.0 * p1 + (1.0 + xi1) / 2.0 * p2;
1921 }
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto norm_sq(const T &a) -> decltype(std::norm(a))

◆ getXFEMqRuleOnSurface()

void XFEM::getXFEMqRuleOnSurface ( std::vector< Point > &  intersection_points,
std::vector< Point > &  quad_pts,
std::vector< Real > &  quad_wts 
) const
virtual

Definition at line 1924 of file XFEM.C.

1927 {
1928  std::size_t nnd_pe = intersection_points.size();
1929  Point xcrd(0.0, 0.0, 0.0);
1930  for (std::size_t i = 0; i < nnd_pe; ++i)
1931  xcrd += intersection_points[i];
1932  xcrd /= nnd_pe;
1933 
1934  quad_pts.resize(nnd_pe);
1935  quad_wts.resize(nnd_pe);
1936 
1937  Real jac = 0.0;
1938 
1939  for (std::size_t j = 0; j < nnd_pe; ++j) // loop all sub-tris
1940  {
1941  std::vector<std::vector<Real>> shape(3, std::vector<Real>(3, 0.0));
1942  std::vector<Point> subtrig_points(3, Point(0.0, 0.0, 0.0)); // sub-trig nodal coords
1943 
1944  int jplus1 = j < nnd_pe - 1 ? j + 1 : 0;
1945  subtrig_points[0] = xcrd;
1946  subtrig_points[1] = intersection_points[j];
1947  subtrig_points[2] = intersection_points[jplus1];
1948 
1949  std::vector<std::vector<Real>> sg2;
1950  Xfem::stdQuadr2D(3, 1, sg2); // get sg2
1951  for (std::size_t l = 0; l < sg2.size(); ++l) // loop all int pts on a sub-trig
1952  {
1953  Xfem::shapeFunc2D(3, sg2[l], subtrig_points, shape, jac, true); // Get shape
1954  std::vector<Real> tsg_line(3, 0.0);
1955  for (std::size_t k = 0; k < 3; ++k) // loop sub-trig nodes
1956  {
1957  tsg_line[0] += shape[k][2] * subtrig_points[k](0);
1958  tsg_line[1] += shape[k][2] * subtrig_points[k](1);
1959  tsg_line[2] += shape[k][2] * subtrig_points[k](2);
1960  }
1961  quad_pts[j + l] = Point(tsg_line[0], tsg_line[1], tsg_line[2]);
1962  quad_wts[j + l] = sg2[l][3] * jac;
1963  }
1964  }
1965 }
void shapeFunc2D(unsigned int nen, std::vector< Real > &ss, std::vector< Point > &xl, std::vector< std::vector< Real >> &shp, Real &xsj, bool natl_flg)
Definition: XFEMFuncs.C:269
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void stdQuadr2D(unsigned int nen, unsigned int iord, std::vector< std::vector< Real >> &sg2)
Definition: XFEMFuncs.C:96
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
static const std::string k
Definition: NS.h:124

◆ getXFEMWeights()

bool XFEM::getXFEMWeights ( MooseArray< Real > &  weights,
const Elem *  elem,
QBase *  qrule,
const MooseArray< Point > &  q_points 
)
overridevirtual

Implements XFEMInterface.

Definition at line 1843 of file XFEM.C.

1847 {
1848  bool have_weights = false;
1849  XFEMCutElem * xfce = nullptr;
1850  if (isElemCut(elem, xfce))
1851  {
1852  mooseAssert(xfce != nullptr, "Must have valid XFEMCutElem object here");
1853  xfce->getWeightMultipliers(weights, qrule, getXFEMQRule(), q_points);
1854  have_weights = true;
1855  }
1856  return have_weights;
1857 }
Xfem::XFEM_QRULE & getXFEMQRule()
Definition: XFEM.C:1813
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1706
void getWeightMultipliers(MooseArray< Real > &weights, QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
Definition: XFEMCutElem.C:63

◆ has_secondary_cut()

bool XFEM::has_secondary_cut ( )
inline

Definition at line 283 of file XFEM.h.

283 { return _has_secondary_cut; }
bool _has_secondary_cut
Definition: XFEM.h:320

◆ healMesh()

bool XFEM::healMesh ( )

Potentially heal the mesh by merging some of the pairs of partial elements cut by XFEM back into single elements if indicated by the cutting objects.

Returns
true if the mesh has been modified due to healing

Definition at line 917 of file XFEM.C.

Referenced by updateHeal().

918 {
919  bool mesh_changed = false;
920 
921  std::set<Node *> nodes_to_delete;
922  std::set<Node *> nodes_to_delete_displaced;
923  std::set<unsigned int> cutelems_to_delete;
924  unsigned int deleted_elem_count = 0;
925  std::vector<std::string> healed_geometric_cuts;
926 
927  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
928  {
929  if (_geometric_cuts[i]->shouldHealMesh())
930  {
931  healed_geometric_cuts.push_back(_geometric_cuts[i]->name());
932  for (auto & it : _sibling_elems[_geometric_cuts[i]->getInterfaceID()])
933  {
934  Elem * elem1 = const_cast<Elem *>(it.first);
935  Elem * elem2 = const_cast<Elem *>(it.second);
936 
937  std::map<unique_id_type, XFEMCutElem *>::iterator cemit =
938  _cut_elem_map.find(elem1->unique_id());
939  if (cemit != _cut_elem_map.end())
940  {
941  const XFEMCutElem * xfce = cemit->second;
942 
943  cutelems_to_delete.insert(elem1->unique_id());
944 
945  for (unsigned int in = 0; in < elem1->n_nodes(); ++in)
946  {
947  Node * e1node = elem1->node_ptr(in);
948  Node * e2node = elem2->node_ptr(in);
949  if (!xfce->isPointPhysical(*e1node) &&
950  e1node != e2node) // This would happen at the crack tip
951  {
952  elem1->set_node(in) = e2node;
953  nodes_to_delete.insert(e1node);
954  }
955  else if (e1node != e2node)
956  nodes_to_delete.insert(e2node);
957  }
958  }
959  else
960  mooseError("Could not find XFEMCutElem for element to be kept in healing");
961 
962  // Store the material properties of the elements to be healed. So that if the element is
963  // immediately re-cut, we can restore the material properties (especially those stateful
964  // ones).
965  std::vector<const Elem *> healed_elems = {elem1, elem2};
966 
967  if (_geometric_cuts[i]->shouldHealMesh())
968  // If the parent element will not be re-cut, then all of its nodes must have the same
969  // CutSubdomainID. Therefore, just query the first node in this parent element to
970  // get its CutSubdomainID.
971  for (auto e : healed_elems)
972  if (elem1->processor_id() == _mesh->processor_id() &&
973  e->processor_id() == _mesh->processor_id())
974  {
975  storeMaterialPropertiesForElement(/*parent_elem = */ elem1, /*child_elem = */ e);
976  // In case the healed element is not re-cut, copy the corresponding material
977  // properties to the parent element now than later.
978  CutSubdomainID parent_gcsid =
979  _geometric_cuts[i]->getCutSubdomainID(elem1->node_ptr(0));
980  CutSubdomainID gcsid = _geom_cut_elems[e]._cut_subdomain_id;
981  if (parent_gcsid == gcsid)
983  }
984 
985  if (_displaced_mesh)
986  {
987  Elem * elem1_displaced = _displaced_mesh->elem_ptr(it.first->id());
988  Elem * elem2_displaced = _displaced_mesh->elem_ptr(it.second->id());
989 
990  std::map<unique_id_type, XFEMCutElem *>::iterator cemit =
991  _cut_elem_map.find(elem1_displaced->unique_id());
992  if (cemit != _cut_elem_map.end())
993  {
994  const XFEMCutElem * xfce = cemit->second;
995 
996  for (unsigned int in = 0; in < elem1_displaced->n_nodes(); ++in)
997  {
998  Node * e1node_displaced = elem1_displaced->node_ptr(in);
999  Node * e2node_displaced = elem2_displaced->node_ptr(in);
1000  if (!xfce->isPointPhysical(*elem1->node_ptr(in)) &&
1001  e1node_displaced != e2node_displaced)
1002  {
1003  elem1_displaced->set_node(in) = e2node_displaced;
1004  nodes_to_delete_displaced.insert(e1node_displaced);
1005  }
1006  else if (e1node_displaced != e2node_displaced)
1007  nodes_to_delete_displaced.insert(e2node_displaced);
1008  }
1009  }
1010  else
1011  mooseError("Could not find XFEMCutElem for element to be kept in healing");
1012 
1013  elem2_displaced->nullify_neighbors();
1014  _displaced_mesh->get_boundary_info().remove(elem2_displaced);
1015  _displaced_mesh->delete_elem(elem2_displaced);
1016  }
1017 
1018  // remove the property storage of deleted element/side
1019  _material_data[0]->eraseProperty(elem2);
1020  _bnd_material_data[0]->eraseProperty(elem2);
1021 
1022  cutelems_to_delete.insert(elem2->unique_id());
1023  elem2->nullify_neighbors();
1024  _mesh->get_boundary_info().remove(elem2);
1025  unsigned int deleted_elem_id = elem2->id();
1026  _mesh->delete_elem(elem2);
1027  if (_debug_output_level > 1)
1028  {
1029  if (deleted_elem_count == 0)
1030  _console << "\n";
1031  _console << "XFEM healing deleted element: " << deleted_elem_id << std::endl;
1032  }
1033  ++deleted_elem_count;
1034  mesh_changed = true;
1035  }
1036  }
1037  }
1038 
1039  for (auto & sit : nodes_to_delete)
1040  {
1041  Node * node_to_delete = sit;
1042  dof_id_type deleted_node_id = node_to_delete->id();
1043  _mesh->get_boundary_info().remove(node_to_delete);
1044  _mesh->delete_node(node_to_delete);
1045  if (_debug_output_level > 1)
1046  _console << "XFEM healing deleted node: " << deleted_node_id << std::endl;
1047  }
1048 
1049  if (_displaced_mesh)
1050  {
1051  for (auto & sit : nodes_to_delete_displaced)
1052  {
1053  Node * node_to_delete_displaced = sit;
1054  _displaced_mesh->get_boundary_info().remove(node_to_delete_displaced);
1055  _displaced_mesh->delete_node(node_to_delete_displaced);
1056  }
1057  }
1058 
1059  for (auto & ced : cutelems_to_delete)
1060  if (_cut_elem_map.find(ced) != _cut_elem_map.end())
1061  {
1062  delete _cut_elem_map.find(ced)->second;
1063  _cut_elem_map.erase(ced);
1064  }
1065 
1066  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1067  if (_geometric_cuts[i]->shouldHealMesh())
1068  _sibling_elems[_geometric_cuts[i]->getInterfaceID()].clear();
1069 
1070  if (_displaced_mesh)
1071  {
1072  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1073  if (_geometric_cuts[i]->shouldHealMesh())
1074  _sibling_displaced_elems[_geometric_cuts[i]->getInterfaceID()].clear();
1075  }
1076 
1077  for (auto & ceh : _crack_tip_elems_to_be_healed)
1078  {
1079  _crack_tip_elems.erase(ceh);
1081  delete _cut_elem_map.find(ceh->unique_id())->second;
1082  _cut_elem_map.erase(ceh->unique_id());
1083  }
1084 
1085  if (!healed_geometric_cuts.empty() && _debug_output_level > 0)
1086  {
1087  _console << "\nXFEM mesh healing complete\n";
1088  _console << "Names of healed geometric cut objects: ";
1089  for (auto geomcut : healed_geometric_cuts)
1090  _console << geomcut << " ";
1091  _console << "\n";
1092  _console << "# deleted nodes: " << nodes_to_delete.size() << "\n";
1093  _console << "# deleted elements: " << deleted_elem_count << "\n";
1094  _console << std::flush;
1095  }
1096 
1097  return mesh_changed;
1098 }
unsigned int _debug_output_level
Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mes...
Definition: XFEM.h:362
void storeMaterialPropertiesForElement(const Elem *parent_elem, const Elem *child_elem)
Helper function to store the material properties of a healed element.
Definition: XFEM.C:2185
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:327
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:332
void mooseError(Args &&... args)
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:333
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:330
std::set< const Elem * > _crack_tip_elems_to_be_healed
Definition: XFEM.h:331
unsigned int CutSubdomainID
Definition: XFEMAppTypes.h:18
std::vector< MaterialData *> _bnd_material_data
void loadMaterialPropertiesForElement(const Elem *elem, const Elem *elem_from, std::unordered_map< const Elem *, Xfem::CutElemInfo > &cached_cei) const
Helper function to store the material properties of a healed element.
Definition: XFEM.C:2237
std::vector< MaterialData *> _material_data
const std::string name
Definition: Setup.h:20
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335
MeshBase * _mesh
MeshBase * _displaced_mesh
const ConsoleStream _console
uint8_t dof_id_type
std::unordered_map< const Elem *, Xfem::CutElemInfo > _geom_cut_elems
All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK...
Definition: XFEM.h:386
bool isPointPhysical(const Point &p) const
Definition: XFEMCutElem.C:193

◆ initCutIntersectionEdge()

bool XFEM::initCutIntersectionEdge ( Point  cut_origin,
RealVectorValue  cut_normal,
Point &  edge_p1,
Point &  edge_p2,
Real dist 
)

Definition at line 898 of file XFEM.C.

Referenced by correctCrackExtensionDirection(), and markCutEdgesByState().

900 {
901  dist = 0.0;
902  bool does_intersect = false;
903  Point origin2p1 = edge_p1 - cut_origin;
904  Real plane2p1 = cut_normal(0) * origin2p1(0) + cut_normal(1) * origin2p1(1);
905  Point origin2p2 = edge_p2 - cut_origin;
906  Real plane2p2 = cut_normal(0) * origin2p2(0) + cut_normal(1) * origin2p2(1);
907 
908  if (plane2p1 * plane2p2 < 0.0)
909  {
910  dist = -plane2p1 / (plane2p2 - plane2p1);
911  does_intersect = true;
912  }
913  return does_intersect;
914 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ initSolution()

void XFEM::initSolution ( const std::vector< std::shared_ptr< NonlinearSystemBase >> &  nl,
AuxiliarySystem aux 
)
overridevirtual

Implements XFEMInterface.

Definition at line 305 of file XFEM.C.

307 {
308  if (nls.size() != 1)
309  mooseError("XFEM does not currently support multiple nonlinear systems");
310 
311  nls[0]->serializeSolution();
312  aux.serializeSolution();
313  NumericVector<Number> & current_solution = *nls[0]->system().current_local_solution;
314  NumericVector<Number> & old_solution = nls[0]->solutionOld();
315  NumericVector<Number> & older_solution = nls[0]->solutionOlder();
316  NumericVector<Number> & current_aux_solution = *aux.system().current_local_solution;
317  NumericVector<Number> & old_aux_solution = aux.solutionOld();
318  NumericVector<Number> & older_aux_solution = aux.solutionOlder();
319 
320  setSolution(*nls[0], _cached_solution, current_solution, old_solution, older_solution);
321  setSolution(
322  aux, _cached_aux_solution, current_aux_solution, old_aux_solution, older_aux_solution);
323 
324  current_solution.close();
325  old_solution.close();
326  older_solution.close();
327  current_aux_solution.close();
328  old_aux_solution.close();
329  older_aux_solution.close();
330 
331  _cached_solution.clear();
332  _cached_aux_solution.clear();
333 }
void setSolution(SystemBase &sys, const std::map< unique_id_type, std::vector< Real >> &stored_solution, NumericVector< Number > &current_solution, NumericVector< Number > &old_solution, NumericVector< Number > &older_solution)
Set the solution for all locally-owned nodes/elements that have stored values.
Definition: XFEM.C:2036
void mooseError(Args &&... args)
NumericVector< Number > & solutionOlder()
std::map< unique_id_type, std::vector< Real > > _cached_aux_solution
Data structure to store the auxiliary solution for nodes/elements affected by XFEM For each node/elem...
Definition: XFEM.h:380
virtual void close()=0
virtual System & system() override
std::map< unique_id_type, std::vector< Real > > _cached_solution
Data structure to store the nonlinear solution for nodes/elements affected by XFEM For each node/elem...
Definition: XFEM.h:371
virtual void serializeSolution()
NumericVector< Number > & solutionOld()

◆ isElemAtCrackTip()

bool XFEM::isElemAtCrackTip ( const Elem *  elem) const

Definition at line 1700 of file XFEM.C.

Referenced by markCutEdgesByGeometry(), and markCutEdgesByState().

1701 {
1702  return (_crack_tip_elems.find(elem) != _crack_tip_elems.end());
1703 }
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:330

◆ isElemCut() [1/2]

bool XFEM::isElemCut ( const Elem *  elem,
XFEMCutElem *&  xfce 
) const

Definition at line 1706 of file XFEM.C.

Referenced by getXFEMFaceWeights(), getXFEMWeights(), and isElemCut().

1707 {
1708  const auto it = _cut_elem_map.find(elem->unique_id());
1709  if (it != _cut_elem_map.end())
1710  {
1711  xfce = it->second;
1712  const EFAElement * EFAelem = xfce->getEFAElement();
1713  if (EFAelem->isPartial()) // exclude the full crack tip elements
1714  return true;
1715  }
1716 
1717  xfce = nullptr;
1718  return false;
1719 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
virtual bool isPartial() const =0
virtual const EFAElement * getEFAElement() const =0

◆ isElemCut() [2/2]

bool XFEM::isElemCut ( const Elem *  elem) const

Definition at line 1722 of file XFEM.C.

1723 {
1724  XFEMCutElem * xfce;
1725  return isElemCut(elem, xfce);
1726 }
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1706

◆ isPointInsidePhysicalDomain()

bool XFEM::isPointInsidePhysicalDomain ( const Elem *  elem,
const Point &  point 
) const

Return true if the point is inside the element physical domain Note: if this element is not cut, return true too.

Definition at line 1648 of file XFEM.C.

Referenced by pickFirstPhysicalNode().

1649 {
1650  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1651  it = _cut_elem_map.find(elem->unique_id());
1652  if (it != _cut_elem_map.end())
1653  {
1654  XFEMCutElem * xfce = it->second;
1655 
1656  if (xfce->isPointPhysical(point))
1657  return true;
1658  }
1659  else
1660  return true;
1661 
1662  return false;
1663 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
bool isPointPhysical(const Point &p) const
Definition: XFEMCutElem.C:193

◆ loadMaterialPropertiesForElement()

void XFEM::loadMaterialPropertiesForElement ( const Elem *  elem,
const Elem *  elem_from,
std::unordered_map< const Elem *, Xfem::CutElemInfo > &  cached_cei 
) const
private

Helper function to store the material properties of a healed element.

Parameters
elemThe cut element to restore material properties to.
elem_fromThe element to copy material properties from.
cached_ceiThe material properties cache to use.

Definition at line 2237 of file XFEM.C.

Referenced by cutMeshWithEFA(), and healMesh().

2241 {
2242  // Restore the element material properties
2243  mooseAssert(cached_cei.count(elem_from) > 0, "XFEM: Unable to find cached material properties.");
2244  Xfem::CutElemInfo & cei = cached_cei[elem_from];
2245 
2246  // Load element material properties from cached properties
2248  cei._elem_material_properties,
2249  _material_data[0]->getMaterialPropertyStorageForXFEM({}));
2250 
2251  // Check if any of the element side need material properties
2252  bool need_boundary_materials = false;
2253  for (unsigned int side = 0; side < elem->n_sides(); ++side)
2254  {
2255  std::vector<boundary_id_type> elem_boundary_ids;
2256  _mesh->get_boundary_info().boundary_ids(elem, side, elem_boundary_ids);
2257  for (auto bdid : elem_boundary_ids)
2259  need_boundary_materials = true;
2260  }
2261 
2262  // Load boundary material properties from cached properties
2263  if (need_boundary_materials)
2265  elem,
2266  cei._bnd_material_properties,
2267  _bnd_material_data[0]->getMaterialPropertyStorageForXFEM({}));
2268 }
bool needBoundaryMaterialOnSide(BoundaryID bnd_id, const THREAD_ID tid)
FEProblemBase * _fe_problem
std::vector< MaterialData *> _bnd_material_data
std::vector< MaterialData *> _material_data
void loadMaterialPropertiesForElementHelper(const Elem *elem, const Xfem::CachedMaterialProperties &cached_props, MaterialPropertyStorage &storage) const
Load the material properties.
Definition: XFEM.C:2213
MeshBase * _mesh
Information about a cut element.
Definition: XFEM.h:57

◆ loadMaterialPropertiesForElementHelper()

void XFEM::loadMaterialPropertiesForElementHelper ( const Elem *  elem,
const Xfem::CachedMaterialProperties cached_props,
MaterialPropertyStorage storage 
) const
private

Load the material properties.

Parameters
props_deserializedThe material properties
props_serializedThe serialized material properties

This does very dirty things and writes back to MOOSE's stateful properties. It should not do this in the future.

Definition at line 2213 of file XFEM.C.

Referenced by loadMaterialPropertiesForElement().

2216 {
2217  if (!storage.hasStatefulProperties())
2218  return;
2219 
2220  for (const auto state : storage.statefulIndexRange())
2221  {
2222  const auto & serialized_props = cached_props[state - 1];
2223  for (const auto & [side, serialized_side_props] : serialized_props)
2224  {
2225  std::istringstream iss;
2226  iss.str(serialized_side_props);
2227  iss.clear();
2228 
2229  // This is very dirty. We should not write to MOOSE's stateful properties.
2230  // Please remove me :(
2231  dataLoad(iss, storage.setProps(elem, side, state), nullptr);
2232  }
2233  }
2234 }
MaterialProperties & setProps(const Elem *elem, unsigned int side, const unsigned int state=0)
IntRange< unsigned int > statefulIndexRange() const
bool hasStatefulProperties() const
void dataLoad(std::istream &stream, FeatureFloodCount::FeatureData &feature, void *context)

◆ markCutEdgesByGeometry()

bool XFEM::markCutEdgesByGeometry ( )

Definition at line 392 of file XFEM.C.

Referenced by markCuts().

393 {
394  bool marked_edges = false;
395  bool marked_nodes = false;
396 
397  for (const auto & gme : _geom_marked_elems_2d)
398  {
399  for (const auto & gmei : gme.second)
400  {
401  EFAElement2D * EFAElem = getEFAElem2D(gme.first);
402 
403  for (unsigned int i = 0; i < gmei._elem_cut_edges.size(); ++i) // mark element edges
404  {
405  if (!EFAElem->isEdgePhantom(
406  gmei._elem_cut_edges[i]._host_side_id)) // must not be phantom edge
407  {
408  _efa_mesh.addElemEdgeIntersection(gme.first->id(),
409  gmei._elem_cut_edges[i]._host_side_id,
410  gmei._elem_cut_edges[i]._distance);
411  marked_edges = true;
412  }
413  }
414 
415  for (unsigned int i = 0; i < gmei._elem_cut_nodes.size(); ++i) // mark element edges
416  {
417  _efa_mesh.addElemNodeIntersection(gme.first->id(), gmei._elem_cut_nodes[i]._host_id);
418  marked_nodes = true;
419  }
420 
421  for (unsigned int i = 0; i < gmei._frag_cut_edges.size();
422  ++i) // MUST DO THIS AFTER MARKING ELEMENT EDGES
423  {
424  if (!EFAElem->getFragment(0)->isSecondaryInteriorEdge(
425  gmei._frag_cut_edges[i]._host_side_id))
426  {
427  if (_efa_mesh.addFragEdgeIntersection(gme.first->id(),
428  gmei._frag_cut_edges[i]._host_side_id,
429  gmei._frag_cut_edges[i]._distance))
430  {
431  marked_edges = true;
432  if (!isElemAtCrackTip(gme.first))
433  _has_secondary_cut = true;
434  }
435  }
436  }
437  }
438  }
439 
440  return marked_edges || marked_nodes;
441 }
void addElemNodeIntersection(unsigned int elemid, unsigned int nodeid)
void addElemEdgeIntersection(unsigned int elemid, unsigned int edgeid, double position)
bool isSecondaryInteriorEdge(unsigned int edge_id) const
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
bool isEdgePhantom(unsigned int edge_id) const
bool _has_secondary_cut
Definition: XFEM.h:320
EFAFragment2D * getFragment(unsigned int frag_id) const
bool addFragEdgeIntersection(unsigned int elemid, unsigned int frag_edge_id, double position)
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo2D > > _geom_marked_elems_2d
Data structure for storing information about all 2D elements to be cut by geometry.
Definition: XFEM.h:344
bool isElemAtCrackTip(const Elem *elem) const
Definition: XFEM.C:1700
EFAElement2D * getEFAElem2D(const Elem *elem)
Get the EFAElement2D object for a specified libMesh element.
Definition: XFEM.C:1746

◆ markCutEdgesByState()

bool XFEM::markCutEdgesByState ( Real  time)

Definition at line 585 of file XFEM.C.

Referenced by markCuts().

586 {
587  bool marked_edges = false;
588  for (std::map<const Elem *, RealVectorValue>::iterator pmeit = _state_marked_elems.begin();
589  pmeit != _state_marked_elems.end();
590  ++pmeit)
591  {
592  const Elem * elem = pmeit->first;
593  RealVectorValue normal = pmeit->second;
594  EFAElement2D * CEMElem = getEFAElem2D(elem);
595 
596  Real volfrac_elem = getPhysicalVolumeFraction(elem);
597  if (volfrac_elem < 0.25)
598  continue;
599 
600  // continue if elem is already cut twice - IMPORTANT
601  if (CEMElem->isFinalCut())
602  continue;
603 
604  // find the first cut edge
605  unsigned int nsides = CEMElem->numEdges();
606  unsigned int orig_cut_side_id = std::numeric_limits<unsigned int>::max();
607  Real orig_cut_distance = -1.0;
608  EFANode * orig_node = nullptr;
609  EFAEdge * orig_edge = nullptr;
610 
611  // crack tip origin coordinates and direction
612  Point crack_tip_origin(0, 0, 0);
613  Point crack_tip_direction(0, 0, 0);
614 
615  if (isElemAtCrackTip(elem)) // crack tip element's crack intiation
616  {
617  orig_cut_side_id = CEMElem->getTipEdgeID();
618  if (orig_cut_side_id < nsides) // valid crack-tip edge found
619  {
620  orig_edge = CEMElem->getEdge(orig_cut_side_id);
621  orig_node = CEMElem->getTipEmbeddedNode();
622  }
623  else
624  mooseError("element ", elem->id(), " has no valid crack-tip edge");
625 
626  // obtain the crack tip origin coordinates and direction.
627  std::map<const Elem *, std::vector<Point>>::iterator ecodm =
629  if (ecodm != _elem_crack_origin_direction_map.end())
630  {
631  crack_tip_origin = (ecodm->second)[0];
632  crack_tip_direction = (ecodm->second)[1];
633  }
634  else
635  mooseError("element ", elem->id(), " cannot find its crack tip origin and direction.");
636  }
637  else
638  {
639  std::map<const Elem *, unsigned int>::iterator mit1;
640  mit1 = _state_marked_elem_sides.find(elem);
641  std::set<const Elem *>::iterator mit2;
642  mit2 = _state_marked_frags.find(elem);
643 
644  if (mit1 != _state_marked_elem_sides.end()) // specified boundary crack initiation
645  {
646  orig_cut_side_id = mit1->second;
647  if (!CEMElem->isEdgePhantom(orig_cut_side_id) &&
648  !CEMElem->getEdge(orig_cut_side_id)->hasIntersection())
649  {
650  orig_cut_distance = 0.5;
651  _efa_mesh.addElemEdgeIntersection(elem->id(), orig_cut_side_id, orig_cut_distance);
652  orig_edge = CEMElem->getEdge(orig_cut_side_id);
653  orig_node = orig_edge->getEmbeddedNode(0);
654  // get a virtual crack tip direction
655  Point elem_center(0.0, 0.0, 0.0);
656  Point edge_center;
657  for (unsigned int i = 0; i < nsides; ++i)
658  {
659  elem_center += getEFANodeCoords(CEMElem->getEdge(i)->getNode(0), CEMElem, elem);
660  elem_center += getEFANodeCoords(CEMElem->getEdge(i)->getNode(1), CEMElem, elem);
661  }
662  elem_center /= nsides * 2.0;
663  edge_center = getEFANodeCoords(orig_edge->getNode(0), CEMElem, elem) +
664  getEFANodeCoords(orig_edge->getNode(1), CEMElem, elem);
665  edge_center /= 2.0;
666  crack_tip_origin = edge_center;
667  crack_tip_direction = elem_center - edge_center;
668  crack_tip_direction /= pow(crack_tip_direction.norm_sq(), 0.5);
669  }
670  else
671  continue; // skip this elem if specified boundary edge is phantom
672  }
673  else if (mit2 != _state_marked_frags.end()) // cut-surface secondary crack initiation
674  {
675  if (CEMElem->numFragments() != 1)
676  mooseError("element ",
677  elem->id(),
678  " flagged for a secondary crack, but has ",
679  CEMElem->numFragments(),
680  " fragments");
681  std::vector<unsigned int> interior_edge_id = CEMElem->getFragment(0)->getInteriorEdgeID();
682  if (interior_edge_id.size() == 1)
683  orig_cut_side_id = interior_edge_id[0];
684  else
685  continue; // skip this elem if more than one interior edges found (i.e. elem's been cut
686  // twice)
687  orig_cut_distance = 0.5;
688  _efa_mesh.addFragEdgeIntersection(elem->id(), orig_cut_side_id, orig_cut_distance);
689  orig_edge = CEMElem->getFragmentEdge(0, orig_cut_side_id);
690  orig_node = orig_edge->getEmbeddedNode(0); // must be an interior embedded node
691  Point elem_center(0.0, 0.0, 0.0);
692  Point edge_center;
693  unsigned int nsides_frag = CEMElem->getFragment(0)->numEdges();
694  for (unsigned int i = 0; i < nsides_frag; ++i)
695  {
696  elem_center +=
697  getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(0), CEMElem, elem);
698  elem_center +=
699  getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(1), CEMElem, elem);
700  }
701  elem_center /= nsides_frag * 2.0;
702  edge_center = getEFANodeCoords(orig_edge->getNode(0), CEMElem, elem) +
703  getEFANodeCoords(orig_edge->getNode(1), CEMElem, elem);
704  edge_center /= 2.0;
705  crack_tip_origin = edge_center;
706  crack_tip_direction = elem_center - edge_center;
707  crack_tip_direction /= pow(crack_tip_direction.norm_sq(), 0.5);
708  }
709  else
710  mooseError("element ",
711  elem->id(),
712  " flagged for state-based growth, but has no edge intersections");
713  }
714 
715  Point cut_origin(0.0, 0.0, 0.0);
716  if (orig_node)
717  cut_origin = getEFANodeCoords(orig_node, CEMElem, elem); // cutting plane origin's coords
718  else
719  mooseError("element ", elem->id(), " does not have valid orig_node");
720 
721  // loop through element edges to add possible second cut points
722  std::vector<Point> edge_ends(2, Point(0.0, 0.0, 0.0));
723  Point edge1(0.0, 0.0, 0.0);
724  Point edge2(0.0, 0.0, 0.0);
725  Point cut_edge_point(0.0, 0.0, 0.0);
726  bool find_compatible_direction = false;
727  unsigned int edge_id_keep = 0;
728  Real distance_keep = 0.0;
729  Point normal_keep(0.0, 0.0, 0.0);
730  Real distance = 0.0;
731  bool edge_cut = false;
732 
733  for (unsigned int i = 0; i < nsides; ++i)
734  {
735  if (!orig_edge->isPartialOverlap(*CEMElem->getEdge(i)))
736  {
737  edge_ends[0] = getEFANodeCoords(CEMElem->getEdge(i)->getNode(0), CEMElem, elem);
738  edge_ends[1] = getEFANodeCoords(CEMElem->getEdge(i)->getNode(1), CEMElem, elem);
740  crack_tip_origin, normal, edge_ends[0], edge_ends[1], distance) &&
741  (!CEMElem->isEdgePhantom(i))))
742  {
743  cut_edge_point = distance * edge_ends[1] + (1.0 - distance) * edge_ends[0];
744  distance_keep = distance;
745  edge_id_keep = i;
746  normal_keep = normal;
747  edge_cut = true;
748  break;
749  }
750  }
751  }
752 
753  Point between_two_cuts = (cut_edge_point - crack_tip_origin);
754  between_two_cuts /= pow(between_two_cuts.norm_sq(), 0.5);
755  Real angle_between_two_cuts = between_two_cuts * crack_tip_direction;
756 
757  if (angle_between_two_cuts > std::cos(45.0 / 180.0 * 3.14159)) // original cut direction is good
758  find_compatible_direction = true;
759 
760  if (!find_compatible_direction && edge_cut)
762  CEMElem,
763  orig_edge,
764  normal,
765  crack_tip_origin,
766  crack_tip_direction,
767  distance_keep,
768  edge_id_keep,
769  normal_keep);
770 
771  if (edge_cut)
772  {
774  _efa_mesh.addElemEdgeIntersection(elem->id(), edge_id_keep, distance_keep);
775  else
776  {
777  Point growth_direction(0.0, 0.0, 0.0);
778 
779  growth_direction(0) = -normal_keep(1);
780  growth_direction(1) = normal_keep(0);
781 
782  if (growth_direction * crack_tip_direction < 1.0e-10)
783  growth_direction *= (-1.0);
784 
785  Real x0 = crack_tip_origin(0);
786  Real y0 = crack_tip_origin(1);
787  Real x1 = x0 + _crack_growth_increment * growth_direction(0);
788  Real y1 = y0 + _crack_growth_increment * growth_direction(1);
789 
790  XFEMCrackGrowthIncrement2DCut geometric_cut(x0, y0, x1, y1, time * 0.9, time * 0.9);
791 
792  for (const auto & elem : _mesh->element_ptr_range())
793  {
794  std::vector<CutEdgeForCrackGrowthIncr> elem_cut_edges;
795  EFAElement2D * CEMElem = getEFAElem2D(elem);
796 
797  // continue if elem has been already cut twice - IMPORTANT
798  if (CEMElem->isFinalCut())
799  continue;
800 
801  // mark cut edges for the element and its fragment
802  geometric_cut.cutElementByCrackGrowthIncrement(elem, elem_cut_edges, time);
803 
804  for (unsigned int i = 0; i < elem_cut_edges.size(); ++i) // mark element edges
805  {
806  if (!CEMElem->isEdgePhantom(
807  elem_cut_edges[i]._host_side_id)) // must not be phantom edge
808  {
810  elem->id(), elem_cut_edges[i]._host_side_id, elem_cut_edges[i]._distance);
811  }
812  }
813  }
814  }
815  }
816  // loop though framgent boundary edges to add possible second cut points
817  // N.B. must do this after marking element edges
818  if (CEMElem->numFragments() > 0 && !edge_cut)
819  {
820  for (unsigned int i = 0; i < CEMElem->getFragment(0)->numEdges(); ++i)
821  {
822  if (!orig_edge->isPartialOverlap(*CEMElem->getFragmentEdge(0, i)))
823  {
824  edge_ends[0] =
825  getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(0), CEMElem, elem);
826  edge_ends[1] =
827  getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(1), CEMElem, elem);
829  crack_tip_origin, normal, edge_ends[0], edge_ends[1], distance) &&
830  (!CEMElem->getFragment(0)->isSecondaryInteriorEdge(i)))
831  {
832  if (_efa_mesh.addFragEdgeIntersection(elem->id(), edge_id_keep, distance_keep))
833  if (!isElemAtCrackTip(elem))
834  _has_secondary_cut = true;
835  break;
836  }
837  }
838  }
839  }
840 
841  marked_edges = true;
842 
843  } // loop over all state_marked_elems
844 
845  return marked_edges;
846 }
void correctCrackExtensionDirection(const Elem *elem, EFAElement2D *CEMElem, EFAEdge *orig_edge, Point normal, Point crack_tip_origin, Point crack_tip_direction, Real &distance_keep, unsigned int &edge_id_keep, Point &normal_keep)
Definition: XFEM.C:444
bool _use_crack_growth_increment
Definition: XFEM.h:324
EFANode * getEmbeddedNode(unsigned int index) const
Definition: EFAEdge.C:332
void addElemEdgeIntersection(unsigned int elemid, unsigned int edgeid, double position)
Real _crack_growth_increment
Definition: XFEM.h:325
bool isSecondaryInteriorEdge(unsigned int edge_id) const
void mooseError(Args &&... args)
bool isPartialOverlap(const EFAEdge &other) const
Definition: EFAEdge.C:59
unsigned int getTipEdgeID() const
unsigned int numEdges() const
EFAEdge * getEdge(unsigned int edge_id) const
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
bool hasIntersection() const
Definition: EFAEdge.C:200
std::vector< unsigned int > getInteriorEdgeID() const
unsigned int numEdges() const
bool isEdgePhantom(unsigned int edge_id) const
virtual unsigned int numFragments() const
Definition: EFAElement2D.C:207
Real distance(const Point &p)
bool _has_secondary_cut
Definition: XFEM.h:320
bool initCutIntersectionEdge(Point cut_origin, RealVectorValue cut_normal, Point &edge_p1, Point &edge_p2, Real &dist)
Definition: XFEM.C:898
EFAFragment2D * getFragment(unsigned int frag_id) const
bool addFragEdgeIntersection(unsigned int elemid, unsigned int frag_edge_id, double position)
virtual bool isFinalCut() const
Definition: EFAElement2D.C:796
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
std::set< const Elem * > _state_marked_frags
Definition: XFEM.h:340
EFANode * getTipEmbeddedNode() const
MeshBase * _mesh
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool isElemAtCrackTip(const Elem *elem) const
Definition: XFEM.C:1700
Real getPhysicalVolumeFraction(const Elem *elem) const
Get the volume fraction of an element that is physical.
Definition: XFEM.C:1628
std::map< const Elem *, unsigned int > _state_marked_elem_sides
Definition: XFEM.h:341
EFANode * getNode(unsigned int index) const
Definition: EFAEdge.C:181
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=nullptr) const
Definition: XFEM.C:1596
std::map< const Elem *, RealVectorValue > _state_marked_elems
Definition: XFEM.h:339
EFAElement2D * getEFAElem2D(const Elem *elem)
Get the EFAElement2D object for a specified libMesh element.
Definition: XFEM.C:1746
EFAEdge * getFragmentEdge(unsigned int frag_id, unsigned int edge_id) const

◆ markCutFacesByGeometry()

bool XFEM::markCutFacesByGeometry ( )

Definition at line 849 of file XFEM.C.

Referenced by markCuts().

850 {
851  bool marked_faces = false;
852 
853  for (const auto & gme : _geom_marked_elems_3d)
854  {
855  for (const auto & gmei : gme.second)
856  {
857  EFAElement3D * EFAElem = getEFAElem3D(gme.first);
858 
859  for (unsigned int i = 0; i < gmei._elem_cut_faces.size(); ++i) // mark element faces
860  {
861  if (!EFAElem->isFacePhantom(gmei._elem_cut_faces[i]._face_id)) // must not be phantom face
862  {
863  _efa_mesh.addElemFaceIntersection(gme.first->id(),
864  gmei._elem_cut_faces[i]._face_id,
865  gmei._elem_cut_faces[i]._face_edge,
866  gmei._elem_cut_faces[i]._position);
867  marked_faces = true;
868  }
869  }
870 
871  for (unsigned int i = 0; i < gmei._frag_cut_faces.size();
872  ++i) // MUST DO THIS AFTER MARKING ELEMENT EDGES
873  {
874  if (!EFAElem->getFragment(0)->isThirdInteriorFace(gmei._frag_cut_faces[i]._face_id))
875  {
876  _efa_mesh.addFragFaceIntersection(gme.first->id(),
877  gmei._frag_cut_faces[i]._face_id,
878  gmei._frag_cut_faces[i]._face_edge,
879  gmei._frag_cut_faces[i]._position);
880  marked_faces = true;
881  }
882  }
883  }
884  }
885 
886  return marked_faces;
887 }
EFAFragment3D * getFragment(unsigned int frag_id) const
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo3D > > _geom_marked_elems_3d
Data structure for storing information about all 3D elements to be cut by geometry.
Definition: XFEM.h:347
bool isFacePhantom(unsigned int face_id) const
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
void addElemFaceIntersection(unsigned int elemid, unsigned int faceid, const std::vector< unsigned int > &edgeid, const std::vector< double > &position)
bool isThirdInteriorFace(unsigned int face_id) const
void addFragFaceIntersection(unsigned int ElemID, unsigned int FragFaceID, const std::vector< unsigned int > &FragFaceEdgeID, const std::vector< double > &position)
EFAElement3D * getEFAElem3D(const Elem *elem)
Get the EFAElement3D object for a specified libMesh element.
Definition: XFEM.C:1758

◆ markCutFacesByState()

bool XFEM::markCutFacesByState ( )

Definition at line 890 of file XFEM.C.

Referenced by markCuts().

891 {
892  bool marked_faces = false;
893  // TODO: need to finish this for 3D problems
894  return marked_faces;
895 }

◆ markCuts()

bool XFEM::markCuts ( Real  time)

Definition at line 375 of file XFEM.C.

Referenced by update().

376 {
377  bool marked_sides = false;
378  if (_mesh->mesh_dimension() == 2)
379  {
380  marked_sides = markCutEdgesByGeometry();
381  marked_sides |= markCutEdgesByState(time);
382  }
383  else if (_mesh->mesh_dimension() == 3)
384  {
385  marked_sides = markCutFacesByGeometry();
386  marked_sides |= markCutFacesByState();
387  }
388  return marked_sides;
389 }
bool markCutEdgesByState(Real time)
Definition: XFEM.C:585
bool markCutEdgesByGeometry()
Definition: XFEM.C:392
MeshBase * _mesh
bool markCutFacesByGeometry()
Definition: XFEM.C:849
bool markCutFacesByState()
Definition: XFEM.C:890

◆ pickFirstPhysicalNode()

const Node * XFEM::pickFirstPhysicalNode ( const Elem *  e,
const Elem *  e0 
) const
private

Return the first node in the provided element that is found to be in the physical domain.

Parameters
eConstant pointer to the child element
e0Constant pointer to the parent element whose nodes will be querried
Returns
A constant pointer to the node

Definition at line 2284 of file XFEM.C.

Referenced by getCutSubdomainID().

2285 {
2286  for (auto i : e0->node_index_range())
2287  if (isPointInsidePhysicalDomain(e, e0->node_ref(i)))
2288  return e0->node_ptr(i);
2289  mooseError("cannot find a physical node in the current element");
2290  return nullptr;
2291 }
void mooseError(Args &&... args)
bool isPointInsidePhysicalDomain(const Elem *elem, const Point &point) const
Return true if the point is inside the element physical domain Note: if this element is not cut...
Definition: XFEM.C:1648

◆ setCrackGrowthMethod()

void XFEM::setCrackGrowthMethod ( bool  use_crack_growth_increment,
Real  crack_growth_increment 
)

Definition at line 1830 of file XFEM.C.

1831 {
1832  _use_crack_growth_increment = use_crack_growth_increment;
1833  _crack_growth_increment = crack_growth_increment;
1834 }
bool _use_crack_growth_increment
Definition: XFEM.h:324
Real _crack_growth_increment
Definition: XFEM.h:325

◆ setDebugOutputLevel()

void XFEM::setDebugOutputLevel ( unsigned int  debug_output_level)

Controls amount of debugging information output.

Parameters
debug_output_levelHow much information to output (see description of _debug_output_level)

Definition at line 1837 of file XFEM.C.

1838 {
1839  _debug_output_level = debug_output_level;
1840 }
unsigned int _debug_output_level
Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mes...
Definition: XFEM.h:362

◆ setSolution()

void XFEM::setSolution ( SystemBase sys,
const std::map< unique_id_type, std::vector< Real >> &  stored_solution,
NumericVector< Number > &  current_solution,
NumericVector< Number > &  old_solution,
NumericVector< Number > &  older_solution 
)
private

Set the solution for all locally-owned nodes/elements that have stored values.

Parameters
sysSystem for which the solution is set
stored_solutionData structure that the stored solution is obtained from
current_solutionCurrent solution vector that will be set
old_solutionOld solution vector that will be set
older_solutionOlder solution vector that will be set

Definition at line 2036 of file XFEM.C.

Referenced by initSolution().

2041 {
2042  for (auto & node : _mesh->local_node_ptr_range())
2043  {
2044  auto mit = stored_solution.find(node->unique_id());
2045  if (mit != stored_solution.end())
2046  {
2047  const std::vector<Real> & stored_node_solution = mit->second;
2048  std::vector<dof_id_type> stored_solution_dofs = getNodeSolutionDofs(node, sys);
2049  setSolutionForDOFs(stored_node_solution,
2050  stored_solution_dofs,
2051  current_solution,
2052  old_solution,
2053  older_solution);
2054  }
2055  }
2056 
2057  for (auto & elem : as_range(_mesh->local_elements_begin(), _mesh->local_elements_end()))
2058  {
2059  auto mit = stored_solution.find(elem->unique_id());
2060  if (mit != stored_solution.end())
2061  {
2062  const std::vector<Real> & stored_elem_solution = mit->second;
2063  std::vector<dof_id_type> stored_solution_dofs = getElementSolutionDofs(elem, sys);
2064  setSolutionForDOFs(stored_elem_solution,
2065  stored_solution_dofs,
2066  current_solution,
2067  old_solution,
2068  older_solution);
2069  }
2070  }
2071 }
void setSolutionForDOFs(const std::vector< Real > &stored_solution, const std::vector< dof_id_type > &stored_solution_dofs, NumericVector< Number > &current_solution, NumericVector< Number > &old_solution, NumericVector< Number > &older_solution)
Set the solution for a set of DOFs.
Definition: XFEM.C:2074
std::vector< dof_id_type > getNodeSolutionDofs(const Node *node, SystemBase &sys) const
Get a vector of the dof indices for all components of all variables associated with a node...
Definition: XFEM.C:2123
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
MeshBase * _mesh
std::vector< dof_id_type > getElementSolutionDofs(const Elem *elem, SystemBase &sys) const
Get a vector of the dof indices for all components of all variables associated with an element...
Definition: XFEM.C:2097

◆ setSolutionForDOFs()

void XFEM::setSolutionForDOFs ( const std::vector< Real > &  stored_solution,
const std::vector< dof_id_type > &  stored_solution_dofs,
NumericVector< Number > &  current_solution,
NumericVector< Number > &  old_solution,
NumericVector< Number > &  older_solution 
)
private

Set the solution for a set of DOFs.

Parameters
stored_solutionStored solution values to set the solution to
stored_solution_dofsDof indices for the entries in stored_solution
current_solutionCurrent solution vector that will be set
old_solutionOld solution vector that will be set
older_solutionOlder solution vector that will be set

Definition at line 2074 of file XFEM.C.

Referenced by setSolution().

2079 {
2080  // Solution vector is stored first for current, then old and older solutions.
2081  // These are the offsets to the beginning of the old and older solutions in the vector.
2082  const auto old_solution_offset = stored_solution_dofs.size();
2083  const auto older_solution_offset = old_solution_offset * 2;
2084 
2085  for (std::size_t i = 0; i < stored_solution_dofs.size(); ++i)
2086  {
2087  current_solution.set(stored_solution_dofs[i], stored_solution[i]);
2088  if (_fe_problem->isTransient())
2089  {
2090  old_solution.set(stored_solution_dofs[i], stored_solution[old_solution_offset + i]);
2091  older_solution.set(stored_solution_dofs[i], stored_solution[older_solution_offset + i]);
2092  }
2093  }
2094 }
FEProblemBase * _fe_problem
virtual void set(const numeric_index_type i, const Number value)=0
virtual bool isTransient() const override

◆ setXFEMQRule()

void XFEM::setXFEMQRule ( std::string &  xfem_qrule)

Definition at line 1819 of file XFEM.C.

1820 {
1821  if (xfem_qrule == "volfrac")
1823  else if (xfem_qrule == "moment_fitting")
1825  else if (xfem_qrule == "direct")
1827 }
Xfem::XFEM_QRULE _XFEM_qrule
Definition: XFEM.h:322

◆ storeCrackTipOriginAndDirection()

void XFEM::storeCrackTipOriginAndDirection ( )

Definition at line 184 of file XFEM.C.

Referenced by update().

185 {
187  std::set<EFAElement *> CrackTipElements = _efa_mesh.getCrackTipElements();
188  std::set<EFAElement *>::iterator sit;
189  for (sit = CrackTipElements.begin(); sit != CrackTipElements.end(); ++sit)
190  {
191  if (_mesh->mesh_dimension() == 2)
192  {
193  EFAElement2D * CEMElem = dynamic_cast<EFAElement2D *>(*sit);
194  EFANode * tip_node = CEMElem->getTipEmbeddedNode();
195  unsigned int cts_id = CEMElem->getCrackTipSplitElementID();
196 
197  Point origin(0, 0, 0);
198  Point direction(0, 0, 0);
199 
200  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
201  it = _cut_elem_map.find(_mesh->elem_ptr(cts_id)->unique_id());
202  if (it != _cut_elem_map.end())
203  {
204  const XFEMCutElem * xfce = it->second;
205  const EFAElement * EFAelem = xfce->getEFAElement();
206  if (EFAelem->isPartial()) // exclude the full crack tip elements
207  {
208  xfce->getCrackTipOriginAndDirection(tip_node->id(), origin, direction);
209  }
210  }
211 
212  std::vector<Point> tip_data;
213  tip_data.push_back(origin);
214  tip_data.push_back(direction);
215  const Elem * elem = _mesh->elem_ptr((*sit)->id());
217  std::pair<const Elem *, std::vector<Point>>(elem, tip_data));
218  }
219  }
220 }
unsigned int getCrackTipSplitElementID() const
Definition: EFAElement2D.C:599
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
virtual bool isPartial() const =0
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
virtual const EFAElement * getEFAElement() const =0
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335
EFANode * getTipEmbeddedNode() const
MeshBase * _mesh
const std::set< EFAElement * > & getCrackTipElements()
unsigned int id() const
Definition: EFANode.C:36
virtual void getCrackTipOriginAndDirection(unsigned tip_id, Point &origin, Point &direction) const =0

◆ storeMaterialPropertiesForElement()

void XFEM::storeMaterialPropertiesForElement ( const Elem *  parent_elem,
const Elem *  child_elem 
)
private

Helper function to store the material properties of a healed element.

Parameters
parent_elemThe parent element
elem1The first child element
elem2The second child element

Definition at line 2185 of file XFEM.C.

Referenced by healMesh().

2186 {
2187  // Set the parent element so that it is consistent post-healing
2188  _geom_cut_elems[child_elem]._parent_elem = parent_elem;
2189 
2190  // Locally store the element material properties
2192  _material_data[0]->getMaterialPropertyStorageForXFEM({}));
2193 
2194  // Locally store the boundary material properties
2195  // First check if any of the side need material properties
2196  bool need_boundary_materials = false;
2197  for (unsigned int side = 0; side < child_elem->n_sides(); ++side)
2198  {
2199  std::vector<boundary_id_type> elem_boundary_ids;
2200  _mesh->get_boundary_info().boundary_ids(child_elem, side, elem_boundary_ids);
2201  for (auto bdid : elem_boundary_ids)
2203  need_boundary_materials = true;
2204  }
2205 
2206  // If boundary material properties are needed for this element, then store them.
2207  if (need_boundary_materials)
2209  child_elem, _bnd_material_data[0]->getMaterialPropertyStorageForXFEM({}));
2210 }
bool needBoundaryMaterialOnSide(BoundaryID bnd_id, const THREAD_ID tid)
FEProblemBase * _fe_problem
std::vector< MaterialData *> _bnd_material_data
std::vector< MaterialData *> _material_data
void storeMaterialPropertiesForElementHelper(const Elem *elem, MaterialPropertyStorage &storage)
Definition: XFEM.C:2167
MeshBase * _mesh
std::unordered_map< const Elem *, Xfem::CutElemInfo > _geom_cut_elems
All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK...
Definition: XFEM.h:386

◆ storeMaterialPropertiesForElementHelper()

void XFEM::storeMaterialPropertiesForElementHelper ( const Elem *  elem,
MaterialPropertyStorage storage 
)
private

Definition at line 2167 of file XFEM.C.

Referenced by storeMaterialPropertiesForElement().

2168 {
2169  for (const auto state : storage.statefulIndexRange())
2170  {
2171  const auto & elem_props = storage.props(state).at(elem);
2172  auto & serialized_props = _geom_cut_elems[elem]._elem_material_properties[state - 1];
2173  serialized_props.clear();
2174  for (const auto & side_props_pair : elem_props)
2175  {
2176  const auto side = side_props_pair.first;
2177  std::ostringstream oss;
2178  dataStore(oss, storage.setProps(elem, side, state), nullptr);
2179  serialized_props[side].assign(oss.str());
2180  }
2181  }
2182 }
MaterialProperties & setProps(const Elem *elem, unsigned int side, const unsigned int state=0)
IntRange< unsigned int > statefulIndexRange() const
const PropsType & props(const unsigned int state=0) const
void dataStore(std::ostream &stream, FeatureFloodCount::FeatureData &feature, void *context)
std::unordered_map< const Elem *, Xfem::CutElemInfo > _geom_cut_elems
All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK...
Definition: XFEM.h:386

◆ storeSolutionForElement()

void XFEM::storeSolutionForElement ( const Elem *  elem_to_store_to,
const Elem *  elem_to_store_from,
SystemBase sys,
std::map< unique_id_type, std::vector< Real >> &  stored_solution,
const NumericVector< Number > &  current_solution,
const NumericVector< Number > &  old_solution,
const NumericVector< Number > &  older_solution 
)
private

Store the solution in stored_solution for a given element.

Parameters
elem_to_store_toElement for which the solution will be stored
elem_to_store_fromElement from which the solution to be stored is obtained
sysSystem from which the solution is stored
stored_solutionData structure that the stored solution is saved to
current_solutionCurrent solution vector that the solution is obtained from
old_solutionOld solution vector that the solution is obtained from
older_solutionOlder solution vector that the solution is obtained from

Definition at line 2002 of file XFEM.C.

Referenced by cutMeshWithEFA().

2009 {
2010  std::vector<dof_id_type> stored_solution_dofs = getElementSolutionDofs(elem_to_store_from, sys);
2011  std::vector<Real> stored_solution_scratch;
2012  // Size for current solution, as well as for old, and older solution only for transient case
2013  std::size_t stored_solution_size =
2014  (_fe_problem->isTransient() ? stored_solution_dofs.size() * 3 : stored_solution_dofs.size());
2015  stored_solution_scratch.reserve(stored_solution_size);
2016 
2017  // Store in the order defined in stored_solution_dofs first for the current, then for old and
2018  // older if applicable
2019  for (auto dof : stored_solution_dofs)
2020  stored_solution_scratch.push_back(current_solution(dof));
2021 
2022  if (_fe_problem->isTransient())
2023  {
2024  for (auto dof : stored_solution_dofs)
2025  stored_solution_scratch.push_back(old_solution(dof));
2026 
2027  for (auto dof : stored_solution_dofs)
2028  stored_solution_scratch.push_back(older_solution(dof));
2029  }
2030 
2031  if (stored_solution_scratch.size() > 0)
2032  stored_solution[elem_to_store_to->unique_id()] = stored_solution_scratch;
2033 }
FEProblemBase * _fe_problem
std::vector< dof_id_type > getElementSolutionDofs(const Elem *elem, SystemBase &sys) const
Get a vector of the dof indices for all components of all variables associated with an element...
Definition: XFEM.C:2097
virtual bool isTransient() const override

◆ storeSolutionForNode()

void XFEM::storeSolutionForNode ( const Node *  node_to_store_to,
const Node *  node_to_store_from,
SystemBase sys,
std::map< unique_id_type, std::vector< Real >> &  stored_solution,
const NumericVector< Number > &  current_solution,
const NumericVector< Number > &  old_solution,
const NumericVector< Number > &  older_solution 
)
private

Store the solution in stored_solution for a given node.

Parameters
node_to_store_toNode for which the solution will be stored
node_to_store_fromNode from which the solution to be stored is obtained
sysSystem from which the solution is stored
stored_solutionData structure that the stored solution is saved to
current_solutionCurrent solution vector that the solution is obtained from
old_solutionOld solution vector that the solution is obtained from
older_solutionOlder solution vector that the solution is obtained from

Definition at line 1968 of file XFEM.C.

Referenced by cutMeshWithEFA().

1975 {
1976  std::vector<dof_id_type> stored_solution_dofs = getNodeSolutionDofs(node_to_store_from, sys);
1977  std::vector<Real> stored_solution_scratch;
1978  // Size for current solution, as well as for old, and older solution only for transient case
1979  std::size_t stored_solution_size =
1980  (_fe_problem->isTransient() ? stored_solution_dofs.size() * 3 : stored_solution_dofs.size());
1981  stored_solution_scratch.reserve(stored_solution_size);
1982 
1983  // Store in the order defined in stored_solution_dofs first for the current, then for old and
1984  // older if applicable
1985  for (auto dof : stored_solution_dofs)
1986  stored_solution_scratch.push_back(current_solution(dof));
1987 
1988  if (_fe_problem->isTransient())
1989  {
1990  for (auto dof : stored_solution_dofs)
1991  stored_solution_scratch.push_back(old_solution(dof));
1992 
1993  for (auto dof : stored_solution_dofs)
1994  stored_solution_scratch.push_back(older_solution(dof));
1995  }
1996 
1997  if (stored_solution_scratch.size() > 0)
1998  stored_solution[node_to_store_to->unique_id()] = stored_solution_scratch;
1999 }
std::vector< dof_id_type > getNodeSolutionDofs(const Node *node, SystemBase &sys) const
Get a vector of the dof indices for all components of all variables associated with a node...
Definition: XFEM.C:2123
FEProblemBase * _fe_problem
virtual bool isTransient() const override

◆ update()

bool XFEM::update ( Real  time,
const std::vector< std::shared_ptr< NonlinearSystemBase >> &  nl,
AuxiliarySystem aux 
)
overridevirtual

Implements XFEMInterface.

Definition at line 260 of file XFEM.C.

263 {
265  mooseError("Use of XFEM with distributed mesh is not yet supported");
266 
267  bool mesh_changed = false;
268 
269  buildEFAMesh();
270 
272 
274 
275  if (markCuts(time))
276  mesh_changed = cutMeshWithEFA(nl, aux);
277 
278  if (mesh_changed)
279  {
280  buildEFAMesh();
282  }
283 
284  if (mesh_changed)
285  {
286  _mesh->allow_renumbering(false);
287  _mesh->skip_partitioning(true);
288  _mesh->prepare_for_use();
289 
290  if (_displaced_mesh)
291  {
292  _displaced_mesh->allow_renumbering(false);
293  _displaced_mesh->skip_partitioning(true);
294  _displaced_mesh->prepare_for_use();
295  }
296  }
297 
300 
301  return mesh_changed;
302 }
void mooseError(Args &&... args)
bool isDistributedMesh() const
bool cutMeshWithEFA(const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux)
Definition: XFEM.C:1101
void clearGeomMarkedElems()
Clear out the list of elements to be marked for cutting.
Definition: XFEM.C:177
FEProblemBase * _fe_problem
virtual void execute(const ExecFlagType &exec_type)
void storeCrackTipOriginAndDirection()
Definition: XFEM.C:184
bool markCuts(Real time)
Definition: XFEM.C:375
MeshBase * _mesh
const ExecFlagType EXEC_XFEM_MARK
Exec flag used to execute MooseObjects while elements are being marked for cutting by XFEM...
Definition: XFEMAppTypes.C:13
MooseMesh * _moose_mesh
MeshBase * _displaced_mesh
void clearStateMarkedElems()
Definition: XFEM.C:149
void buildEFAMesh()
Definition: XFEM.C:336

◆ updateHeal()

bool XFEM::updateHeal ( )
overridevirtual

Implements XFEMInterface.

Definition at line 223 of file XFEM.C.

224 {
225  bool mesh_changed = false;
226 
227  mesh_changed = healMesh();
228 
229  if (mesh_changed)
230  buildEFAMesh();
231 
232  if (mesh_changed)
233  {
234  _mesh->update_parallel_id_counts();
235  MeshCommunication().make_elems_parallel_consistent(*_mesh);
236  MeshCommunication().make_nodes_parallel_consistent(*_mesh);
237  // _mesh->find_neighbors();
238  // _mesh->contract();
239  _mesh->allow_renumbering(false);
240  _mesh->skip_partitioning(true);
241  _mesh->prepare_for_use();
242 
243  if (_displaced_mesh)
244  {
245  _displaced_mesh->update_parallel_id_counts();
246  MeshCommunication().make_elems_parallel_consistent(*_displaced_mesh);
247  MeshCommunication().make_nodes_parallel_consistent(*_displaced_mesh);
248  _displaced_mesh->allow_renumbering(false);
249  _displaced_mesh->skip_partitioning(true);
250  _displaced_mesh->prepare_for_use();
251  }
252  }
253 
254  _geom_marker_id_elems.clear();
255 
256  return mesh_changed;
257 }
bool healMesh()
Potentially heal the mesh by merging some of the pairs of partial elements cut by XFEM back into sing...
Definition: XFEM.C:917
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
Data structure for storing the elements cut by specific geometric cutters.
Definition: XFEM.h:350
MeshBase * _mesh
MeshBase * _displaced_mesh
void buildEFAMesh()
Definition: XFEM.C:336

Member Data Documentation

◆ _cached_aux_solution

std::map<unique_id_type, std::vector<Real> > XFEM::_cached_aux_solution
private

Data structure to store the auxiliary solution for nodes/elements affected by XFEM For each node/element, this is stored as a vector that contains all components of all applicable variables in an order defined by getElementSolutionDofs() or getNodeSolutionDofs().

This vector first contains the current solution in that order, followed by the old and older solutions, also in that same order.

Definition at line 380 of file XFEM.h.

Referenced by cutMeshWithEFA(), and initSolution().

◆ _cached_solution

std::map<unique_id_type, std::vector<Real> > XFEM::_cached_solution
private

Data structure to store the nonlinear solution for nodes/elements affected by XFEM For each node/element, this is stored as a vector that contains all components of all applicable variables in an order defined by getElementSolutionDofs() or getNodeSolutionDofs().

This vector first contains the current solution in that order, followed by the old and older solutions, also in that same order.

Definition at line 371 of file XFEM.h.

Referenced by cutMeshWithEFA(), and initSolution().

◆ _crack_growth_increment

Real XFEM::_crack_growth_increment
private

Definition at line 325 of file XFEM.h.

Referenced by markCutEdgesByState(), and setCrackGrowthMethod().

◆ _crack_tip_elems

std::set<const Elem *> XFEM::_crack_tip_elems
private

Definition at line 330 of file XFEM.h.

Referenced by cutMeshWithEFA(), healMesh(), and isElemAtCrackTip().

◆ _crack_tip_elems_to_be_healed

std::set<const Elem *> XFEM::_crack_tip_elems_to_be_healed
private

Definition at line 331 of file XFEM.h.

Referenced by cutMeshWithEFA(), and healMesh().

◆ _cut_elem_map

std::map<unique_id_type, XFEMCutElem *> XFEM::_cut_elem_map
private

◆ _debug_output_level

unsigned int XFEM::_debug_output_level
private

Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mesh 3: Full dump of element fragment algorithm mesh.

Definition at line 362 of file XFEM.h.

Referenced by cutMeshWithEFA(), healMesh(), and setDebugOutputLevel().

◆ _efa_mesh

ElementFragmentAlgorithm XFEM::_efa_mesh
private

◆ _elem_crack_origin_direction_map

std::map<const Elem *, std::vector<Point> > XFEM::_elem_crack_origin_direction_map
private

◆ _geom_cut_elems

std::unordered_map<const Elem *, Xfem::CutElemInfo> XFEM::_geom_cut_elems
private

All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK.

This data structure is updated everytime a new cut element is created.

Definition at line 386 of file XFEM.h.

Referenced by cutMeshWithEFA(), healMesh(), storeMaterialPropertiesForElement(), and storeMaterialPropertiesForElementHelper().

◆ _geom_marked_elems_2d

std::map<const Elem *, std::vector<Xfem::GeomMarkedElemInfo2D> > XFEM::_geom_marked_elems_2d
private

Data structure for storing information about all 2D elements to be cut by geometry.

Definition at line 344 of file XFEM.h.

Referenced by addGeomMarkedElem2D(), clearGeomMarkedElems(), and markCutEdgesByGeometry().

◆ _geom_marked_elems_3d

std::map<const Elem *, std::vector<Xfem::GeomMarkedElemInfo3D> > XFEM::_geom_marked_elems_3d
private

Data structure for storing information about all 3D elements to be cut by geometry.

Definition at line 347 of file XFEM.h.

Referenced by addGeomMarkedElem3D(), clearGeomMarkedElems(), and markCutFacesByGeometry().

◆ _geom_marker_id_elems

std::map<unsigned int, std::set<unsigned int> > XFEM::_geom_marker_id_elems
private

Data structure for storing the elements cut by specific geometric cutters.

Definition at line 350 of file XFEM.h.

Referenced by addGeomMarkedElem2D(), addGeomMarkedElem3D(), cutMeshWithEFA(), getGeometricCutForElem(), and updateHeal().

◆ _geom_marker_id_map

std::map<const GeometricCutUserObject *, unsigned int> XFEM::_geom_marker_id_map
private

Data structure for storing the GeommetricCutUserObjects and their corresponding id.

Definition at line 353 of file XFEM.h.

Referenced by addGeometricCut(), and getGeometricCutID().

◆ _geometric_cuts

std::vector<const GeometricCutUserObject *> XFEM::_geometric_cuts
private

Definition at line 327 of file XFEM.h.

Referenced by addGeometricCut(), cutMeshWithEFA(), getGeometricCutForElem(), and healMesh().

◆ _has_secondary_cut

bool XFEM::_has_secondary_cut
private

Definition at line 320 of file XFEM.h.

Referenced by has_secondary_cut(), markCutEdgesByGeometry(), markCutEdgesByState(), and XFEM().

◆ _old_geom_cut_elems

std::unordered_map<const Elem *, Xfem::CutElemInfo> XFEM::_old_geom_cut_elems
private

All geometrically cut elements and their CutElemInfo before the current execution of XFEM_MARK.

Definition at line 392 of file XFEM.h.

Referenced by cutMeshWithEFA().

◆ _sibling_displaced_elems

std::map<unsigned int, ElementPairLocator::ElementPairList> XFEM::_sibling_displaced_elems
private

Definition at line 333 of file XFEM.h.

Referenced by cutMeshWithEFA(), getXFEMDisplacedCutElemPairs(), and healMesh().

◆ _sibling_elems

std::map<unsigned int, ElementPairLocator::ElementPairList> XFEM::_sibling_elems
private

Definition at line 332 of file XFEM.h.

Referenced by cutMeshWithEFA(), getXFEMCutElemPairs(), and healMesh().

◆ _state_marked_elem_sides

std::map<const Elem *, unsigned int> XFEM::_state_marked_elem_sides
private

Definition at line 341 of file XFEM.h.

Referenced by addStateMarkedElem(), clearStateMarkedElems(), and markCutEdgesByState().

◆ _state_marked_elems

std::map<const Elem *, RealVectorValue> XFEM::_state_marked_elems
private

Definition at line 339 of file XFEM.h.

Referenced by addStateMarkedElem(), clearStateMarkedElems(), and markCutEdgesByState().

◆ _state_marked_frags

std::set<const Elem *> XFEM::_state_marked_frags
private

Definition at line 340 of file XFEM.h.

Referenced by addStateMarkedFrag(), clearStateMarkedElems(), and markCutEdgesByState().

◆ _use_crack_growth_increment

bool XFEM::_use_crack_growth_increment
private

Definition at line 324 of file XFEM.h.

Referenced by markCutEdgesByState(), and setCrackGrowthMethod().

◆ _XFEM_qrule

Xfem::XFEM_QRULE XFEM::_XFEM_qrule
private

Definition at line 322 of file XFEM.h.

Referenced by getXFEMQRule(), and setXFEMQRule().


The documentation for this class was generated from the following files: