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

#include <CircleCutUserObject.h>

Inheritance diagram for CircleCutUserObject:
[legend]

Public Member Functions

 CircleCutUserObject (const InputParameters &parameters)
 
virtual const std::vector< Point > getCrackFrontPoints (unsigned int num_crack_front_points) const override
 get a set of points along a crack front from a XFEM GeometricCutUserObject More...
 
virtual bool cutElementByGeometry (const Elem *elem, std::vector< Xfem::CutEdge > &cut_edges, std::vector< Xfem::CutNode > &cut_nodes, Real time) const override
 Check to see whether a specified 2D element should be cut based on geometric conditions. More...
 
virtual bool cutElementByGeometry (const Elem *elem, std::vector< Xfem::CutFace > &cut_faces, Real time) const override
 Check to see whether a specified 3D element should be cut based on geometric conditions. More...
 
virtual bool cutFragmentByGeometry (std::vector< std::vector< Point >> &frag_edges, std::vector< Xfem::CutEdge > &cut_edges, Real time) const override
 Check to see whether a fragment of a 2D element should be cut based on geometric conditions. More...
 
virtual bool cutFragmentByGeometry (std::vector< std::vector< Point >> &frag_faces, std::vector< Xfem::CutFace > &cut_faces, Real time) const override
 Check to see whether a fragment of a 3D element should be cut based on geometric conditions. More...
 
virtual void initialize () override
 
virtual void execute () override
 
virtual void threadJoin (const UserObject &y) override
 
virtual void finalize () override
 
unsigned int getInterfaceID () const
 Get the interface ID for this cutting object. More...
 
void setInterfaceID (unsigned int interface_id)
 Set the interface ID for this cutting object. More...
 
bool shouldHealMesh () const
 Should the elements cut by this cutting object be healed in the current time step? More...
 

Protected Member Functions

virtual bool intersectWithEdge (const Point &p1, const Point &p2, Point &pint) const
 
bool isInsideEdge (const Point &p1, const Point &p2, const Point &p) const
 
Real getRelativePosition (const Point &p1, const Point &p2, const Point &p) const
 
void serialize (std::string &serialized_buffer)
 Methods to pack/unpack the _marked_elems_2d and _marked_elems_3d data into a structure suitable for parallel communication. More...
 
void deserialize (std::vector< std::string > &serialized_buffers)
 

Protected Attributes

std::vector< Real > _cut_data
 
Point _center
 
Point _normal
 
std::shared_ptr< XFEM_xfem
 Pointer to the XFEM controller object. More...
 
unsigned int _interface_id
 Associated interface id. More...
 
bool _heal_always
 Heal the mesh. More...
 
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
 Containers with information about all 2D and 3D elements marked for cutting by this object. More...
 
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d
 

Private Member Functions

virtual bool isInsideCutPlane (Point p) const override
 

Private Attributes

std::vector< Point > _vertices
 
Real _radius
 
Real _angle
 

Detailed Description

Definition at line 21 of file CircleCutUserObject.h.

Constructor & Destructor Documentation

◆ CircleCutUserObject()

CircleCutUserObject::CircleCutUserObject ( const InputParameters &  parameters)

Definition at line 36 of file CircleCutUserObject.C.

