www.mooseframework.org
Public Member Functions | 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, NonlinearSystemBase &nl, AuxiliarySystem &aux) override
 
virtual void initSolution (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 (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=NULL) 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)
 
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::ElementPairList * getXFEMCutElemPairs (unsigned int interface_id)
 Get the list of cut element pairs corresponding to a given interface ID. More...
 
virtual const ElementPairLocator::ElementPairList * getXFEMDisplacedCutElemPairs (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
 

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_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. More...
 
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. 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
 
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...
 

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 62 of file XFEM.h.

Constructor & Destructor Documentation

◆ XFEM()

XFEM::XFEM ( const InputParameters &  params)
explicit

Definition at line 35 of file XFEM.C.

35  : XFEMInterface(params), _efa_mesh(Moose::out)
36 {
37 #ifndef LIBMESH_ENABLE_UNIQUE_ID
38  mooseError("MOOSE requires unique ids to be enabled in libmesh (configure with "
39  "--enable-unique-id) to use XFEM!");
40 #endif
41  _has_secondary_cut = false;
42 }
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:288
bool _has_secondary_cut
Definition: XFEM.h:253

◆ ~XFEM()

XFEM::~XFEM ( )

Definition at line 44 of file XFEM.C.

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

Member Function Documentation

◆ addGeometricCut()

void XFEM::addGeometricCut ( GeometricCutUserObject geometric_cut)

Definition at line 53 of file XFEM.C.

54 {
55  _geometric_cuts.push_back(geometric_cut);
56 
57  geometric_cut->setInterfaceID(_geometric_cuts.size() - 1);
58 
59  _geom_marker_id_map[geometric_cut] = _geometric_cuts.size() - 1;
60 }
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:260
std::map< const GeometricCutUserObject *, unsigned int > _geom_marker_id_map
Data structure for storing the GeommetricCutUserObjects and their corresponding id.
Definition: XFEM.h:286
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 161 of file XFEM.C.

164 {
165  Elem * elem = _mesh->elem(elem_id);
166  _geom_marked_elems_2d[elem].push_back(geom_info);
167  _geom_marker_id_elems[interface_id].insert(elem_id);
168 }
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:277
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:283

◆ 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 171 of file XFEM.C.

174 {
175  Elem * elem = _mesh->elem(elem_id);
176  _geom_marked_elems_3d[elem].push_back(geom_info);
177  _geom_marker_id_elems[interface_id].insert(elem_id);
178 }
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:280
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:283

◆ 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(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 }
std::map< const Elem *, RealVectorValue > _state_marked_elems
Definition: XFEM.h:272

◆ 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(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  {
130  mooseError(" ERROR: side of element ", elem->id(), " already marked for crack initiation.");
131  exit(1);
132  }
133  _state_marked_elem_sides[elem] = marked_side;
134 }
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:274

◆ addStateMarkedFrag()

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

Definition at line 137 of file XFEM.C.

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

◆ buildEFAMesh()

void XFEM::buildEFAMesh ( )

Definition at line 339 of file XFEM.C.

Referenced by update().

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

Referenced by update().

182 {
183  _geom_marked_elems_2d.clear();
184  _geom_marked_elems_3d.clear();
185 }
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:280
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:277

◆ clearStateMarkedElems()

void XFEM::clearStateMarkedElems ( )

Definition at line 153 of file XFEM.C.

Referenced by update().

154 {
155  _state_marked_elems.clear();
156  _state_marked_frags.clear();
157  _state_marked_elem_sides.clear();
158 }
std::set< const Elem * > _state_marked_frags
Definition: XFEM.h:273
std::map< const Elem *, unsigned int > _state_marked_elem_sides
Definition: XFEM.h:274
std::map< const Elem *, RealVectorValue > _state_marked_elems
Definition: XFEM.h:272

◆ 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 446 of file XFEM.C.

Referenced by markCutEdgesByState().

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

◆ cutMeshWithEFA()

bool XFEM::cutMeshWithEFA ( NonlinearSystemBase &  nl,
AuxiliarySystem &  aux 
)

Definition at line 1056 of file XFEM.C.

Referenced by update().

1057 {
1058  std::map<unsigned int, Node *> efa_id_to_new_node;
1059  std::map<unsigned int, Node *> efa_id_to_new_node2;
1060  std::map<unsigned int, Elem *> efa_id_to_new_elem;
1061  _cached_solution.clear();
1062  _cached_aux_solution.clear();
1063 
1065  // DEBUG
1066  //_efa_mesh.printMesh();
1067 
1069  // DEBUG
1070  //_efa_mesh.printMesh();
1071 
1072  const std::vector<EFANode *> new_nodes = _efa_mesh.getNewNodes();
1073  const std::vector<EFAElement *> new_elements = _efa_mesh.getChildElements();
1074  const std::vector<EFAElement *> delete_elements = _efa_mesh.getParentElements();
1075 
1076  bool mesh_changed = (new_nodes.size() + new_elements.size() + delete_elements.size() > 0);
1077 
1078  // Prepare to cache solution on DOFs modified by XFEM
1079  if (mesh_changed)
1080  {
1081  nl.serializeSolution();
1082  aux.serializeSolution();
1083  }
1084  NumericVector<Number> & current_solution = *nl.system().current_local_solution;
1085  NumericVector<Number> & old_solution = nl.solutionOld();
1086  NumericVector<Number> & older_solution = nl.solutionOlder();
1087  NumericVector<Number> & current_aux_solution = *aux.system().current_local_solution;
1088  NumericVector<Number> & old_aux_solution = aux.solutionOld();
1089  NumericVector<Number> & older_aux_solution = aux.solutionOlder();
1090 
1091  std::map<Node *, Node *> new_nodes_to_parents;
1092 
1093  // Add new nodes
1094  for (unsigned int i = 0; i < new_nodes.size(); ++i)
1095  {
1096  unsigned int new_node_id = new_nodes[i]->id();
1097  unsigned int parent_id = new_nodes[i]->parent()->id();
1098 
1099  Node * parent_node = _mesh->node_ptr(parent_id);
1100  Node * new_node = Node::build(*parent_node, _mesh->n_nodes()).release();
1101  _mesh->add_node(new_node);
1102 
1103  new_nodes_to_parents[new_node] = parent_node;
1104 
1105  new_node->set_n_systems(parent_node->n_systems());
1106  efa_id_to_new_node.insert(std::make_pair(new_node_id, new_node));
1107  _console << "XFEM added new node: " << new_node->id() << "\n";
1108  if (_displaced_mesh)
1109  {
1110  const Node * parent_node2 = _displaced_mesh->node_ptr(parent_id);
1111  Node * new_node2 = Node::build(*parent_node2, _displaced_mesh->n_nodes()).release();
1112  _displaced_mesh->add_node(new_node2);
1113 
1114  new_node2->set_n_systems(parent_node2->n_systems());
1115  efa_id_to_new_node2.insert(std::make_pair(new_node_id, new_node2));
1116  }
1117  }
1118 
1119  // Add new elements
1120  std::map<unsigned int, std::vector<const Elem *>> temporary_parent_children_map;
1121 
1122  for (unsigned int i = 0; i < new_elements.size(); ++i)
1123  {
1124  unsigned int parent_id = new_elements[i]->getParent()->id();
1125  unsigned int efa_child_id = new_elements[i]->id();
1126 
1127  Elem * parent_elem = _mesh->elem(parent_id);
1128  Elem * libmesh_elem = Elem::build(parent_elem->type()).release();
1129 
1130  for (unsigned int m = 0; m < _geometric_cuts.size(); ++m)
1131  {
1132  for (auto & it : _sibling_elems[_geometric_cuts[m]->getInterfaceID()])
1133  {
1134  if (parent_elem == it.first)
1135  it.first = libmesh_elem;
1136  else if (parent_elem == it.second)
1137  it.second = libmesh_elem;
1138  }
1139  }
1140 
1141  // parent has at least two children
1142  if (new_elements[i]->getParent()->numChildren() > 1)
1143  temporary_parent_children_map[parent_elem->id()].push_back(libmesh_elem);
1144 
1145  Elem * parent_elem2 = NULL;
1146  Elem * libmesh_elem2 = NULL;
1147  if (_displaced_mesh)
1148  {
1149  parent_elem2 = _displaced_mesh->elem(parent_id);
1150  libmesh_elem2 = Elem::build(parent_elem2->type()).release();
1151 
1152  for (unsigned int m = 0; m < _geometric_cuts.size(); ++m)
1153  {
1154  for (auto & it : _sibling_displaced_elems[_geometric_cuts[m]->getInterfaceID()])
1155  {
1156  if (parent_elem2 == it.first)
1157  it.first = libmesh_elem2;
1158  else if (parent_elem2 == it.second)
1159  it.second = libmesh_elem2;
1160  }
1161  }
1162  }
1163 
1164  for (unsigned int j = 0; j < new_elements[i]->numNodes(); ++j)
1165  {
1166  unsigned int node_id = new_elements[i]->getNode(j)->id();
1167  Node * libmesh_node;
1168 
1169  std::map<unsigned int, Node *>::iterator nit = efa_id_to_new_node.find(node_id);
1170  if (nit != efa_id_to_new_node.end())
1171  libmesh_node = nit->second;
1172  else
1173  libmesh_node = _mesh->node_ptr(node_id);
1174 
1175  if (libmesh_node->processor_id() == DofObject::invalid_processor_id)
1176  libmesh_node->processor_id() = parent_elem->processor_id();
1177 
1178  libmesh_elem->set_node(j) = libmesh_node;
1179 
1180  // Store solution for all nodes affected by XFEM (even existing nodes)
1181  if (parent_elem->is_semilocal(_mesh->processor_id()))
1182  {
1183  Node * solution_node = libmesh_node; // Node from which to store solution
1184  if (new_nodes_to_parents.find(libmesh_node) != new_nodes_to_parents.end())
1185  solution_node = new_nodes_to_parents[libmesh_node];
1186 
1187  if ((_moose_mesh->isSemiLocal(solution_node)) ||
1188  (solution_node->processor_id() == _mesh->processor_id()))
1189  {
1190  storeSolutionForNode(libmesh_node,
1191  solution_node,
1192  nl,
1194  current_solution,
1195  old_solution,
1196  older_solution);
1197  storeSolutionForNode(libmesh_node,
1198  solution_node,
1199  aux,
1201  current_aux_solution,
1202  old_aux_solution,
1203  older_aux_solution);
1204  }
1205  }
1206 
1207  Node * parent_node = parent_elem->get_node(j);
1208  std::vector<boundary_id_type> parent_node_boundary_ids =
1209  _mesh->boundary_info->boundary_ids(parent_node);
1210  _mesh->boundary_info->add_node(libmesh_node, parent_node_boundary_ids);
1211 
1212  if (_displaced_mesh)
1213  {
1214  std::map<unsigned int, Node *>::iterator nit2 = efa_id_to_new_node2.find(node_id);
1215  if (nit2 != efa_id_to_new_node2.end())
1216  libmesh_node = nit2->second;
1217  else
1218  libmesh_node = _displaced_mesh->node_ptr(node_id);
1219 
1220  if (libmesh_node->processor_id() == DofObject::invalid_processor_id)
1221  libmesh_node->processor_id() = parent_elem2->processor_id();
1222 
1223  libmesh_elem2->set_node(j) = libmesh_node;
1224 
1225  parent_node = parent_elem2->get_node(j);
1226  parent_node_boundary_ids.clear();
1227  parent_node_boundary_ids = _displaced_mesh->boundary_info->boundary_ids(parent_node);
1228  _displaced_mesh->boundary_info->add_node(libmesh_node, parent_node_boundary_ids);
1229  }
1230  }
1231 
1232  libmesh_elem->set_p_level(parent_elem->p_level());
1233  libmesh_elem->set_p_refinement_flag(parent_elem->p_refinement_flag());
1234  _mesh->add_elem(libmesh_elem);
1235  libmesh_elem->set_n_systems(parent_elem->n_systems());
1236  libmesh_elem->subdomain_id() = parent_elem->subdomain_id();
1237  libmesh_elem->processor_id() = parent_elem->processor_id();
1238 
1239  // TODO: The 0 here is the thread ID. Need to sort out how to do this correctly
1240  // TODO: Also need to copy neighbor material data
1241  if (parent_elem->processor_id() == _mesh->processor_id())
1242  {
1243  (*_material_data)[0]->copy(*libmesh_elem, *parent_elem, 0);
1244  for (unsigned int side = 0; side < parent_elem->n_sides(); ++side)
1245  {
1246  std::vector<boundary_id_type> parent_elem_boundary_ids =
1247  _mesh->boundary_info->boundary_ids(parent_elem, side);
1248  std::vector<boundary_id_type>::iterator it_bd = parent_elem_boundary_ids.begin();
1249  for (; it_bd != parent_elem_boundary_ids.end(); ++it_bd)
1250  {
1251  if (_fe_problem->needBoundaryMaterialOnSide(*it_bd, 0))
1252  (*_bnd_material_data)[0]->copy(*libmesh_elem, *parent_elem, side);
1253  }
1254  }
1255 
1256  // Store solution for all elements affected by XFEM
1257  storeSolutionForElement(libmesh_elem,
1258  parent_elem,
1259  nl,
1261  current_solution,
1262  old_solution,
1263  older_solution);
1264  storeSolutionForElement(libmesh_elem,
1265  parent_elem,
1266  aux,
1268  current_aux_solution,
1269  old_aux_solution,
1270  older_aux_solution);
1271  }
1272 
1273  // The crack tip origin map is stored before cut, thus the elem should be updated with new
1274  // element.
1275  std::map<const Elem *, std::vector<Point>>::iterator mit =
1276  _elem_crack_origin_direction_map.find(parent_elem);
1277  if (mit != _elem_crack_origin_direction_map.end())
1278  {
1279  std::vector<Point> crack_data = _elem_crack_origin_direction_map[parent_elem];
1281  _elem_crack_origin_direction_map[libmesh_elem] = crack_data;
1282  }
1283 
1284  _console << "XFEM added new element: " << libmesh_elem->id() << "\n";
1285 
1286  XFEMCutElem * xfce = NULL;
1287  if (_mesh->mesh_dimension() == 2)
1288  {
1289  EFAElement2D * new_efa_elem2d = dynamic_cast<EFAElement2D *>(new_elements[i]);
1290  if (!new_efa_elem2d)
1291  mooseError("EFAelem is not of EFAelement2D type");
1292  xfce = new XFEMCutElem2D(libmesh_elem,
1293  new_efa_elem2d,
1294  _fe_problem->assembly(0).qRule()->n_points(),
1295  libmesh_elem->n_sides());
1296  }
1297  else if (_mesh->mesh_dimension() == 3)
1298  {
1299  EFAElement3D * new_efa_elem3d = dynamic_cast<EFAElement3D *>(new_elements[i]);
1300  if (!new_efa_elem3d)
1301  mooseError("EFAelem is not of EFAelement3D type");
1302  xfce = new XFEMCutElem3D(libmesh_elem,
1303  new_efa_elem3d,
1304  _fe_problem->assembly(0).qRule()->n_points(),
1305  libmesh_elem->n_sides());
1306  }
1307  _cut_elem_map.insert(std::pair<unique_id_type, XFEMCutElem *>(libmesh_elem->unique_id(), xfce));
1308  efa_id_to_new_elem.insert(std::make_pair(efa_child_id, libmesh_elem));
1309 
1310  if (_displaced_mesh)
1311  {
1312  libmesh_elem2->set_p_level(parent_elem2->p_level());
1313  libmesh_elem2->set_p_refinement_flag(parent_elem2->p_refinement_flag());
1314  _displaced_mesh->add_elem(libmesh_elem2);
1315  libmesh_elem2->set_n_systems(parent_elem2->n_systems());
1316  libmesh_elem2->subdomain_id() = parent_elem2->subdomain_id();
1317  libmesh_elem2->processor_id() = parent_elem2->processor_id();
1318  }
1319 
1320  unsigned int n_sides = parent_elem->n_sides();
1321  for (unsigned int side = 0; side < n_sides; ++side)
1322  {
1323  std::vector<boundary_id_type> parent_elem_boundary_ids =
1324  _mesh->boundary_info->boundary_ids(parent_elem, side);
1325  _mesh->boundary_info->add_side(libmesh_elem, side, parent_elem_boundary_ids);
1326  }
1327  if (_displaced_mesh)
1328  {
1329  n_sides = parent_elem2->n_sides();
1330  for (unsigned int side = 0; side < n_sides; ++side)
1331  {
1332  std::vector<boundary_id_type> parent_elem_boundary_ids =
1333  _displaced_mesh->boundary_info->boundary_ids(parent_elem2, side);
1334  _displaced_mesh->boundary_info->add_side(libmesh_elem2, side, parent_elem_boundary_ids);
1335  }
1336  }
1337 
1338  unsigned int n_edges = parent_elem->n_edges();
1339  for (unsigned int edge = 0; edge < n_edges; ++edge)
1340  {
1341  std::vector<boundary_id_type> parent_elem_boundary_ids =
1342  _mesh->boundary_info->edge_boundary_ids(parent_elem, edge);
1343  _mesh->boundary_info->add_edge(libmesh_elem, edge, parent_elem_boundary_ids);
1344  }
1345  if (_displaced_mesh)
1346  {
1347  n_edges = parent_elem2->n_edges();
1348  for (unsigned int edge = 0; edge < n_edges; ++edge)
1349  {
1350  std::vector<boundary_id_type> parent_elem_boundary_ids =
1351  _displaced_mesh->boundary_info->edge_boundary_ids(parent_elem2, edge);
1352  _displaced_mesh->boundary_info->add_edge(libmesh_elem2, edge, parent_elem_boundary_ids);
1353  }
1354  }
1355  }
1356 
1357  // delete elements
1358  for (std::size_t i = 0; i < delete_elements.size(); ++i)
1359  {
1360  Elem * elem_to_delete = _mesh->elem(delete_elements[i]->id());
1361 
1362  // delete the XFEMCutElem object for any elements that are to be deleted
1363  std::map<unique_id_type, XFEMCutElem *>::iterator cemit =
1364  _cut_elem_map.find(elem_to_delete->unique_id());
1365  if (cemit != _cut_elem_map.end())
1366  {
1367  delete cemit->second;
1368  _cut_elem_map.erase(cemit);
1369  }
1370 
1371  elem_to_delete->nullify_neighbors();
1372  _mesh->boundary_info->remove(elem_to_delete);
1373  unsigned int deleted_elem_id = elem_to_delete->id();
1374  _mesh->delete_elem(elem_to_delete);
1375  _console << "XFEM deleted element: " << deleted_elem_id << "\n";
1376 
1377  if (_displaced_mesh)
1378  {
1379  Elem * elem_to_delete2 = _displaced_mesh->elem(delete_elements[i]->id());
1380  elem_to_delete2->nullify_neighbors();
1381  _displaced_mesh->boundary_info->remove(elem_to_delete2);
1382  _displaced_mesh->delete_elem(elem_to_delete2);
1383  }
1384  }
1385 
1386  for (std::map<unsigned int, std::vector<const Elem *>>::iterator it =
1387  temporary_parent_children_map.begin();
1388  it != temporary_parent_children_map.end();
1389  ++it)
1390  {
1391  std::vector<const Elem *> & sibling_elem_vec = it->second;
1392  // TODO: for cut-node case, how to find the sibling elements?
1393  // if (sibling_elem_vec.size() != 2)
1394  // mooseError("Must have exactly 2 sibling elements");
1395 
1396  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1397  for (auto const & elem_id : _geom_marker_id_elems[_geometric_cuts[i]->getInterfaceID()])
1398  if (it->first == elem_id)
1399  _sibling_elems[_geometric_cuts[i]->getInterfaceID()].push_back(
1400  std::make_pair(sibling_elem_vec[0], sibling_elem_vec[1]));
1401  }
1402 
1403  // add sibling elems on displaced mesh
1404  if (_displaced_mesh)
1405  {
1406  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1407  {
1408  for (auto & se : _sibling_elems[_geometric_cuts[i]->getInterfaceID()])
1409  {
1410  Elem * elem = _displaced_mesh->elem(se.first->id());
1411  Elem * elem_pair = _displaced_mesh->elem(se.second->id());
1412  _sibling_displaced_elems[_geometric_cuts[i]->getInterfaceID()].push_back(
1413  std::make_pair(elem, elem_pair));
1414  }
1415  }
1416  }
1417 
1418  // clear the temporary map
1419  temporary_parent_children_map.clear();
1420 
1421  // Store information about crack tip elements
1422  if (mesh_changed)
1423  {
1424  _crack_tip_elems.clear();
1426  const std::set<EFAElement *> CrackTipElements = _efa_mesh.getCrackTipElements();
1427  std::set<EFAElement *>::const_iterator sit;
1428  for (sit = CrackTipElements.begin(); sit != CrackTipElements.end(); ++sit)
1429  {
1430  unsigned int eid = (*sit)->id();
1431  Elem * crack_tip_elem;
1432  std::map<unsigned int, Elem *>::iterator eit = efa_id_to_new_elem.find(eid);
1433  if (eit != efa_id_to_new_elem.end())
1434  crack_tip_elem = eit->second;
1435  else
1436  crack_tip_elem = _mesh->elem(eid);
1437  _crack_tip_elems.insert(crack_tip_elem);
1438 
1439  // Store the crack tip elements which are going to be healed
1440  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1441  {
1442  if (_geometric_cuts[i]->shouldHealMesh())
1443  {
1444  for (auto const & mie : _geom_marker_id_elems[_geometric_cuts[i]->getInterfaceID()])
1445  if ((*sit)->getParent() != nullptr)
1446  {
1447  if ((*sit)->getParent()->id() == mie)
1448  _crack_tip_elems_to_be_healed.insert(crack_tip_elem);
1449  }
1450  }
1451  }
1452  }
1453  }
1454  _console << std::flush;
1455 
1456  // store virtual nodes
1457  // store cut edge info
1458  return mesh_changed;
1459 }
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:260
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:265
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:266
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
const std::vector< EFAElement * > & getChildElements()
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:263
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:288
std::set< const Elem * > _crack_tip_elems_to_be_healed
Definition: XFEM.h:264
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:306
void updateTopology(bool mergeUncutVirtualEdges=true)
const std::vector< EFAElement * > & getParentElements()
const std::vector< EFANode * > & getNewNodes()
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:1863
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:1829
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:268
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:283
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:297
const std::set< EFAElement * > & getCrackTipElements()

◆ getCrackTipOrigin()

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

Definition at line 63 of file XFEM.C.

65 {
66  elem_id_crack_tip.clear();
67  crack_front_points.clear();
68  crack_front_points.resize(_elem_crack_origin_direction_map.size());
69 
70  unsigned int crack_tip_index = 0;
71  // This map is used to sort the order in _elem_crack_origin_direction_map such that every process
72  // has same order
73  std::map<unsigned int, const Elem *> elem_id_map;
74 
75  int m = -1;
76  for (std::map<const Elem *, std::vector<Point>>::iterator mit1 =
79  ++mit1)
80  {
81  unsigned int elem_id = mit1->first->id();
82  if (elem_id > 999999)
83  {
84  elem_id_map[m] = mit1->first;
85  m--;
86  }
87  else
88  {
89  elem_id_map[elem_id] = mit1->first;
90  }
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:268

◆ getCrackTipOriginMap()

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

Definition at line 247 of file XFEM.h.

248  {
250  }
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:268

◆ 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 1532 of file XFEM.C.

1535 {
1536  Real comp = 0.0;
1537  Point planedata(0.0, 0.0, 0.0);
1538  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1539  it = _cut_elem_map.find(elem->unique_id());
1540  if (it != _cut_elem_map.end())
1541  {
1542  const XFEMCutElem * xfce = it->second;
1543  const EFAElement * EFAelem = xfce->getEFAElement();
1544  if (EFAelem->isPartial()) // exclude the full crack tip elements
1545  {
1546  if ((unsigned int)quantity < 3)
1547  {
1548  unsigned int index = (unsigned int)quantity;
1549  planedata = xfce->getCutPlaneOrigin(plane_id, _displaced_mesh);
1550  comp = planedata(index);
1551  }
1552  else if ((unsigned int)quantity < 6)
1553  {
1554  unsigned int index = (unsigned int)quantity - 3;
1555  planedata = xfce->getCutPlaneNormal(plane_id, _displaced_mesh);
1556  comp = planedata(index);
1557  }
1558  else
1559  mooseError("In get_cut_plane index out of range");
1560  }
1561  }
1562  return comp;
1563 }
virtual Point getCutPlaneOrigin(unsigned int plane_id, MeshBase *displaced_mesh=NULL) const =0
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
virtual bool isPartial() const =0
virtual const EFAElement * getEFAElement() const =0
virtual Point getCutPlaneNormal(unsigned int plane_id, MeshBase *displaced_mesh=NULL) const =0

◆ 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 1613 of file XFEM.C.

Referenced by markCutEdgesByGeometry(), and markCutEdgesByState().

1614 {
1615  EFAElement * EFAelem = _efa_mesh.getElemByID(elem->id());
1616  EFAElement2D * EFAelem2D = dynamic_cast<EFAElement2D *>(EFAelem);
1617 
1618  if (!EFAelem2D)
1619  mooseError("EFAelem is not of EFAelement2D type");
1620 
1621  return EFAelem2D;
1622 }
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:288
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 1625 of file XFEM.C.

Referenced by markCutFacesByGeometry().

1626 {
1627  EFAElement * EFAelem = _efa_mesh.getElemByID(elem->id());
1628  EFAElement3D * EFAelem3D = dynamic_cast<EFAElement3D *>(EFAelem);
1629 
1630  if (!EFAelem3D)
1631  mooseError("EFAelem is not of EFAelement3D type");
1632 
1633  return EFAelem3D;
1634 }
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:288
EFAElement * getElemByID(unsigned int id)

◆ getEFANodeCoords()

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

Definition at line 1462 of file XFEM.C.

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

1466 {
1467  Point node_coor(0.0, 0.0, 0.0);
1468  std::vector<EFANode *> master_nodes;
1469  std::vector<Point> master_points;
1470  std::vector<double> master_weights;
1471 
1472  CEMElem->getMasterInfo(CEMnode, master_nodes, master_weights);
1473  for (std::size_t i = 0; i < master_nodes.size(); ++i)
1474  {
1475  if (master_nodes[i]->category() == EFANode::N_CATEGORY_PERMANENT)
1476  {
1477  unsigned int local_node_id = CEMElem->getLocalNodeIndex(master_nodes[i]);
1478  Node * node = elem->get_node(local_node_id);
1479  if (displaced_mesh)
1480  node = displaced_mesh->node_ptr(node->id());
1481  Point node_p((*node)(0), (*node)(1), (*node)(2));
1482  master_points.push_back(node_p);
1483  }
1484  else
1485  mooseError("master nodes must be permanent");
1486  }
1487  for (std::size_t i = 0; i < master_nodes.size(); ++i)
1488  node_coor += master_weights[i] * master_points[i];
1489 
1490  return node_coor;
1491 }
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 1958 of file XFEM.C.

Referenced by setSolution(), and storeSolutionForElement().

1959 {
1960  SubdomainID sid = elem->subdomain_id();
1961  const std::vector<MooseVariableFEBase *> & vars = sys.getVariables(0);
1962  std::vector<dof_id_type> solution_dofs;
1963  solution_dofs.reserve(vars.size()); // just an approximation
1964  for (auto var : vars)
1965  {
1966  if (!var->isNodal())
1967  {
1968  const std::set<SubdomainID> & var_subdomains = sys.getSubdomainsForVar(var->number());
1969  if (var_subdomains.empty() || var_subdomains.find(sid) != var_subdomains.end())
1970  {
1971  unsigned int n_comp = elem->n_comp(sys.number(), var->number());
1972  for (unsigned int icomp = 0; icomp < n_comp; ++icomp)
1973  {
1974  dof_id_type elem_dof = elem->dof_number(sys.number(), var->number(), icomp);
1975  solution_dofs.push_back(elem_dof);
1976  }
1977  }
1978  }
1979  }
1980  return solution_dofs;
1981 }

◆ getFragmentEdges()

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

Definition at line 1637 of file XFEM.C.

1640 {
1641  // N.B. CEMElem here has global EFAnode
1642  frag_edges.clear();
1643  if (CEMElem->numFragments() > 0)
1644  {
1645  if (CEMElem->numFragments() > 1)
1646  mooseError("element ", elem->id(), " has more than one fragment at this point");
1647  for (unsigned int i = 0; i < CEMElem->getFragment(0)->numEdges(); ++i)
1648  {
1649  std::vector<Point> p_line(2, Point(0.0, 0.0, 0.0));
1650  p_line[0] = getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(0), CEMElem, elem);
1651  p_line[1] = getEFANodeCoords(CEMElem->getFragmentEdge(0, i)->getNode(1), CEMElem, elem);
1652  frag_edges.push_back(p_line);
1653  }
1654  }
1655 }
unsigned int numEdges() const
virtual unsigned int numFragments() const
Definition: EFAElement2D.C:202
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=NULL) const
Definition: XFEM.C:1462
EFAFragment2D * getFragment(unsigned int frag_id) const
EFANode * getNode(unsigned int index) const
Definition: EFAEdge.C:179
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 1596 of file XFEM.C.

1599 {
1600  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1601  it = _cut_elem_map.find(elem->unique_id());
1602  if (it != _cut_elem_map.end())
1603  {
1604  const XFEMCutElem * xfce = it->second;
1605  if (displaced_mesh)
1606  xfce->getFragmentFaces(frag_faces, _displaced_mesh);
1607  else
1608  xfce->getFragmentFaces(frag_faces);
1609  }
1610 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
virtual void getFragmentFaces(std::vector< std::vector< Point >> &frag_faces, MeshBase *displaced_mesh=NULL) const =0

◆ getFragmentFaces() [2/2]

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

Definition at line 1658 of file XFEM.C.

1661 {
1662  // N.B. CEMElem here has global EFAnode
1663  frag_faces.clear();
1664  if (CEMElem->numFragments() > 0)
1665  {
1666  if (CEMElem->numFragments() > 1)
1667  mooseError("element ", elem->id(), " has more than one fragment at this point");
1668  for (unsigned int i = 0; i < CEMElem->getFragment(0)->numFaces(); ++i)
1669  {
1670  unsigned int num_face_nodes = CEMElem->getFragmentFace(0, i)->numNodes();
1671  std::vector<Point> p_line(num_face_nodes, Point(0.0, 0.0, 0.0));
1672  for (unsigned int j = 0; j < num_face_nodes; ++j)
1673  p_line[j] = getEFANodeCoords(CEMElem->getFragmentFace(0, i)->getNode(j), CEMElem, elem);
1674  frag_faces.push_back(p_line);
1675  }
1676  }
1677 }
EFAFragment3D * getFragment(unsigned int frag_id) const
virtual unsigned int numFragments() const
Definition: EFAElement3D.C:269
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=NULL) const
Definition: XFEM.C:1462
EFANode * getNode(unsigned int node_id) const
Definition: EFAFace.C:99
unsigned int numNodes() const
Definition: EFAFace.C:87
EFAFace * getFragmentFace(unsigned int frag_id, unsigned int face_id) const
unsigned int numFaces() const

◆ 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 210 of file XFEM.h.

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

◆ 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 1984 of file XFEM.C.

Referenced by setSolution(), and storeSolutionForNode().

1985 {
1986  const std::set<SubdomainID> & sids = _moose_mesh->getNodeBlockIds(*node);
1987  const std::vector<MooseVariableFEBase *> & vars = sys.getVariables(0);
1988  std::vector<dof_id_type> solution_dofs;
1989  solution_dofs.reserve(vars.size()); // just an approximation
1990  for (auto var : vars)
1991  {
1992  if (var->isNodal())
1993  {
1994  const std::set<SubdomainID> & var_subdomains = sys.getSubdomainsForVar(var->number());
1995  std::set<SubdomainID> intersect;
1996  set_intersection(var_subdomains.begin(),
1997  var_subdomains.end(),
1998  sids.begin(),
1999  sids.end(),
2000  std::inserter(intersect, intersect.begin()));
2001  if (var_subdomains.empty() || !intersect.empty())
2002  {
2003  unsigned int n_comp = node->n_comp(sys.number(), var->number());
2004  for (unsigned int icomp = 0; icomp < n_comp; ++icomp)
2005  {
2006  dof_id_type node_dof = node->dof_number(sys.number(), var->number(), icomp);
2007  solution_dofs.push_back(node_dof);
2008  }
2009  }
2010  }
2011  }
2012  return solution_dofs;
2013 }

◆ getPhysicalVolumeFraction()

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

Get the volume fraction of an element that is physical.

Definition at line 1494 of file XFEM.C.

Referenced by markCutEdgesByState().

1495 {
1496  Real phys_volfrac = 1.0;
1497  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1498  it = _cut_elem_map.find(elem->unique_id());
1499  if (it != _cut_elem_map.end())
1500  {
1501  XFEMCutElem * xfce = it->second;
1502  const EFAElement * EFAelem = xfce->getEFAElement();
1503  if (EFAelem->isPartial())
1504  { // exclude the full crack tip elements
1506  phys_volfrac = xfce->getPhysicalVolumeFraction();
1507  }
1508  }
1509 
1510  return phys_volfrac;
1511 }
Real getPhysicalVolumeFraction() const
Returns the volume fraction of the element fragment.
Definition: XFEMCutElem.C:49
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
virtual bool isPartial() const =0
virtual const EFAElement * getEFAElement() const =0
virtual void computePhysicalVolumeFraction()=0
Computes the volume fraction of the element fragment.

◆ 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 188 of file XFEM.h.

189  {
190  return &_sibling_elems[interface_id];
191  }
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:265

◆ 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 200 of file XFEM.h.

201  {
202  return &_sibling_displaced_elems[interface_id];
203  }
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:266

◆ getXFEMFaceWeights()

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

Definition at line 1721 of file XFEM.C.

1726 {
1727  bool have_weights = false;
1728  XFEMCutElem * xfce = NULL;
1729  if (isElemCut(elem, xfce))
1730  {
1731  mooseAssert(xfce != NULL, "Must have valid XFEMCutElem object here");
1732  xfce->getFaceWeightMultipliers(weights, qrule, getXFEMQRule(), q_points, side);
1733  have_weights = true;
1734  }
1735  return have_weights;
1736 }
void getFaceWeightMultipliers(MooseArray< Real > &face_weights, QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points, unsigned int side)
Definition: XFEMCutElem.C:75
Xfem::XFEM_QRULE & getXFEMQRule()
Definition: XFEM.C:1680
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1572

◆ 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 1739 of file XFEM.C.

1744 {
1745  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1746  it = _cut_elem_map.find(elem->unique_id());
1747  if (it != _cut_elem_map.end())
1748  {
1749  const XFEMCutElem * xfce = it->second;
1750  if (displaced_mesh)
1751  xfce->getIntersectionInfo(plane_id, normal, intersectionPoints, _displaced_mesh);
1752  else
1753  xfce->getIntersectionInfo(plane_id, normal, intersectionPoints);
1754  }
1755 }
virtual void getIntersectionInfo(unsigned int plane_id, Point &normal, std::vector< Point > &intersectionPoints, MeshBase *displaced_mesh=NULL) const =0
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262

◆ getXFEMQRule()

Xfem::XFEM_QRULE & XFEM::getXFEMQRule ( )

Definition at line 1680 of file XFEM.C.

Referenced by getXFEMFaceWeights(), and getXFEMWeights().

1681 {
1682  return _XFEM_qrule;
1683 }
Xfem::XFEM_QRULE _XFEM_qrule
Definition: XFEM.h:255

◆ getXFEMqRuleOnLine()

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

Definition at line 1758 of file XFEM.C.

1761 {
1762  Point p1 = intersection_points[0];
1763  Point p2 = intersection_points[1];
1764 
1765  // number of quadrature points
1766  std::size_t num_qpoints = 2;
1767 
1768  // quadrature coordinates
1769  Real xi0 = -std::sqrt(1.0 / 3.0);
1770  Real xi1 = std::sqrt(1.0 / 3.0);
1771 
1772  quad_wts.resize(num_qpoints);
1773  quad_pts.resize(num_qpoints);
1774 
1775  Real integ_jacobian = pow((p1 - p2).size_sq(), 0.5) * 0.5;
1776 
1777  quad_wts[0] = 1.0 * integ_jacobian;
1778  quad_wts[1] = 1.0 * integ_jacobian;
1779 
1780  quad_pts[0] = (1.0 - xi0) / 2.0 * p1 + (1.0 + xi0) / 2.0 * p2;
1781  quad_pts[1] = (1.0 - xi1) / 2.0 * p1 + (1.0 + xi1) / 2.0 * p2;
1782 }
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ getXFEMqRuleOnSurface()

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

Definition at line 1785 of file XFEM.C.

1788 {
1789  std::size_t nnd_pe = intersection_points.size();
1790  Point xcrd(0.0, 0.0, 0.0);
1791  for (std::size_t i = 0; i < nnd_pe; ++i)
1792  xcrd += intersection_points[i];
1793  xcrd /= nnd_pe;
1794 
1795  quad_pts.resize(nnd_pe);
1796  quad_wts.resize(nnd_pe);
1797 
1798  Real jac = 0.0;
1799 
1800  for (std::size_t j = 0; j < nnd_pe; ++j) // loop all sub-tris
1801  {
1802  std::vector<std::vector<Real>> shape(3, std::vector<Real>(3, 0.0));
1803  std::vector<Point> subtrig_points(3, Point(0.0, 0.0, 0.0)); // sub-trig nodal coords
1804 
1805  int jplus1 = j < nnd_pe - 1 ? j + 1 : 0;
1806  subtrig_points[0] = xcrd;
1807  subtrig_points[1] = intersection_points[j];
1808  subtrig_points[2] = intersection_points[jplus1];
1809 
1810  std::vector<std::vector<Real>> sg2;
1811  Xfem::stdQuadr2D(3, 1, sg2); // get sg2
1812  for (std::size_t l = 0; l < sg2.size(); ++l) // loop all int pts on a sub-trig
1813  {
1814  Xfem::shapeFunc2D(3, sg2[l], subtrig_points, shape, jac, true); // Get shape
1815  std::vector<Real> tsg_line(3, 0.0);
1816  for (std::size_t k = 0; k < 3; ++k) // loop sub-trig nodes
1817  {
1818  tsg_line[0] += shape[k][2] * subtrig_points[k](0);
1819  tsg_line[1] += shape[k][2] * subtrig_points[k](1);
1820  tsg_line[2] += shape[k][2] * subtrig_points[k](2);
1821  }
1822  quad_pts[j + l] = Point(tsg_line[0], tsg_line[1], tsg_line[2]);
1823  quad_wts[j + l] = sg2[l][3] * jac;
1824  }
1825  }
1826 }
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
void stdQuadr2D(unsigned int nen, unsigned int iord, std::vector< std::vector< Real >> &sg2)
Definition: XFEMFuncs.C:96

◆ getXFEMWeights()

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

Definition at line 1704 of file XFEM.C.

1708 {
1709  bool have_weights = false;
1710  XFEMCutElem * xfce = NULL;
1711  if (isElemCut(elem, xfce))
1712  {
1713  mooseAssert(xfce != NULL, "Must have valid XFEMCutElem object here");
1714  xfce->getWeightMultipliers(weights, qrule, getXFEMQRule(), q_points);
1715  have_weights = true;
1716  }
1717  return have_weights;
1718 }
Xfem::XFEM_QRULE & getXFEMQRule()
Definition: XFEM.C:1680
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1572
void getWeightMultipliers(MooseArray< Real > &weights, QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
Definition: XFEMCutElem.C:61

◆ has_secondary_cut()

bool XFEM::has_secondary_cut ( )
inline

Definition at line 226 of file XFEM.h.

226 { return _has_secondary_cut; }
bool _has_secondary_cut
Definition: XFEM.h:253

◆ 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 919 of file XFEM.C.

Referenced by updateHeal().

920 {
921  bool mesh_changed = false;
922 
923  std::set<Node *> nodes_to_delete;
924  std::set<Node *> nodes_to_delete_displaced;
925  std::set<unsigned int> cutelems_to_delete;
926 
927  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
928  {
929  if (_geometric_cuts[i]->shouldHealMesh())
930  {
931  for (auto & it : _sibling_elems[_geometric_cuts[i]->getInterfaceID()])
932  {
933  Elem * elem1 = const_cast<Elem *>(it.first);
934  Elem * elem2 = const_cast<Elem *>(it.second);
935 
936  std::map<unique_id_type, XFEMCutElem *>::iterator cemit =
937  _cut_elem_map.find(elem1->unique_id());
938  if (cemit != _cut_elem_map.end())
939  {
940  const XFEMCutElem * xfce = cemit->second;
941 
942  cutelems_to_delete.insert(elem1->unique_id());
943 
944  for (unsigned int in = 0; in < elem1->n_nodes(); ++in)
945  {
946  Node * e1node = elem1->get_node(in);
947  Node * e2node = elem2->get_node(in);
948  if (!xfce->isPointPhysical(*e1node) &&
949  e1node != e2node) // This would happen at the crack tip
950  {
951  elem1->set_node(in) = e2node;
952  nodes_to_delete.insert(e1node);
953  }
954  else if (e1node != e2node)
955  nodes_to_delete.insert(e2node);
956  }
957  }
958  else
959  mooseError("Could not find XFEMCutElem for element to be kept in healing");
960 
961  if (_displaced_mesh)
962  {
963  Elem * elem1_displaced = _displaced_mesh->elem(it.first->id());
964  Elem * elem2_displaced = _displaced_mesh->elem(it.second->id());
965 
966  std::map<unique_id_type, XFEMCutElem *>::iterator cemit =
967  _cut_elem_map.find(elem1_displaced->unique_id());
968  if (cemit != _cut_elem_map.end())
969  {
970  const XFEMCutElem * xfce = cemit->second;
971 
972  for (unsigned int in = 0; in < elem1_displaced->n_nodes(); ++in)
973  {
974  Node * e1node_displaced = elem1_displaced->get_node(in);
975  Node * e2node_displaced = elem2_displaced->get_node(in);
976  if (!xfce->isPointPhysical(*e1node_displaced) &&
977  e1node_displaced != e2node_displaced) // This would happen at the crack tip
978  {
979  elem1_displaced->set_node(in) = e2node_displaced;
980  nodes_to_delete_displaced.insert(e1node_displaced);
981  }
982  else if (e1node_displaced != e2node_displaced)
983  nodes_to_delete_displaced.insert(e2node_displaced);
984  }
985  }
986  else
987  mooseError("Could not find XFEMCutElem for element to be kept in healing");
988 
989  elem2_displaced->nullify_neighbors();
990  _displaced_mesh->boundary_info->remove(elem2_displaced);
991  _displaced_mesh->delete_elem(elem2_displaced);
992  }
993 
994  cutelems_to_delete.insert(elem2->unique_id());
995  elem2->nullify_neighbors();
996  _mesh->boundary_info->remove(elem2);
997  unsigned int deleted_elem_id = elem2->id();
998  _mesh->delete_elem(elem2);
999  _console << "XFEM healing deleted element: " << deleted_elem_id << "\n";
1000  mesh_changed = true;
1001  }
1002  }
1003  }
1004 
1005  for (auto & sit : nodes_to_delete)
1006  {
1007  Node * node_to_delete = sit;
1008  dof_id_type deleted_node_id = node_to_delete->id();
1009  _mesh->boundary_info->remove(node_to_delete);
1010  _mesh->delete_node(node_to_delete);
1011  _console << "XFEM healing deleted node: " << deleted_node_id << "\n";
1012  }
1013 
1014  _console << std::flush;
1015 
1016  if (_displaced_mesh)
1017  {
1018  for (auto & sit : nodes_to_delete_displaced)
1019  {
1020  Node * node_to_delete = sit;
1021  _displaced_mesh->boundary_info->remove(node_to_delete);
1022  _displaced_mesh->delete_node(node_to_delete);
1023  }
1024  }
1025 
1026  for (auto & ced : cutelems_to_delete)
1027  if (_cut_elem_map.find(ced) != _cut_elem_map.end())
1028  {
1029  delete _cut_elem_map.find(ced)->second;
1030  _cut_elem_map.erase(ced);
1031  }
1032 
1033  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1034  if (_geometric_cuts[i]->shouldHealMesh())
1035  _sibling_elems[_geometric_cuts[i]->getInterfaceID()].clear();
1036 
1037  if (_displaced_mesh)
1038  {
1039  for (unsigned int i = 0; i < _geometric_cuts.size(); ++i)
1040  if (_geometric_cuts[i]->shouldHealMesh())
1041  _sibling_displaced_elems[_geometric_cuts[i]->getInterfaceID()].clear();
1042  }
1043 
1044  for (auto & ceh : _crack_tip_elems_to_be_healed)
1045  {
1046  _crack_tip_elems.erase(ceh);
1048  delete _cut_elem_map.find(ceh->unique_id())->second;
1049  _cut_elem_map.erase(ceh->unique_id());
1050  }
1051 
1052  return mesh_changed;
1053 }
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:260
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:265
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:266
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:263
std::set< const Elem * > _crack_tip_elems_to_be_healed
Definition: XFEM.h:264
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:268
bool isPointPhysical(const Point &p) const
Definition: XFEMCutElem.C:191

◆ initCutIntersectionEdge()

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

Definition at line 900 of file XFEM.C.

Referenced by correctCrackExtensionDirection(), and markCutEdgesByState().

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

◆ initSolution()

void XFEM::initSolution ( NonlinearSystemBase &  nl,
AuxiliarySystem &  aux 
)
overridevirtual

Definition at line 312 of file XFEM.C.

313 {
314  nl.serializeSolution();
315  aux.serializeSolution();
316  NumericVector<Number> & current_solution = *nl.system().current_local_solution;
317  NumericVector<Number> & old_solution = nl.solutionOld();
318  NumericVector<Number> & older_solution = nl.solutionOlder();
319  NumericVector<Number> & current_aux_solution = *aux.system().current_local_solution;
320  NumericVector<Number> & old_aux_solution = aux.solutionOld();
321  NumericVector<Number> & older_aux_solution = aux.solutionOlder();
322 
323  setSolution(nl, _cached_solution, current_solution, old_solution, older_solution);
324  setSolution(
325  aux, _cached_aux_solution, current_aux_solution, old_aux_solution, older_aux_solution);
326 
327  current_solution.close();
328  old_solution.close();
329  older_solution.close();
330  current_aux_solution.close();
331  old_aux_solution.close();
332  older_aux_solution.close();
333 
334  _cached_solution.clear();
335  _cached_aux_solution.clear();
336 }
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:1897
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:306
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:297

◆ isElemAtCrackTip()

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

Definition at line 1566 of file XFEM.C.

Referenced by markCutEdgesByGeometry(), and markCutEdgesByState().

1567 {
1568  return (_crack_tip_elems.find(elem) != _crack_tip_elems.end());
1569 }
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:263

◆ isElemCut() [1/2]

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

Definition at line 1572 of file XFEM.C.

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

1573 {
1574  xfce = NULL;
1575  bool is_cut = false;
1576  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1577  it = _cut_elem_map.find(elem->unique_id());
1578  if (it != _cut_elem_map.end())
1579  {
1580  xfce = it->second;
1581  const EFAElement * EFAelem = xfce->getEFAElement();
1582  if (EFAelem->isPartial()) // exclude the full crack tip elements
1583  is_cut = true;
1584  }
1585  return is_cut;
1586 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
virtual bool isPartial() const =0
virtual const EFAElement * getEFAElement() const =0

◆ isElemCut() [2/2]

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

Definition at line 1589 of file XFEM.C.

1590 {
1591  XFEMCutElem * xfce = NULL;
1592  return isElemCut(elem, xfce);
1593 }
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1572

◆ 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 1514 of file XFEM.C.

1515 {
1516  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
1517  it = _cut_elem_map.find(elem->unique_id());
1518  if (it != _cut_elem_map.end())
1519  {
1520  XFEMCutElem * xfce = it->second;
1521 
1522  if (xfce->isPointPhysical(point))
1523  return true;
1524  }
1525  else
1526  return true;
1527 
1528  return false;
1529 }
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
bool isPointPhysical(const Point &p) const
Definition: XFEMCutElem.C:191

◆ markCutEdgesByGeometry()

bool XFEM::markCutEdgesByGeometry ( )

Definition at line 394 of file XFEM.C.

Referenced by markCuts().

395 {
396  bool marked_edges = false;
397  bool marked_nodes = false;
398 
399  for (const auto & gme : _geom_marked_elems_2d)
400  {
401  for (const auto & gmei : gme.second)
402  {
403  EFAElement2D * EFAElem = getEFAElem2D(gme.first);
404 
405  for (unsigned int i = 0; i < gmei._elem_cut_edges.size(); ++i) // mark element edges
406  {
407  if (!EFAElem->isEdgePhantom(
408  gmei._elem_cut_edges[i]._host_side_id)) // must not be phantom edge
409  {
410  _efa_mesh.addElemEdgeIntersection(gme.first->id(),
411  gmei._elem_cut_edges[i]._host_side_id,
412  gmei._elem_cut_edges[i]._distance);
413  marked_edges = true;
414  }
415  }
416 
417  for (unsigned int i = 0; i < gmei._elem_cut_nodes.size(); ++i) // mark element edges
418  {
419  _efa_mesh.addElemNodeIntersection(gme.first->id(), gmei._elem_cut_nodes[i]._host_id);
420  marked_nodes = true;
421  }
422 
423  for (unsigned int i = 0; i < gmei._frag_cut_edges.size();
424  ++i) // MUST DO THIS AFTER MARKING ELEMENT EDGES
425  {
426  if (!EFAElem->getFragment(0)->isSecondaryInteriorEdge(
427  gmei._frag_cut_edges[i]._host_side_id))
428  {
429  if (_efa_mesh.addFragEdgeIntersection(gme.first->id(),
430  gmei._frag_cut_edges[i]._host_side_id,
431  gmei._frag_cut_edges[i]._distance))
432  {
433  marked_edges = true;
434  if (!isElemAtCrackTip(gme.first))
435  _has_secondary_cut = true;
436  }
437  }
438  }
439  }
440  }
441 
442  return marked_edges || marked_nodes;
443 }
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:288
bool isEdgePhantom(unsigned int edge_id) const
bool _has_secondary_cut
Definition: XFEM.h:253
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:277
bool isElemAtCrackTip(const Elem *elem) const
Definition: XFEM.C:1566
EFAElement2D * getEFAElem2D(const Elem *elem)
Get the EFAElement2D object for a specified libMesh element.
Definition: XFEM.C:1613

◆ markCutEdgesByState()

bool XFEM::markCutEdgesByState ( Real  time)

Definition at line 587 of file XFEM.C.

Referenced by markCuts().

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

◆ markCutFacesByGeometry()

bool XFEM::markCutFacesByGeometry ( )

Definition at line 851 of file XFEM.C.

Referenced by markCuts().

852 {
853  bool marked_faces = false;
854 
855  for (const auto & gme : _geom_marked_elems_3d)
856  {
857  for (const auto & gmei : gme.second)
858  {
859  EFAElement3D * EFAElem = getEFAElem3D(gme.first);
860 
861  for (unsigned int i = 0; i < gmei._elem_cut_faces.size(); ++i) // mark element faces
862  {
863  if (!EFAElem->isFacePhantom(gmei._elem_cut_faces[i]._face_id)) // must not be phantom face
864  {
865  _efa_mesh.addElemFaceIntersection(gme.first->id(),
866  gmei._elem_cut_faces[i]._face_id,
867  gmei._elem_cut_faces[i]._face_edge,
868  gmei._elem_cut_faces[i]._position);
869  marked_faces = true;
870  }
871  }
872 
873  for (unsigned int i = 0; i < gmei._frag_cut_faces.size();
874  ++i) // MUST DO THIS AFTER MARKING ELEMENT EDGES
875  {
876  if (!EFAElem->getFragment(0)->isThirdInteriorFace(gmei._frag_cut_faces[i]._face_id))
877  {
878  _efa_mesh.addFragFaceIntersection(gme.first->id(),
879  gmei._frag_cut_faces[i]._face_id,
880  gmei._frag_cut_faces[i]._face_edge,
881  gmei._frag_cut_faces[i]._position);
882  marked_faces = true;
883  }
884  }
885  }
886  }
887 
888  return marked_faces;
889 }
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:280
bool isFacePhantom(unsigned int face_id) const
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:288
void addFragFaceIntersection(unsigned int ElemID, unsigned int FragFaceID, std::vector< unsigned int > FragFaceEdgeID, std::vector< double > position)
void addElemFaceIntersection(unsigned int elemid, unsigned int faceid, std::vector< unsigned int > edgeid, std::vector< double > position)
bool isThirdInteriorFace(unsigned int face_id) const
EFAElement3D * getEFAElem3D(const Elem *elem)
Get the EFAElement3D object for a specified libMesh element.
Definition: XFEM.C:1625

◆ markCutFacesByState()

bool XFEM::markCutFacesByState ( )

Definition at line 892 of file XFEM.C.

Referenced by markCuts().

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

◆ markCuts()

bool XFEM::markCuts ( Real  time)

Definition at line 377 of file XFEM.C.

Referenced by update().

378 {
379  bool marked_sides = false;
380  if (_mesh->mesh_dimension() == 2)
381  {
382  marked_sides = markCutEdgesByGeometry();
383  marked_sides |= markCutEdgesByState(time);
384  }
385  else if (_mesh->mesh_dimension() == 3)
386  {
387  marked_sides = markCutFacesByGeometry();
388  marked_sides |= markCutFacesByState();
389  }
390  return marked_sides;
391 }
bool markCutEdgesByState(Real time)
Definition: XFEM.C:587
bool markCutEdgesByGeometry()
Definition: XFEM.C:394
bool markCutFacesByGeometry()
Definition: XFEM.C:851
bool markCutFacesByState()
Definition: XFEM.C:892

◆ setCrackGrowthMethod()

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

Definition at line 1697 of file XFEM.C.

1698 {
1699  _use_crack_growth_increment = use_crack_growth_increment;
1700  _crack_growth_increment = crack_growth_increment;
1701 }
bool _use_crack_growth_increment
Definition: XFEM.h:257
Real _crack_growth_increment
Definition: XFEM.h:258

◆ 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 1897 of file XFEM.C.

Referenced by initSolution().

1902 {
1903  for (auto & node : _mesh->local_node_ptr_range())
1904  {
1905  auto mit = stored_solution.find(node->unique_id());
1906  if (mit != stored_solution.end())
1907  {
1908  const std::vector<Real> & stored_node_solution = mit->second;
1909  std::vector<dof_id_type> stored_solution_dofs = getNodeSolutionDofs(node, sys);
1910  setSolutionForDOFs(stored_node_solution,
1911  stored_solution_dofs,
1912  current_solution,
1913  old_solution,
1914  older_solution);
1915  }
1916  }
1917 
1918  for (auto & elem : as_range(_mesh->local_elements_begin(), _mesh->local_elements_end()))
1919  {
1920  auto mit = stored_solution.find(elem->unique_id());
1921  if (mit != stored_solution.end())
1922  {
1923  const std::vector<Real> & stored_elem_solution = mit->second;
1924  std::vector<dof_id_type> stored_solution_dofs = getElementSolutionDofs(elem, sys);
1925  setSolutionForDOFs(stored_elem_solution,
1926  stored_solution_dofs,
1927  current_solution,
1928  old_solution,
1929  older_solution);
1930  }
1931  }
1932 }
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:1935
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:1984
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:1958

◆ 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 1935 of file XFEM.C.

Referenced by setSolution().

1940 {
1941  // Solution vector is stored first for current, then old and older solutions.
1942  // These are the offsets to the beginning of the old and older solutions in the vector.
1943  const auto old_solution_offset = stored_solution_dofs.size();
1944  const auto older_solution_offset = old_solution_offset * 2;
1945 
1946  for (std::size_t i = 0; i < stored_solution_dofs.size(); ++i)
1947  {
1948  current_solution.set(stored_solution_dofs[i], stored_solution[i]);
1949  if (_fe_problem->isTransient())
1950  {
1951  old_solution.set(stored_solution_dofs[i], stored_solution[old_solution_offset + i]);
1952  older_solution.set(stored_solution_dofs[i], stored_solution[older_solution_offset + i]);
1953  }
1954  }
1955 }

◆ setXFEMQRule()

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

Definition at line 1686 of file XFEM.C.

1687 {
1688  if (xfem_qrule == "volfrac")
1690  else if (xfem_qrule == "moment_fitting")
1692  else if (xfem_qrule == "direct")
1694 }
Xfem::XFEM_QRULE _XFEM_qrule
Definition: XFEM.h:255

◆ storeCrackTipOriginAndDirection()

void XFEM::storeCrackTipOriginAndDirection ( )

Definition at line 188 of file XFEM.C.

Referenced by update().

189 {
191  std::set<EFAElement *> CrackTipElements = _efa_mesh.getCrackTipElements();
192  std::set<EFAElement *>::iterator sit;
193  for (sit = CrackTipElements.begin(); sit != CrackTipElements.end(); ++sit)
194  {
195  if (_mesh->mesh_dimension() == 2)
196  {
197  EFAElement2D * CEMElem = dynamic_cast<EFAElement2D *>(*sit);
198  EFANode * tip_node = CEMElem->getTipEmbeddedNode();
199  unsigned int cts_id = CEMElem->getCrackTipSplitElementID();
200 
201  Point origin(0, 0, 0);
202  Point direction(0, 0, 0);
203 
204  std::map<unique_id_type, XFEMCutElem *>::const_iterator it;
205  it = _cut_elem_map.find(_mesh->elem(cts_id)->unique_id());
206  if (it != _cut_elem_map.end())
207  {
208  const XFEMCutElem * xfce = it->second;
209  const EFAElement * EFAelem = xfce->getEFAElement();
210  if (EFAelem->isPartial()) // exclude the full crack tip elements
211  {
212  xfce->getCrackTipOriginAndDirection(tip_node->id(), origin, direction);
213  }
214  }
215 
216  std::vector<Point> tip_data;
217  tip_data.push_back(origin);
218  tip_data.push_back(direction);
219  const Elem * elem = _mesh->elem((*sit)->id());
221  std::pair<const Elem *, std::vector<Point>>(elem, tip_data));
222  }
223  }
224 }
unsigned int getCrackTipSplitElementID() const
Definition: EFAElement2D.C:594
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:262
virtual bool isPartial() const =0
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:288
virtual const EFAElement * getEFAElement() const =0
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:268
EFANode * getTipEmbeddedNode() const
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

◆ 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 1863 of file XFEM.C.

Referenced by cutMeshWithEFA().

1870 {
1871  std::vector<dof_id_type> stored_solution_dofs = getElementSolutionDofs(elem_to_store_from, sys);
1872  std::vector<Real> stored_solution_scratch;
1873  // Size for current solution, as well as for old, and older solution only for transient case
1874  std::size_t stored_solution_size =
1875  (_fe_problem->isTransient() ? stored_solution_dofs.size() * 3 : stored_solution_dofs.size());
1876  stored_solution_scratch.reserve(stored_solution_size);
1877 
1878  // Store in the order defined in stored_solution_dofs first for the current, then for old and
1879  // older if applicable
1880  for (auto dof : stored_solution_dofs)
1881  stored_solution_scratch.push_back(current_solution(dof));
1882 
1883  if (_fe_problem->isTransient())
1884  {
1885  for (auto dof : stored_solution_dofs)
1886  stored_solution_scratch.push_back(old_solution(dof));
1887 
1888  for (auto dof : stored_solution_dofs)
1889  stored_solution_scratch.push_back(older_solution(dof));
1890  }
1891 
1892  if (stored_solution_scratch.size() > 0)
1893  stored_solution[elem_to_store_to->unique_id()] = stored_solution_scratch;
1894 }
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:1958

◆ 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 1829 of file XFEM.C.

Referenced by cutMeshWithEFA().

1836 {
1837  std::vector<dof_id_type> stored_solution_dofs = getNodeSolutionDofs(node_to_store_from, sys);
1838  std::vector<Real> stored_solution_scratch;
1839  // Size for current solution, as well as for old, and older solution only for transient case
1840  std::size_t stored_solution_size =
1841  (_fe_problem->isTransient() ? stored_solution_dofs.size() * 3 : stored_solution_dofs.size());
1842  stored_solution_scratch.reserve(stored_solution_size);
1843 
1844  // Store in the order defined in stored_solution_dofs first for the current, then for old and
1845  // older if applicable
1846  for (auto dof : stored_solution_dofs)
1847  stored_solution_scratch.push_back(current_solution(dof));
1848 
1849  if (_fe_problem->isTransient())
1850  {
1851  for (auto dof : stored_solution_dofs)
1852  stored_solution_scratch.push_back(old_solution(dof));
1853 
1854  for (auto dof : stored_solution_dofs)
1855  stored_solution_scratch.push_back(older_solution(dof));
1856  }
1857 
1858  if (stored_solution_scratch.size() > 0)
1859  stored_solution[node_to_store_to->unique_id()] = stored_solution_scratch;
1860 }
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:1984

◆ update()

bool XFEM::update ( Real  time,
NonlinearSystemBase &  nl,
AuxiliarySystem &  aux 
)
overridevirtual

Definition at line 264 of file XFEM.C.

265 {
266  if (_moose_mesh->isDistributedMesh())
267  mooseError("Use of XFEM with distributed mesh is not yet supported");
268 
269  bool mesh_changed = false;
270 
271  buildEFAMesh();
272 
273  _fe_problem->execute(EXEC_XFEM_MARK);
274 
276 
277  if (markCuts(time))
278  mesh_changed = cutMeshWithEFA(nl, aux);
279 
280  if (mesh_changed)
281  {
282  buildEFAMesh();
284  }
285 
286  if (mesh_changed)
287  {
288  // _mesh->find_neighbors();
289  // _mesh->contract();
290  _mesh->allow_renumbering(false);
291  _mesh->skip_partitioning(true);
292  _mesh->prepare_for_use();
293  // _mesh->prepare_for_use(true,true); //doing this preserves the numbering, but generates
294  // warning
295 
296  if (_displaced_mesh)
297  {
298  _displaced_mesh->allow_renumbering(false);
299  _displaced_mesh->skip_partitioning(true);
300  _displaced_mesh->prepare_for_use();
301  // _displaced_mesh->prepare_for_use(true,true);
302  }
303  }
304 
307 
308  return mesh_changed;
309 }
void clearGeomMarkedElems()
Clear out the list of elements to be marked for cutting.
Definition: XFEM.C:181
void storeCrackTipOriginAndDirection()
Definition: XFEM.C:188
bool markCuts(Real time)
Definition: XFEM.C:377
const ExecFlagType EXEC_XFEM_MARK
Exec flag used to execute MooseObjects while elements are being marked for cutting by XFEM...
bool cutMeshWithEFA(NonlinearSystemBase &nl, AuxiliarySystem &aux)
Definition: XFEM.C:1056
void clearStateMarkedElems()
Definition: XFEM.C:153
void buildEFAMesh()
Definition: XFEM.C:339

◆ updateHeal()

bool XFEM::updateHeal ( )
overridevirtual

Definition at line 227 of file XFEM.C.

228 {
229  bool mesh_changed = false;
230 
231  mesh_changed = healMesh();
232 
233  if (mesh_changed)
234  {
235  _mesh->update_parallel_id_counts();
236  MeshCommunication().make_elems_parallel_consistent(*_mesh);
237  MeshCommunication().make_nodes_parallel_consistent(*_mesh);
238  // _mesh->find_neighbors();
239  // _mesh->contract();
240  _mesh->allow_renumbering(false);
241  _mesh->skip_partitioning(true);
242  _mesh->prepare_for_use();
243  // _mesh->prepare_for_use(true,true); //doing this preserves the numbering, but generates
244  // warning
245 
246  if (_displaced_mesh)
247  {
248  _displaced_mesh->update_parallel_id_counts();
249  MeshCommunication().make_elems_parallel_consistent(*_displaced_mesh);
250  MeshCommunication().make_nodes_parallel_consistent(*_displaced_mesh);
251  _displaced_mesh->allow_renumbering(false);
252  _displaced_mesh->skip_partitioning(true);
253  _displaced_mesh->prepare_for_use();
254  // _displaced_mesh->prepare_for_use(true,true);
255  }
256  }
257 
258  _geom_marker_id_elems.clear();
259 
260  return mesh_changed;
261 }
bool healMesh()
Potentially heal the mesh by merging some of the pairs of partial elements cut by XFEM back into sing...
Definition: XFEM.C:919
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:283

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 306 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 297 of file XFEM.h.

Referenced by cutMeshWithEFA(), and initSolution().

◆ _crack_growth_increment

Real XFEM::_crack_growth_increment
private

Definition at line 258 of file XFEM.h.

Referenced by markCutEdgesByState(), and setCrackGrowthMethod().

◆ _crack_tip_elems

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

Definition at line 263 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 264 of file XFEM.h.

Referenced by cutMeshWithEFA(), and healMesh().

◆ _cut_elem_map

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

◆ _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_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 277 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 280 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 283 of file XFEM.h.

Referenced by addGeomMarkedElem2D(), addGeomMarkedElem3D(), cutMeshWithEFA(), 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 286 of file XFEM.h.

Referenced by addGeometricCut(), and getGeometricCutID().

◆ _geometric_cuts

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

Definition at line 260 of file XFEM.h.

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

◆ _has_secondary_cut

bool XFEM::_has_secondary_cut
private

Definition at line 253 of file XFEM.h.

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

◆ _sibling_displaced_elems

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

Definition at line 266 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 265 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 274 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 272 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 273 of file XFEM.h.

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

◆ _use_crack_growth_increment

bool XFEM::_use_crack_growth_increment
private

Definition at line 257 of file XFEM.h.

Referenced by markCutEdgesByState(), and setCrackGrowthMethod().

◆ _XFEM_qrule

Xfem::XFEM_QRULE XFEM::_XFEM_qrule
private

Definition at line 255 of file XFEM.h.

Referenced by getXFEMQRule(), and setXFEMQRule().


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