www.mooseframework.org
XFEM.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "ElementPairLocator.h"
14 #include "XFEMInterface.h"
16 
17 #include "libmesh/vector_value.h"
18 #include "libmesh/quadrature.h"
19 
20 #include "GeometricCutUserObject.h"
21 
22 // Forward declarations
23 class SystemBase;
24 
25 namespace Xfem
26 {
28 {
35 };
36 
38 {
42 };
43 
50 typedef std::array<std::unordered_map<unsigned int, std::string>, 2> CachedMaterialProperties;
51 
58 {
59  const Elem * _parent_elem;
64  // TODO: add neighbor material properties
65  // CachedMaterialProperties _neighbor_material_properties;
66 
68  : _parent_elem(nullptr),
69  _geometric_cut(nullptr),
70  _cut_subdomain_id(std::numeric_limits<CutSubdomainID>::max())
71  {
72  }
73 
74  CutElemInfo(const Elem * parent_elem,
75  const GeometricCutUserObject * geometric_cut,
76  CutSubdomainID cut_subdomain_id)
77  : _parent_elem(parent_elem), _geometric_cut(geometric_cut), _cut_subdomain_id(cut_subdomain_id)
78  {
79  }
80 
81  // A new child element is said to be previously healed if an entry in the old CutElemInfo has the
82  // same parent element ID, the same cut, AND the same cut subdomain ID.
83  bool match(const CutElemInfo & rhs)
84  {
85  return _parent_elem == rhs._parent_elem && _geometric_cut == rhs._geometric_cut &&
87  }
88 };
89 } // namespace Xfem
90 
91 class XFEMCutElem;
93 class EFANode;
94 class EFAEdge;
95 class EFAElement;
96 class EFAElement2D;
97 class EFAElement3D;
98 
105 // ------------------------------------------------------------
106 // XFEM class definition
107 class XFEM : public XFEMInterface
108 {
109 public:
110  explicit XFEM(const InputParameters & params);
111 
112  ~XFEM();
113 
114  void addGeometricCut(GeometricCutUserObject * geometric_cut);
115 
116  void addStateMarkedElem(unsigned int elem_id, RealVectorValue & normal);
117  void addStateMarkedElem(unsigned int elem_id, RealVectorValue & normal, unsigned int marked_side);
118  void addStateMarkedFrag(unsigned int elem_id, RealVectorValue & normal);
119 
120  void clearStateMarkedElems();
121 
128  void addGeomMarkedElem2D(const unsigned int elem_id,
129  const Xfem::GeomMarkedElemInfo2D geom_info,
130  const unsigned int interface_id);
131 
138  void addGeomMarkedElem3D(const unsigned int elem_id,
139  const Xfem::GeomMarkedElemInfo3D geom_info,
140  const unsigned int interface_id);
141 
145  void clearGeomMarkedElems();
146 
147  virtual bool update(Real time,
148  const std::vector<std::shared_ptr<NonlinearSystemBase>> & nl,
149  AuxiliarySystem & aux) override;
150 
151  virtual void initSolution(const std::vector<std::shared_ptr<NonlinearSystemBase>> & nl,
152  AuxiliarySystem & aux) override;
153 
154  void buildEFAMesh();
155  bool markCuts(Real time);
156  bool markCutEdgesByGeometry();
157  bool markCutEdgesByState(Real time);
158  bool markCutFacesByGeometry();
159  bool markCutFacesByState();
161  Point cut_origin, RealVectorValue cut_normal, Point & edge_p1, Point & edge_p2, Real & dist);
162  bool cutMeshWithEFA(const std::vector<std::shared_ptr<NonlinearSystemBase>> & nl,
163  AuxiliarySystem & aux);
164 
171  bool healMesh();
172 
173  virtual bool updateHeal() override;
174  Point getEFANodeCoords(EFANode * CEMnode,
175  EFAElement * CEMElem,
176  const Elem * elem,
177  MeshBase * displaced_mesh = nullptr) const;
178 
182  Real getPhysicalVolumeFraction(const Elem * elem) const;
183 
188  bool isPointInsidePhysicalDomain(const Elem * elem, const Point & point) const;
189 
193  Real getCutPlane(const Elem * elem,
194  const Xfem::XFEM_CUTPLANE_QUANTITY quantity,
195  unsigned int plane_id) const;
196 
197  bool isElemAtCrackTip(const Elem * elem) const;
198  bool isElemCut(const Elem * elem, XFEMCutElem *& xfce) const;
199  bool isElemCut(const Elem * elem) const;
200  void getFragmentFaces(const Elem * elem,
201  std::vector<std::vector<Point>> & frag_faces,
202  bool displaced_mesh = false) const;
204 
205  void correctCrackExtensionDirection(const Elem * elem,
206  EFAElement2D * CEMElem,
207  EFAEdge * orig_edge,
208  Point normal,
209  Point crack_tip_origin,
210  Point crack_tip_direction,
211  Real & distance_keep,
212  unsigned int & edge_id_keep,
213  Point & normal_keep);
214 
215  void getCrackTipOrigin(std::map<unsigned int, const Elem *> & elem_id_crack_tip,
216  std::vector<Point> & crack_front_points);
217 
219  void setXFEMQRule(std::string & xfem_qrule);
220  void setCrackGrowthMethod(bool use_crack_growth_increment, Real crack_growth_increment);
221 
227  void setDebugOutputLevel(unsigned int debug_output_level);
228 
229  virtual bool getXFEMWeights(MooseArray<Real> & weights,
230  const Elem * elem,
231  QBase * qrule,
232  const MooseArray<Point> & q_points) override;
233  virtual bool getXFEMFaceWeights(MooseArray<Real> & weights,
234  const Elem * elem,
235  QBase * qrule,
236  const MooseArray<Point> & q_points,
237  unsigned int side) override;
238 
245  virtual const ElementPairLocator::ElementPairList * getXFEMCutElemPairs(unsigned int interface_id)
246  {
247  return &_sibling_elems[interface_id];
248  }
249 
257  getXFEMDisplacedCutElemPairs(unsigned int interface_id)
258  {
259  return &_sibling_displaced_elems[interface_id];
260  }
261 
267  virtual unsigned int getGeometricCutID(const GeometricCutUserObject * gcu)
268  {
269  return _geom_marker_id_map[gcu];
270  };
271 
272  virtual void getXFEMIntersectionInfo(const Elem * elem,
273  unsigned int plane_id,
274  Point & normal,
275  std::vector<Point> & intersectionPoints,
276  bool displaced_mesh = false) const;
277  virtual void getXFEMqRuleOnLine(std::vector<Point> & intersection_points,
278  std::vector<Point> & quad_pts,
279  std::vector<Real> & quad_wts) const;
280  virtual void getXFEMqRuleOnSurface(std::vector<Point> & intersection_points,
281  std::vector<Point> & quad_pts,
282  std::vector<Real> & quad_wts) const;
284 
289  EFAElement2D * getEFAElem2D(const Elem * elem);
290 
295  EFAElement3D * getEFAElem3D(const Elem * elem);
296 
297  void getFragmentEdges(const Elem * elem,
298  EFAElement2D * CEMElem,
299  std::vector<std::vector<Point>> & frag_edges) const;
300  void getFragmentFaces(const Elem * elem,
301  EFAElement3D * CEMElem,
302  std::vector<std::vector<Point>> & frag_faces) const;
303 
304  const std::map<const Elem *, std::vector<Point>> & getCrackTipOriginMap() const
305  {
307  }
308 
316  const Elem * cut_elem,
317  const Elem * parent_elem = nullptr) const;
318 
319 private:
321 
323 
326 
327  std::vector<const GeometricCutUserObject *> _geometric_cuts;
328 
329  std::map<unique_id_type, XFEMCutElem *> _cut_elem_map;
330  std::set<const Elem *> _crack_tip_elems;
331  std::set<const Elem *> _crack_tip_elems_to_be_healed;
332  std::map<unsigned int, ElementPairLocator::ElementPairList> _sibling_elems;
333  std::map<unsigned int, ElementPairLocator::ElementPairList> _sibling_displaced_elems;
334 
335  std::map<const Elem *, std::vector<Point>> _elem_crack_origin_direction_map;
336 
337  // std::map<const Elem*, Point> _crack_propagation_direction_map;
338 
339  std::map<const Elem *, RealVectorValue> _state_marked_elems;
340  std::set<const Elem *> _state_marked_frags;
341  std::map<const Elem *, unsigned int> _state_marked_elem_sides;
342 
344  std::map<const Elem *, std::vector<Xfem::GeomMarkedElemInfo2D>> _geom_marked_elems_2d;
345 
347  std::map<const Elem *, std::vector<Xfem::GeomMarkedElemInfo3D>> _geom_marked_elems_3d;
348 
350  std::map<unsigned int, std::set<unsigned int>> _geom_marker_id_elems;
351 
353  std::map<const GeometricCutUserObject *, unsigned int> _geom_marker_id_map;
354 
356 
362  unsigned int _debug_output_level;
363 
371  std::map<unique_id_type, std::vector<Real>> _cached_solution;
372 
380  std::map<unique_id_type, std::vector<Real>> _cached_aux_solution;
381 
386  std::unordered_map<const Elem *, Xfem::CutElemInfo> _geom_cut_elems;
387 
392  std::unordered_map<const Elem *, Xfem::CutElemInfo> _old_geom_cut_elems;
393 
404  void storeSolutionForNode(const Node * node_to_store_to,
405  const Node * node_to_store_from,
406  SystemBase & sys,
407  std::map<unique_id_type, std::vector<Real>> & stored_solution,
408  const NumericVector<Number> & current_solution,
409  const NumericVector<Number> & old_solution,
410  const NumericVector<Number> & older_solution);
411 
422  void storeSolutionForElement(const Elem * elem_to_store_to,
423  const Elem * elem_to_store_from,
424  SystemBase & sys,
425  std::map<unique_id_type, std::vector<Real>> & stored_solution,
426  const NumericVector<Number> & current_solution,
427  const NumericVector<Number> & old_solution,
428  const NumericVector<Number> & older_solution);
429 
438  void setSolution(SystemBase & sys,
439  const std::map<unique_id_type, std::vector<Real>> & stored_solution,
440  NumericVector<Number> & current_solution,
441  NumericVector<Number> & old_solution,
442  NumericVector<Number> & older_solution);
443 
452  void setSolutionForDOFs(const std::vector<Real> & stored_solution,
453  const std::vector<dof_id_type> & stored_solution_dofs,
454  NumericVector<Number> & current_solution,
455  NumericVector<Number> & old_solution,
456  NumericVector<Number> & older_solution);
457 
464  std::vector<dof_id_type> getElementSolutionDofs(const Elem * elem, SystemBase & sys) const;
465 
472  std::vector<dof_id_type> getNodeSolutionDofs(const Node * node, SystemBase & sys) const;
473 
479  const GeometricCutUserObject * getGeometricCutForElem(const Elem * elem) const;
480 
481  void storeMaterialPropertiesForElementHelper(const Elem * elem,
482  MaterialPropertyStorage & storage);
483 
490  void storeMaterialPropertiesForElement(const Elem * parent_elem, const Elem * child_elem);
491 
500  void loadMaterialPropertiesForElementHelper(const Elem * elem,
501  const Xfem::CachedMaterialProperties & cached_props,
502  MaterialPropertyStorage & storage) const;
503 
511  const Elem * elem,
512  const Elem * elem_from,
513  std::unordered_map<const Elem *, Xfem::CutElemInfo> & cached_cei) const;
514 
521  const Node * pickFirstPhysicalNode(const Elem * e, const Elem * e0) const;
522 };
~XFEM()
Definition: XFEM.C:46
void getCrackTipOrigin(std::map< unsigned int, const Elem *> &elem_id_crack_tip, std::vector< Point > &crack_front_points)
Definition: XFEM.C:65
void getFragmentFaces(const Elem *elem, std::vector< std::vector< Point >> &frag_faces, bool displaced_mesh=false) const
Definition: XFEM.C:1729
void correctCrackExtensionDirection(const Elem *elem, EFAElement2D *CEMElem, EFAEdge *orig_edge, Point normal, Point crack_tip_origin, Point crack_tip_direction, Real &distance_keep, unsigned int &edge_id_keep, Point &normal_keep)
Definition: XFEM.C:444
bool _use_crack_growth_increment
Definition: XFEM.h:324
unsigned int _debug_output_level
Controls amount of debugging output information 0: None 1: Summary 2: Details on modifications to mes...
Definition: XFEM.h:362
void setSolutionForDOFs(const std::vector< Real > &stored_solution, const std::vector< dof_id_type > &stored_solution_dofs, NumericVector< Number > &current_solution, NumericVector< Number > &old_solution, NumericVector< Number > &older_solution)
Set the solution for a set of DOFs.
Definition: XFEM.C:2074
const GeometricCutUserObject * getGeometricCutForElem(const Elem *elem) const
Get the GeometricCutUserObject associated with an element.
Definition: XFEM.C:2155
bool markCutEdgesByState(Real time)
Definition: XFEM.C:585
void storeMaterialPropertiesForElement(const Elem *parent_elem, const Elem *child_elem)
Helper function to store the material properties of a healed element.
Definition: XFEM.C:2185
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo3D > > _geom_marked_elems_3d
Data structure for storing information about all 3D elements to be cut by geometry.
Definition: XFEM.h:347
CutSubdomainID getCutSubdomainID(const GeometricCutUserObject *gcuo, const Elem *cut_elem, const Elem *parent_elem=nullptr) const
Determine which cut subdomain the element belongs to relative to the cut.
Definition: XFEM.C:2271
Data structure describing geometrically described cut through 3D element.
XFEM_QRULE
Definition: XFEM.h:37
std::vector< const GeometricCutUserObject * > _geometric_cuts
Definition: XFEM.h:327
void setSolution(SystemBase &sys, const std::map< unique_id_type, std::vector< Real >> &stored_solution, NumericVector< Number > &current_solution, NumericVector< Number > &old_solution, NumericVector< Number > &older_solution)
Set the solution for all locally-owned nodes/elements that have stored values.
Definition: XFEM.C:2036
virtual bool update(Real time, const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux) override
Definition: XFEM.C:260
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_elems
Definition: XFEM.h:332
Real _crack_growth_increment
Definition: XFEM.h:325
virtual bool getXFEMFaceWeights(MooseArray< Real > &weights, const Elem *elem, QBase *qrule, const MooseArray< Point > &q_points, unsigned int side) override
Definition: XFEM.C: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:2123
bool isPointInsidePhysicalDomain(const Elem *elem, const Point &point) const
Return true if the point is inside the element physical domain Note: if this element is not cut...
Definition: XFEM.C:1648
void getFragmentEdges(const Elem *elem, EFAElement2D *CEMElem, std::vector< std::vector< Point >> &frag_edges) const
Definition: XFEM.C:1770
bool match(const CutElemInfo &rhs)
Definition: XFEM.h:83
bool cutMeshWithEFA(const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux)
Definition: XFEM.C:1101
const std::map< const Elem *, std::vector< Point > > & getCrackTipOriginMap() const
Definition: XFEM.h:304
std::map< unsigned int, ElementPairLocator::ElementPairList > _sibling_displaced_elems
Definition: XFEM.h:333
void clearGeomMarkedElems()
Clear out the list of elements to be marked for cutting.
Definition: XFEM.C:177
std::map< unique_id_type, XFEMCutElem * > _cut_elem_map
Definition: XFEM.h:329
std::set< const Elem * > _crack_tip_elems
Definition: XFEM.h:330
CutSubdomainID _cut_subdomain_id
Definition: XFEM.h:61
bool healMesh()
Potentially heal the mesh by merging some of the pairs of partial elements cut by XFEM back into sing...
Definition: XFEM.C:917
XFEM(const InputParameters &params)
Definition: XFEM.C:36
ElementFragmentAlgorithm _efa_mesh
Definition: XFEM.h:355
std::array< std::unordered_map< unsigned int, std::string >, 2 > CachedMaterialProperties
Convenient typedef for local storage of stateful material properties.
Definition: XFEM.h:50
const GeometricCutUserObject * _geometric_cut
Definition: XFEM.h:60
CachedMaterialProperties _bnd_material_properties
Definition: XFEM.h:63
This is the XFEM class.
Definition: XFEM.h:107
CutElemInfo(const Elem *parent_elem, const GeometricCutUserObject *geometric_cut, CutSubdomainID cut_subdomain_id)
Definition: XFEM.h:74
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...
Definition: XFEM.h:257
std::set< const Elem * > _crack_tip_elems_to_be_healed
Definition: XFEM.h:331
bool _has_secondary_cut
Definition: XFEM.h:320
void addGeometricCut(GeometricCutUserObject *geometric_cut)
Definition: XFEM.C:55
auto max(const L &left, const R &right)
std::map< unique_id_type, std::vector< Real > > _cached_aux_solution
Data structure to store the auxiliary solution for nodes/elements affected by XFEM For each node/elem...
Definition: XFEM.h:380
bool markCutEdgesByGeometry()
Definition: XFEM.C:392
unsigned int CutSubdomainID
Definition: XFEMAppTypes.h:18
const Elem * _parent_elem
Definition: XFEM.h:59
virtual const ElementPairLocator::ElementPairList * getXFEMCutElemPairs(unsigned int interface_id)
Get the list of cut element pairs corresponding to a given interface ID.
Definition: XFEM.h:245
bool initCutIntersectionEdge(Point cut_origin, RealVectorValue cut_normal, Point &edge_p1, Point &edge_p2, Real &dist)
Definition: XFEM.C:898
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.
Definition: XFEM.C:167
void storeCrackTipOriginAndDirection()
Definition: XFEM.C:184
const Node * pickFirstPhysicalNode(const Elem *e, const Elem *e0) const
Return the first node in the provided element that is found to be in the physical domain...
Definition: XFEM.C:2284
void loadMaterialPropertiesForElement(const Elem *elem, const Elem *elem_from, std::unordered_map< const Elem *, Xfem::CutElemInfo > &cached_cei) const
Helper function to store the material properties of a healed element.
Definition: XFEM.C:2237
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.
Definition: XFEM.C:1666
bool markCuts(Real time)
Definition: XFEM.C:375
Definition: XFEM.h:25
std::map< const GeometricCutUserObject *, unsigned int > _geom_marker_id_map
Data structure for storing the GeommetricCutUserObjects and their corresponding id.
Definition: XFEM.h:353
void 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.
Definition: XFEM.C:157
void storeSolutionForElement(const Elem *elem_to_store_to, const Elem *elem_to_store_from, SystemBase &sys, std::map< unique_id_type, std::vector< Real >> &stored_solution, const NumericVector< Number > &current_solution, const NumericVector< Number > &old_solution, const NumericVector< Number > &older_solution)
Store the solution in stored_solution for a given element.
Definition: XFEM.C:2002
Xfem::XFEM_QRULE _XFEM_qrule
Definition: XFEM.h:322
void storeSolutionForNode(const Node *node_to_store_to, const Node *node_to_store_from, SystemBase &sys, std::map< unique_id_type, std::vector< Real >> &stored_solution, const NumericVector< Number > &current_solution, const NumericVector< Number > &old_solution, const NumericVector< Number > &older_solution)
Store the solution in stored_solution for a given node.
Definition: XFEM.C:1968
std::list< std::pair< const Elem *, const Elem *> > ElementPairList
virtual void getXFEMqRuleOnLine(std::vector< Point > &intersection_points, std::vector< Point > &quad_pts, std::vector< Real > &quad_wts) const
Definition: XFEM.C:1897
Data structure describing geometrically described cut through 2D element.
void loadMaterialPropertiesForElementHelper(const Elem *elem, const Xfem::CachedMaterialProperties &cached_props, MaterialPropertyStorage &storage) const
Load the material properties.
Definition: XFEM.C:2213
void setDebugOutputLevel(unsigned int debug_output_level)
Controls amount of debugging information output.
Definition: XFEM.C:1837
std::map< const Elem *, std::vector< Point > > _elem_crack_origin_direction_map
Definition: XFEM.h:335
virtual void getXFEMIntersectionInfo(const Elem *elem, unsigned int plane_id, Point &normal, std::vector< Point > &intersectionPoints, bool displaced_mesh=false) const
Definition: XFEM.C:1878
std::map< const Elem *, std::vector< Xfem::GeomMarkedElemInfo2D > > _geom_marked_elems_2d
Data structure for storing information about all 2D elements to be cut by geometry.
Definition: XFEM.h:344
std::map< unsigned int, std::set< unsigned int > > _geom_marker_id_elems
Data structure for storing the elements cut by specific geometric cutters.
Definition: XFEM.h:350
virtual unsigned int getGeometricCutID(const GeometricCutUserObject *gcu)
Get the interface ID corresponding to a given GeometricCutUserObject.
Definition: XFEM.h:267
XFEM_CUTPLANE_QUANTITY
Definition: XFEM.h:27
virtual void initSolution(const std::vector< std::shared_ptr< NonlinearSystemBase >> &nl, AuxiliarySystem &aux) override
Definition: XFEM.C:305
std::set< const Elem * > _state_marked_frags
Definition: XFEM.h:340
void storeMaterialPropertiesForElementHelper(const Elem *elem, MaterialPropertyStorage &storage)
Definition: XFEM.C:2167
EFAElement3D * getEFAElem3D(const Elem *elem)
Get the EFAElement3D object for a specified libMesh element.
Definition: XFEM.C:1758
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool isElemAtCrackTip(const Elem *elem) const
Definition: XFEM.C:1700
bool has_secondary_cut()
Definition: XFEM.h:283
CachedMaterialProperties _elem_material_properties
Definition: XFEM.h:62
std::vector< dof_id_type > getElementSolutionDofs(const Elem *elem, SystemBase &sys) const
Get a vector of the dof indices for all components of all variables associated with an element...
Definition: XFEM.C:2097
Xfem::XFEM_QRULE & getXFEMQRule()
Definition: XFEM.C:1813
Real getPhysicalVolumeFraction(const Elem *elem) const
Get the volume fraction of an element that is physical.
Definition: XFEM.C:1628
virtual bool getXFEMWeights(MooseArray< Real > &weights, const Elem *elem, QBase *qrule, const MooseArray< Point > &q_points) override
Definition: XFEM.C:1843
virtual void getXFEMqRuleOnSurface(std::vector< Point > &intersection_points, std::vector< Point > &quad_pts, std::vector< Real > &quad_wts) const
Definition: XFEM.C:1924
bool markCutFacesByGeometry()
Definition: XFEM.C:849
void setCrackGrowthMethod(bool use_crack_growth_increment, Real crack_growth_increment)
Definition: XFEM.C:1830
void addStateMarkedElem(unsigned int elem_id, RealVectorValue &normal)
Definition: XFEM.C:111
std::map< unique_id_type, std::vector< Real > > _cached_solution
Data structure to store the nonlinear solution for nodes/elements affected by XFEM For each node/elem...
Definition: XFEM.h:371
bool isElemCut(const Elem *elem, XFEMCutElem *&xfce) const
Definition: XFEM.C:1706
virtual bool updateHeal() override
Definition: XFEM.C:223
void clearStateMarkedElems()
Definition: XFEM.C:149
void setXFEMQRule(std::string &xfem_qrule)
Definition: XFEM.C:1819
void buildEFAMesh()
Definition: XFEM.C:336
std::map< const Elem *, unsigned int > _state_marked_elem_sides
Definition: XFEM.h:341
Information about a cut element.
Definition: XFEM.h:57
uint8_t unique_id_type
Point getEFANodeCoords(EFANode *CEMnode, EFAElement *CEMElem, const Elem *elem, MeshBase *displaced_mesh=nullptr) const
Definition: XFEM.C:1596
std::map< const Elem *, RealVectorValue > _state_marked_elems
Definition: XFEM.h:339
EFAElement2D * getEFAElem2D(const Elem *elem)
Get the EFAElement2D object for a specified libMesh element.
Definition: XFEM.C:1746
bool markCutFacesByState()
Definition: XFEM.C:890
void addStateMarkedFrag(unsigned int elem_id, RealVectorValue &normal)
Definition: XFEM.C:135
std::unordered_map< const Elem *, Xfem::CutElemInfo > _geom_cut_elems
All geometrically cut elements and their CutElemInfo during the current execution of XFEM_MARK...
Definition: XFEM.h:386
std::unordered_map< const Elem *, Xfem::CutElemInfo > _old_geom_cut_elems
All geometrically cut elements and their CutElemInfo before the current execution of XFEM_MARK...
Definition: XFEM.h:392