37  : GeometricCut3DUserObject(parameters), _cut_data(getParam<std::vector<Real>>("cut_data"))
38 {
39  // Set up constant parameters
40  const int cut_data_len = 9;
41 
42  // Throw error if length of cut_data is incorrect
43  if (_cut_data.size() != cut_data_len)
44  mooseError("Length of CircleCutUserObject cut_data must be 9");
45 
46  // Assign cut_data to vars used to construct cuts
47  _center = Point(_cut_data[0], _cut_data[1], _cut_data[2]);
48  _vertices.push_back(Point(_cut_data[3], _cut_data[4], _cut_data[5]));
49  _vertices.push_back(Point(_cut_data[6], _cut_data[7], _cut_data[8]));
50 
51  std::pair<Point, Point> rays = std::make_pair(_vertices[0] - _center, _vertices[1] - _center);
52 
53  _normal = rays.first.cross(rays.second);
55 
56  std::pair<Real, Real> ray_radii =
57  std::make_pair(std::sqrt(rays.first.norm_sq()), std::sqrt(rays.second.norm_sq()));
58 
59  if (std::abs(ray_radii.first - ray_radii.second) > 1e-10)
60  mooseError("CircleCutUserObject only works for a circular cut");
61 
62  _radius = 0.5 * (ray_radii.first + ray_radii.second);
63  _angle = std::acos((rays.first * rays.second) / (ray_radii.first * ray_radii.second));
64 }
GeometricCut3DUserObject(const InputParameters &parameters)
std::vector< Real > _cut_data
std::vector< Point > _vertices
void normalizePoint(Point &p)
Definition: XFEMFuncs.C:621

Member Function Documentation

◆ cutElementByGeometry() [1/2]

bool GeometricCut3DUserObject::cutElementByGeometry ( const Elem *  elem,
std::vector< Xfem::CutEdge > &  cut_edges,
std::vector< Xfem::CutNode > &  cut_nodes,
Real  time 
) const
overridevirtualinherited

Check to see whether a specified 2D element should be cut based on geometric conditions.

Parameters
elemPointer to the libMesh element to be considered for cutting
cut_edgesData structure filled with information about edges to be cut
cut_nodesData structure filled with information about nodes to be cut
timeCurrent simulation time
Returns
bool true if element is to be cut

Implements GeometricCutUserObject.

Definition at line 39 of file GeometricCut3DUserObject.C.

43 {
44  mooseError("Invalid method: must use vector of element faces for 3D mesh cutting");
45  return false;
46 }

◆ cutElementByGeometry() [2/2]

bool GeometricCut3DUserObject::cutElementByGeometry ( const Elem *  elem,
std::vector< Xfem::CutFace > &  cut_faces,
Real  time 
) const
overridevirtualinherited

Check to see whether a specified 3D element should be cut based on geometric conditions.

Parameters
elemPointer to the libMesh element to be considered for cutting
cut_facesData structure filled with information about edges to be cut
timeCurrent simulation time
Returns
bool true if element is to be cut

Implements GeometricCutUserObject.

Definition at line 49 of file GeometricCut3DUserObject.C.

53 {
54  bool cut_elem = false;
55 
56  for (unsigned int i = 0; i < elem->n_sides(); ++i)
57  {
58  // This returns the lowest-order type of side.
59  std::unique_ptr<Elem> curr_side = elem->side(i);
60  if (curr_side->dim() != 2)
61  mooseError("In cutElementByGeometry dimension of side must be 2, but it is ",
62  curr_side->dim());
63  unsigned int n_edges = curr_side->n_sides();
64 
65  std::vector<unsigned int> cut_edges;
66  std::vector<Real> cut_pos;
67 
68  for (unsigned int j = 0; j < n_edges; j++)
69  {
70  // This returns the lowest-order type of side.
71  std::unique_ptr<Elem> curr_edge = curr_side->side(j);
72  if (curr_edge->type() != EDGE2)
73  mooseError("In cutElementByGeometry face edge must be EDGE2, but type is: ",
74  libMesh::Utility::enum_to_string(curr_edge->type()),
75  " base element type is: ",
76  libMesh::Utility::enum_to_string(elem->type()));
77  Node * node1 = curr_edge->get_node(0);
78  Node * node2 = curr_edge->get_node(1);
79 
80  Point intersection;
81  if (intersectWithEdge(*node1, *node2, intersection))
82  {
83  cut_edges.push_back(j);
84  cut_pos.push_back(getRelativePosition(*node1, *node2, intersection));
85  }
86  }
87 
88  if (cut_edges.size() == 2)
89  {
90  cut_elem = true;
91  Xfem::CutFace mycut;
92  mycut._face_id = i;
93  mycut._face_edge.push_back(cut_edges[0]);
94  mycut._face_edge.push_back(cut_edges[1]);
95  mycut._position.push_back(cut_pos[0]);
96  mycut._position.push_back(cut_pos[1]);
97  cut_faces.push_back(mycut);
98  }
99  }
100 
101  return cut_elem;
102 }
std::vector< Real > _position
Fractional distance along the cut edges where the cut is located.
virtual bool intersectWithEdge(const Point &p1, const Point &p2, Point &pint) const
Data structure defining a cut through a face.
Real getRelativePosition(const Point &p1, const Point &p2, const Point &p) const
std::vector< unsigned int > _face_edge
IDs of all cut faces.
unsigned int _face_id
ID of the cut face.

◆ cutFragmentByGeometry() [1/2]

bool GeometricCut3DUserObject::cutFragmentByGeometry ( std::vector< std::vector< Point >> &  frag_edges,
std::vector< Xfem::CutEdge > &  cut_edges,
Real  time 
) const
overridevirtualinherited

Check to see whether a fragment of a 2D element should be cut based on geometric conditions.

Parameters
frag_edgesData structure defining the current fragment to be considered
cut_edgesData structure filled with information about fragment edges to be cut
timeCurrent simulation time
Returns
bool true if fragment is to be cut

Implements GeometricCutUserObject.

Definition at line 105 of file GeometricCut3DUserObject.C.

108 {
109  mooseError("Invalid method: must use vector of element faces for 3D mesh cutting");
110  return false;
111 }

◆ cutFragmentByGeometry() [2/2]

bool GeometricCut3DUserObject::cutFragmentByGeometry ( std::vector< std::vector< Point >> &  frag_faces,
std::vector< Xfem::CutFace > &  cut_faces,
Real  time 
) const
overridevirtualinherited

Check to see whether a fragment of a 3D element should be cut based on geometric conditions.

Parameters
frag_facesData structure defining the current fragment to be considered
cut_facesData structure filled with information about fragment faces to be cut
timeCurrent simulation time
Returns
bool true if fragment is to be cut

Implements GeometricCutUserObject.

Definition at line 114 of file GeometricCut3DUserObject.C.

117 {
118  // TODO: Need this for branching in 3D
119  mooseError("cutFragmentByGeometry not yet implemented for 3D mesh cutting");
120  return false;
121 }

◆ deserialize()

void GeometricCutUserObject::deserialize ( std::vector< std::string > &  serialized_buffers)
protectedinherited

Definition at line 222 of file GeometricCutUserObject.C.

Referenced by GeometricCutUserObject::finalize().

223 {
224  mooseAssert(serialized_buffers.size() == _app.n_processors(),
225  "Unexpected size of serialized_buffers: " << serialized_buffers.size());
226 
227  // The input string stream used for deserialization
228  std::istringstream iss;
229 
230  // Loop over all datastructures for all processors to perfrom the gather operation
231  for (unsigned int rank = 0; rank < serialized_buffers.size(); ++rank)
232  {
233  // skip the current processor (its data is already in the structures)
234  if (rank == processor_id())
235  continue;
236 
237  // populate the stream with a new buffer and reset stream state
238  iss.clear();
239  iss.str(serialized_buffers[rank]);
240 
241  // Load the communicated data into temporary structures
242  std::map<unsigned int, std::vector<Xfem::GeomMarkedElemInfo2D>> other_marked_elems_2d;
243  std::map<unsigned int, std::vector<Xfem::GeomMarkedElemInfo3D>> other_marked_elems_3d;
244  dataLoad(iss, other_marked_elems_2d, this);
245  dataLoad(iss, other_marked_elems_3d, this);
246 
247  // merge the data in with the current processor's data
248  _marked_elems_2d.insert(other_marked_elems_2d.begin(), other_marked_elems_2d.end());
249  _marked_elems_3d.insert(other_marked_elems_3d.begin(), other_marked_elems_3d.end());
250  }
251 }
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
Containers with information about all 2D and 3D elements marked for cutting by this object...
void dataLoad(std::istream &stream, Xfem::CutFace &cf, void *context)
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d

◆ execute()

void GeometricCutUserObject::execute ( )
overridevirtualinherited

Reimplemented in MovingLineSegmentCutSetUserObject.

Definition at line 67 of file GeometricCutUserObject.C.

Referenced by MovingLineSegmentCutSetUserObject::execute().

68 {
69  if (_current_elem->dim() == 2)
70  {
71  std::vector<Xfem::CutEdge> elem_cut_edges;
72  std::vector<Xfem::CutNode> elem_cut_nodes;
73  std::vector<Xfem::CutEdge> frag_cut_edges;
74  std::vector<std::vector<Point>> frag_edges;
75 
76  EFAElement2D * EFAElem = _xfem->getEFAElem2D(_current_elem);
77 
78  // Don't cut again if elem has been already cut twice
79  if (!EFAElem->isFinalCut())
80  {
81  // get fragment edges
82  _xfem->getFragmentEdges(_current_elem, EFAElem, frag_edges);
83 
84  // mark cut edges for the element and its fragment
85  bool cut = cutElementByGeometry(_current_elem, elem_cut_edges, elem_cut_nodes, _t);
86  if (EFAElem->numFragments() > 0)
87  cut |= cutFragmentByGeometry(frag_edges, frag_cut_edges, _t);
88 
89  if (cut)
90  {
92  gmei2d._elem_cut_edges = elem_cut_edges;
93  gmei2d._elem_cut_nodes = elem_cut_nodes;
94  gmei2d._frag_cut_edges = frag_cut_edges;
95  gmei2d._frag_edges = frag_edges;
96  _marked_elems_2d[_current_elem->id()].push_back(gmei2d);
97  }
98  }
99  }
100  else if (_current_elem->dim() == 3)
101  {
102  std::vector<Xfem::CutFace> elem_cut_faces;
103  std::vector<Xfem::CutFace> frag_cut_faces;
104  std::vector<std::vector<Point>> frag_faces;
105 
106  EFAElement3D * EFAElem = _xfem->getEFAElem3D(_current_elem);
107 
108  // Don't cut again if elem has been already cut twice
109  if (!EFAElem->isFinalCut())
110  {
111  // get fragment edges
112  _xfem->getFragmentFaces(_current_elem, EFAElem, frag_faces);
113 
114  // mark cut faces for the element and its fragment
115  bool cut = cutElementByGeometry(_current_elem, elem_cut_faces, _t);
116  // TODO: This would be done for branching, which is not yet supported in 3D
117  // if (EFAElem->numFragments() > 0)
118  // cut |= cutFragmentByGeometry(frag_faces, frag_cut_faces, _t);
119 
120  if (cut)
121  {
123  gmei3d._elem_cut_faces = elem_cut_faces;
124  gmei3d._frag_cut_faces = frag_cut_faces;
125  gmei3d._frag_faces = frag_faces;
126  _marked_elems_3d[_current_elem->id()].push_back(gmei3d);
127  }
128  }
129  }
130 }
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
Containers with information about all 2D and 3D elements marked for cutting by this object...
Data structure describing geometrically described cut through 3D element.
virtual bool cutFragmentByGeometry(std::vector< std::vector< Point >> &frag_edges, std::vector< Xfem::CutEdge > &cut_edges, Real time) const =0
Check to see whether a fragment of a 2D element should be cut based on geometric conditions.
virtual unsigned int numFragments() const
Definition: EFAElement2D.C:202
virtual bool cutElementByGeometry(const Elem *elem, std::vector< Xfem::CutEdge > &cut_edges, std::vector< Xfem::CutNode > &cut_nodes, Real time) const =0
Check to see whether a specified 2D element should be cut based on geometric conditions.
std::vector< CutNode > _elem_cut_nodes
Container for data about all cut nodes in this element.
std::vector< CutEdge > _elem_cut_edges
Container for data about all cut edges in this element.
std::vector< std::vector< Point > > _frag_edges
Container for data about all cut edges in cut fragments in this element.
virtual bool isFinalCut() const
Definition: EFAElement2D.C:791
std::vector< CutEdge > _frag_cut_edges
Container for data about all cut fragments in this element.
Data structure describing geometrically described cut through 2D element.
std::shared_ptr< XFEM > _xfem
Pointer to the XFEM controller object.
std::vector< CutFace > _frag_cut_faces
Container for data about all faces this element&#39;s fragment.
std::vector< std::vector< Point > > _frag_faces
Container for data about all cut faces in cut fragments in this element.
std::vector< CutFace > _elem_cut_faces
Container for data about all cut faces in this element.
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d
virtual bool isFinalCut() const
Definition: EFAElement3D.C:838

◆ finalize()

void GeometricCutUserObject::finalize ( )
overridevirtualinherited

Reimplemented in MovingLineSegmentCutSetUserObject.

Definition at line 254 of file GeometricCutUserObject.C.

Referenced by MovingLineSegmentCutSetUserObject::finalize().

255 {
256  // for single processor runs we do not need to do anything here
257  if (_app.n_processors() > 1)
258  {
259  // create send buffer
260  std::string send_buffer;
261 
262  // create byte buffers for the streams received from all processors
263  std::vector<std::string> recv_buffers;
264 
265  // pack the complex datastructures into the string stream
266  serialize(send_buffer);
267 
268  // broadcast serialized data to and receive from all processors
269  _communicator.allgather(send_buffer, recv_buffers);
270 
271  // unpack the received data and merge it into the local data structures
272  deserialize(recv_buffers);
273  }
274 
275  for (const auto & it : _marked_elems_2d)
276  for (const auto & gmei : it.second)
277  _xfem->addGeomMarkedElem2D(it.first, gmei, _interface_id);
278 
279  for (const auto & it : _marked_elems_3d)
280  for (const auto & gmei : it.second)
281  _xfem->addGeomMarkedElem3D(it.first, gmei, _interface_id);
282 
283  _marked_elems_2d.clear();
284  _marked_elems_3d.clear();
285 }
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
Containers with information about all 2D and 3D elements marked for cutting by this object...
unsigned int _interface_id
Associated interface id.
void deserialize(std::vector< std::string > &serialized_buffers)
void serialize(std::string &serialized_buffer)
Methods to pack/unpack the _marked_elems_2d and _marked_elems_3d data into a structure suitable for p...
std::shared_ptr< XFEM > _xfem
Pointer to the XFEM controller object.
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d

◆ getCrackFrontPoints()

const std::vector< Point > CircleCutUserObject::getCrackFrontPoints ( unsigned int  int) const
overridevirtual

get a set of points along a crack front from a XFEM GeometricCutUserObject

Returns
A vector which contains all crack front points

Implements CrackFrontPointsProvider.

Definition at line 76 of file CircleCutUserObject.C.

77 {
78  std::vector<Point> crack_front_points(number_crack_front_points);
79  Point v1 = _vertices[0] - _center;
80  Point v2 = _normal.cross(v1);
81  v1 /= v1.norm();
82  v2 /= v2.norm();
83  // parametric circle in 3D: center + r * cos(theta) * v1 + r * sin(theta) * v2
84  for (unsigned int i = 0; i < number_crack_front_points; ++i)
85  {
86  Real theta = 2.0 * libMesh::pi / number_crack_front_points * i;
87  crack_front_points[i] =
88  _center + _radius * std::cos(theta) * v1 + _radius * std::sin(theta) * v2;
89  }
90 
91  return crack_front_points;
92 }
std::vector< Point > _vertices

◆ getInterfaceID()

unsigned int GeometricCutUserObject::getInterfaceID ( ) const
inlineinherited

Get the interface ID for this cutting object.

Returns
the interface ID

Definition at line 173 of file GeometricCutUserObject.h.

173 { return _interface_id; };
unsigned int _interface_id
Associated interface id.

◆ getRelativePosition()

Real GeometricCut3DUserObject::getRelativePosition ( const Point &  p1,
const Point &  p2,
const Point &  p 
) const
protectedinherited

Definition at line 155 of file GeometricCut3DUserObject.C.

158 {
159  // get the relative position of p from p1
160  Real full_len = (p2 - p1).norm();
161  Real len_p1_p = (p - p1).norm();
162  return len_p1_p / full_len;
163 }

◆ initialize()

void GeometricCutUserObject::initialize ( )
overridevirtualinherited

Reimplemented in MeshCut3DUserObject, and MovingLineSegmentCutSetUserObject.

Definition at line 60 of file GeometricCutUserObject.C.

61 {
62  _marked_elems_2d.clear();
63  _marked_elems_3d.clear();
64 }
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
Containers with information about all 2D and 3D elements marked for cutting by this object...
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d

◆ intersectWithEdge()

bool GeometricCut3DUserObject::intersectWithEdge ( const Point &  p1,
const Point &  p2,
Point &  pint 
) const
protectedvirtualinherited

Definition at line 124 of file GeometricCut3DUserObject.C.

125 {
126  bool has_intersection = false;
127  double plane_point[3] = {_center(0), _center(1), _center(2)};
128  double plane_normal[3] = {_normal(0), _normal(1), _normal(2)};
129  double edge_point1[3] = {p1(0), p1(1), p1(2)};
130  double edge_point2[3] = {p2(0), p2(1), p2(2)};
131  double cut_point[3] = {0.0, 0.0, 0.0};
132 
134  plane_point, plane_normal, edge_point1, edge_point2, cut_point) == 1)
135  {
136  Point temp_p(cut_point[0], cut_point[1], cut_point[2]);
137  if (isInsideCutPlane(temp_p) && isInsideEdge(p1, p2, temp_p))
138  {
139  pint = temp_p;
140  has_intersection = true;
141  }
142  }
143  return has_intersection;
144 }
virtual bool isInsideCutPlane(Point p) const =0
int plane_normal_line_exp_int_3d(double pp[3], double normal[3], double p1[3], double p2[3], double pint[3])
Definition: XFEMFuncs.C:403
bool isInsideEdge(const Point &p1, const Point &p2, const Point &p) const

◆ isInsideCutPlane()

bool CircleCutUserObject::isInsideCutPlane ( Point  p) const
overrideprivatevirtual

Implements GeometricCut3DUserObject.

Definition at line 67 of file CircleCutUserObject.C.

68 {
69  Point ray = p - _center;
70  if (std::abs(ray * _normal) < 1e-15 && std::sqrt(ray.norm_sq()) < _radius)
71  return true;
72  return false;
73 }

◆ isInsideEdge()

bool GeometricCut3DUserObject::isInsideEdge ( const Point &  p1,
const Point &  p2,
const Point &  p 
) const
protectedinherited

Definition at line 147 of file GeometricCut3DUserObject.C.

Referenced by GeometricCut3DUserObject::intersectWithEdge().

148 {
149  Real dotp1 = (p1 - p) * (p2 - p1);
150  Real dotp2 = (p2 - p) * (p2 - p1);
151  return (dotp1 * dotp2 <= 0.0);
152 }

◆ serialize()

void GeometricCutUserObject::serialize ( std::string &  serialized_buffer)
protectedinherited

Methods to pack/unpack the _marked_elems_2d and _marked_elems_3d data into a structure suitable for parallel communication.

Definition at line 209 of file GeometricCutUserObject.C.

Referenced by GeometricCutUserObject::finalize().

210 {
211  // stream for serializing the _marked_elems_2d and _marked_elems_3d data structures to a byte
212  // stream
213  std::ostringstream oss;
214  dataStore(oss, _marked_elems_2d, this);
215  dataStore(oss, _marked_elems_3d, this);
216 
217  // Populate the passed in string pointer with the string stream's buffer contents
218  serialized_buffer.assign(oss.str());
219 }
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
Containers with information about all 2D and 3D elements marked for cutting by this object...
void dataStore(std::ostream &stream, Xfem::CutFace &cf, void *context)
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d

◆ setInterfaceID()

void GeometricCutUserObject::setInterfaceID ( unsigned int  interface_id)
inlineinherited

Set the interface ID for this cutting object.

Parameters
theinterface ID

Definition at line 179 of file GeometricCutUserObject.h.

Referenced by XFEM::addGeometricCut().

179 { _interface_id = interface_id; };
unsigned int _interface_id
Associated interface id.

◆ shouldHealMesh()

bool GeometricCutUserObject::shouldHealMesh ( ) const
inlineinherited

Should the elements cut by this cutting object be healed in the current time step?

Returns
true if the cut element should be healed

Definition at line 186 of file GeometricCutUserObject.h.

186 { return _heal_always; };
bool _heal_always
Heal the mesh.

◆ threadJoin()

void GeometricCutUserObject::threadJoin ( const UserObject &  y)
overridevirtualinherited

Definition at line 133 of file GeometricCutUserObject.C.

134 {
135  const GeometricCutUserObject & gcuo = dynamic_cast<const GeometricCutUserObject &>(y);
136 
137  for (const auto & it : gcuo._marked_elems_2d)
138  {
139  mooseAssert(_marked_elems_2d.find(it.first) == _marked_elems_2d.end(),
140  "Element already inserted in map from a different thread");
141  _marked_elems_2d[it.first] = it.second;
142  }
143  for (const auto & it : gcuo._marked_elems_3d)
144  {
145  mooseAssert(_marked_elems_3d.find(it.first) == _marked_elems_3d.end(),
146  "Element already inserted in map from a different thread");
147  _marked_elems_3d[it.first] = it.second;
148  }
149 }
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo2D > > _marked_elems_2d
Containers with information about all 2D and 3D elements marked for cutting by this object...
std::map< unsigned int, std::vector< Xfem::GeomMarkedElemInfo3D > > _marked_elems_3d

Member Data Documentation

◆ _angle

Real CircleCutUserObject::_angle
private

Definition at line 35 of file CircleCutUserObject.h.

Referenced by CircleCutUserObject().

◆ _center

Point GeometricCut3DUserObject::_center
protectedinherited

◆ _cut_data

std::vector<Real> CircleCutUserObject::_cut_data
protected

Definition at line 30 of file CircleCutUserObject.h.

Referenced by CircleCutUserObject().

◆ _heal_always

bool GeometricCutUserObject::_heal_always
protectedinherited

Heal the mesh.

Definition at line 196 of file GeometricCutUserObject.h.

◆ _interface_id

unsigned int GeometricCutUserObject::_interface_id
protectedinherited

Associated interface id.

Definition at line 193 of file GeometricCutUserObject.h.

Referenced by GeometricCutUserObject::finalize(), and GeometricCutUserObject::GeometricCutUserObject().

◆ _marked_elems_2d

std::map<unsigned int, std::vector<Xfem::GeomMarkedElemInfo2D> > GeometricCutUserObject::_marked_elems_2d
protectedinherited

◆ _marked_elems_3d

std::map<unsigned int, std::vector<Xfem::GeomMarkedElemInfo3D> > GeometricCutUserObject::_marked_elems_3d
protectedinherited

◆ _normal

Point GeometricCut3DUserObject::_normal
protectedinherited

◆ _radius

Real CircleCutUserObject::_radius
private

◆ _vertices

std::vector<Point> CircleCutUserObject::_vertices
private

Definition at line 33 of file CircleCutUserObject.h.

Referenced by CircleCutUserObject(), and getCrackFrontPoints().

◆ _xfem

std::shared_ptr<XFEM> GeometricCutUserObject::_xfem
protectedinherited

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