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

#include <XFEMCutElem3D.h>

Inheritance diagram for XFEMCutElem3D:
[legend]

Public Member Functions

 XFEMCutElem3D (Elem *elem, const EFAElement3D *const CEMelem, unsigned int n_qpoints, unsigned int n_sides)
 Constructor initializes XFEMCutElem3D object. More...
 
 ~XFEMCutElem3D ()
 
virtual void computePhysicalVolumeFraction ()
 Computes the volume fraction of the element fragment. More...
 
virtual void computePhysicalFaceAreaFraction (unsigned int side)
 Computes the surface area fraction of the element side. More...
 
virtual void computeMomentFittingWeights ()
 
virtual Point getCutPlaneOrigin (unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const
 
virtual Point getCutPlaneNormal (unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const
 
virtual void getCrackTipOriginAndDirection (unsigned tip_id, Point &origin, Point &direction) const
 
virtual void getFragmentFaces (std::vector< std::vector< Point >> &frag_faces, MeshBase *displaced_mesh=nullptr) const
 
virtual const EFAElementgetEFAElement () const
 
virtual unsigned int numCutPlanes () const
 
virtual void getIntersectionInfo (unsigned int plane_id, Point &normal, std::vector< Point > &intersectionPoints, MeshBase *displaced_mesh=nullptr) const
 
void setQuadraturePointsAndWeights (const std::vector< Point > &qp_points, const std::vector< Real > &qp_weights)
 
Real getPhysicalVolumeFraction () const
 Returns the volume fraction of the element fragment. More...
 
Real getPhysicalFaceAreaFraction (unsigned int side) const
 Returns the surface area fraction of the element side. More...
 
Real getMomentFittingWeight (unsigned int i_qp) const
 
void getWeightMultipliers (MooseArray< Real > &weights, QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
 
void getFaceWeightMultipliers (MooseArray< Real > &face_weights, QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points, unsigned int side)
 
void computeXFEMWeights (QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
 Computes integration weights for the cut element. More...
 
void computeXFEMFaceWeights (QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points, unsigned int side)
 Computes face integration weights for the cut element side. More...
 
bool isPointPhysical (const Point &p) const
 

Protected Attributes

unsigned int _n_nodes
 
unsigned int _n_qpoints
 
unsigned int _n_sides
 
std::vector< Node * > _nodes
 
std::vector< Point_qp_points
 
std::vector< Real_qp_weights
 
Real _elem_volume
 
std::vector< Real_elem_side_area
 
Real _physical_volfrac
 
std::vector< Real_physical_areafrac
 
bool _have_weights
 
std::vector< bool > _have_face_weights
 
std::vector< Real_new_weights
 quadrature weights from volume fraction and moment fitting More...
 
std::vector< std::vector< Real > > _new_face_weights
 face quadrature weights from surface area fraction More...
 

Private Member Functions

virtual Point getNodeCoordinates (EFANode *node, MeshBase *displaced_mesh=nullptr) const
 

Private Attributes

EFAElement3D _efa_elem3d
 

Detailed Description

Definition at line 24 of file XFEMCutElem3D.h.

Constructor & Destructor Documentation

◆ XFEMCutElem3D()

XFEMCutElem3D::XFEMCutElem3D ( Elem *  elem,
const EFAElement3D *const  CEMelem,
unsigned int  n_qpoints,
unsigned int  n_sides 
)

Constructor initializes XFEMCutElem3D object.

Parameters
elemThe element on which XFEMCutElem3D is built
CEMelemThe EFAFragment3D object that belongs to XFEMCutElem3D
n_qpointsThe number of quadrature points
n_sidesThe number of sides which the element has

Definition at line 23 of file XFEMCutElem3D.C.

27  : XFEMCutElem(elem, n_qpoints, n_sides), _efa_elem3d(CEMelem, true)
28 {
31 }
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
virtual void computeMomentFittingWeights()
virtual void computePhysicalVolumeFraction()
Computes the volume fraction of the element fragment.
Definition: XFEMCutElem3D.C:63
XFEMCutElem(Elem *elem, unsigned int n_qpoints, unsigned int n_sides)
Constructor initializes XFEMCutElem object.
Definition: XFEMCutElem.C:22

◆ ~XFEMCutElem3D()

XFEMCutElem3D::~XFEMCutElem3D ( )

Definition at line 33 of file XFEMCutElem3D.C.

33 {}

Member Function Documentation

◆ computeMomentFittingWeights()

void XFEMCutElem3D::computeMomentFittingWeights ( )
virtual

Implements XFEMCutElem.

Definition at line 152 of file XFEMCutElem3D.C.

Referenced by XFEMCutElem3D().

153 {
154  // TODO: 3D moment fitting method nod coded yet - use volume fraction for now
156 }
unsigned int _n_qpoints
Definition: XFEMCutElem.h:43
Real _physical_volfrac
Definition: XFEMCutElem.h:50
std::vector< Real > _new_weights
quadrature weights from volume fraction and moment fitting
Definition: XFEMCutElem.h:55

◆ computePhysicalFaceAreaFraction()

void XFEMCutElem3D::computePhysicalFaceAreaFraction ( unsigned int  side)
virtual

Computes the surface area fraction of the element side.

Parameters
sideThe side of the element

find a fragment surface which is covered by element side

Implements XFEMCutElem.

Definition at line 108 of file XFEMCutElem3D.C.

109 {
110  Real frag_surf = 0.0;
111 
112  std::vector<std::vector<unsigned int>> frag_face_indices;
113  std::vector<EFANode *> frag_nodes;
114  _efa_elem3d.getFragment(0)->getNodeInfo(frag_face_indices, frag_nodes);
115  int face_num = frag_face_indices.size();
116 
117  EFAFace * efa_face = _efa_elem3d.getFace(side);
118  bool contains_all = true;
119 
121  for (int i = 0; i < face_num; ++i)
122  {
123  contains_all = true;
124  for (unsigned int j = 0; j < frag_face_indices[i].size(); ++j)
125  {
126  EFANode * efa_node = frag_nodes[frag_face_indices[i][j]];
127  if (!efa_face->containsNode(efa_node))
128  {
129  contains_all = false;
130  break;
131  }
132  }
133 
134  if (contains_all)
135  {
136  for (unsigned int j = 0; j < frag_face_indices[i].size(); ++j)
137  {
138  unsigned int m = ((j + 1) == frag_face_indices[i].size()) ? 0 : j + 1;
139  Point edge_p1 = getNodeCoordinates(frag_nodes[frag_face_indices[i][j]]);
140  Point edge_p2 = getNodeCoordinates(frag_nodes[frag_face_indices[i][m]]);
141 
142  frag_surf += 0.5 * (edge_p1(0) - edge_p2(0)) * (edge_p1(1) + edge_p2(1));
143  }
144  _physical_areafrac[side] = std::abs(frag_surf) / _elem_side_area[side];
145  return;
146  }
147  }
148  _physical_areafrac[side] = 1.0;
149 }
EFAFragment3D * getFragment(unsigned int frag_id) const
EFAFace * getFace(unsigned int face_id) const
std::vector< Real > _elem_side_area
Definition: XFEMCutElem.h:49
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
bool containsNode(const EFANode *node) const
Definition: EFAFace.C:362
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
void getNodeInfo(std::vector< std::vector< unsigned int >> &face_node_indices, std::vector< EFANode *> &nodes) const
virtual Point getNodeCoordinates(EFANode *node, MeshBase *displaced_mesh=nullptr) const
Definition: XFEMCutElem3D.C:36
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
std::vector< Real > _physical_areafrac
Definition: XFEMCutElem.h:51

◆ computePhysicalVolumeFraction()

void XFEMCutElem3D::computePhysicalVolumeFraction ( )
virtual

Computes the volume fraction of the element fragment.

Implements XFEMCutElem.

Definition at line 63 of file XFEMCutElem3D.C.

Referenced by XFEMCutElem3D().

64 {
65  Real frag_vol = 0.0;
66 
67  // collect fragment info needed by polyhedron_volume_3d()
68  std::vector<std::vector<unsigned int>> frag_face_indices;
69  std::vector<EFANode *> frag_nodes;
70  _efa_elem3d.getFragment(0)->getNodeInfo(frag_face_indices, frag_nodes);
71  int face_num = frag_face_indices.size();
72  int node_num = frag_nodes.size();
73 
74  int order_max = 0;
75  int * order = new int[face_num];
76  for (int i = 0; i < face_num; ++i)
77  {
78  if (frag_face_indices[i].size() > (unsigned int)order_max)
79  order_max = frag_face_indices[i].size();
80  order[i] = frag_face_indices[i].size();
81  }
82 
83  double * coord = new double[3 * node_num];
84  for (unsigned int i = 0; i < frag_nodes.size(); ++i)
85  {
86  Point p = getNodeCoordinates(frag_nodes[i]);
87  coord[3 * i + 0] = p(0);
88  coord[3 * i + 1] = p(1);
89  coord[3 * i + 2] = p(2);
90  }
91 
92  int * node = new int[face_num * order_max];
93  Xfem::i4vec_zero(face_num * order_max, node);
94  for (int i = 0; i < face_num; ++i)
95  for (unsigned int j = 0; j < frag_face_indices[i].size(); ++j)
96  node[order_max * i + j] = frag_face_indices[i][j];
97 
98  // compute fragment volume and volume fraction
99  frag_vol = Xfem::polyhedron_volume_3d(coord, order_max, face_num, node, node_num, order);
100  _physical_volfrac = frag_vol / _elem_volume;
101 
102  delete[] order;
103  delete[] coord;
104  delete[] node;
105 }
EFAFragment3D * getFragment(unsigned int frag_id) const
double polyhedron_volume_3d(double coord[], int order_max, int face_num, int node[], int node_num, int order[])
Definition: XFEMFuncs.C:489
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
Real _physical_volfrac
Definition: XFEMCutElem.h:50
void getNodeInfo(std::vector< std::vector< unsigned int >> &face_node_indices, std::vector< EFANode *> &nodes) const
void i4vec_zero(int n, int a[])
Definition: XFEMFuncs.C:586
virtual Point getNodeCoordinates(EFANode *node, MeshBase *displaced_mesh=nullptr) const
Definition: XFEMCutElem3D.C:36
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real _elem_volume
Definition: XFEMCutElem.h:48
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")

◆ computeXFEMFaceWeights()

void XFEMCutElem::computeXFEMFaceWeights ( QBase *  qrule,
Xfem::XFEM_QRULE  xfem_qrule,
const MooseArray< Point > &  q_points,
unsigned int  side 
)
inherited

Computes face integration weights for the cut element side.

Parameters
qruleThe standard MOOSE face quadrature rule
xfem_qruleThe integration scheme for the cut element (We use surface area fraction only)
q_pointsThe quadrature points for the element side
sideThe side of the element

Definition at line 92 of file XFEMCutElem.C.

Referenced by XFEMCutElem::getFaceWeightMultipliers().

96 {
97  _new_face_weights[side].clear();
98  _new_face_weights[side].resize(qrule->n_points(), 1.0);
99 
101  Real surffrac = getPhysicalFaceAreaFraction(side);
102  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
103  _new_face_weights[side][qp] = surffrac;
104 
105  _have_face_weights[side] = true;
106 }
std::vector< std::vector< Real > > _new_face_weights
face quadrature weights from surface area fraction
Definition: XFEMCutElem.h:57
std::vector< bool > _have_face_weights
Definition: XFEMCutElem.h:53
virtual void computePhysicalFaceAreaFraction(unsigned int side)=0
Computes the surface area fraction of the element side.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real getPhysicalFaceAreaFraction(unsigned int side) const
Returns the surface area fraction of the element side.
Definition: XFEMCutElem.C:57

◆ computeXFEMWeights()

void XFEMCutElem::computeXFEMWeights ( QBase *  qrule,
Xfem::XFEM_QRULE  xfem_qrule,
const MooseArray< Point > &  q_points 
)
inherited

Computes integration weights for the cut element.

Parameters
qruleThe standard MOOSE quadrature rule
xfem_qruleThe integration scheme for the cut element
q_pointsThe quadrature points for the element

Definition at line 109 of file XFEMCutElem.C.

Referenced by XFEMCutElem::getWeightMultipliers().

112 {
113  _new_weights.clear();
114 
115  _new_weights.resize(qrule->n_points(), 1.0);
116 
117  switch (xfem_qrule)
118  {
119  case Xfem::VOLFRAC:
120  {
122  Real volfrac = getPhysicalVolumeFraction();
123  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
124  _new_weights[qp] = volfrac;
125  break;
126  }
128  {
129  // These are the coordinates in parametric coordinates
130  _qp_points = qrule->get_points();
131  _qp_weights = qrule->get_weights();
132 
134 
135  // Blend weights from moment fitting and volume fraction to avoid negative weights
136  Real alpha = 1.0;
137  if (*std::min_element(_new_weights.begin(), _new_weights.end()) < 0.0)
138  {
139  // One or more of the weights computed by moment fitting is negative.
140  // Blend moment and volume fraction weights to keep them nonnegative.
141  // Find the largest value of alpha that will keep all of the weights nonnegative.
142  for (unsigned int i = 0; i < _n_qpoints; ++i)
143  {
144  const Real denominator = _physical_volfrac - _new_weights[i];
145  if (denominator > 0.0) // Negative values would give a negative value for alpha, which
146  // must be between 0 and 1.
147  {
148  const Real alpha_i = _physical_volfrac / denominator;
149  if (alpha_i < alpha)
150  alpha = alpha_i;
151  }
152  }
153  for (unsigned int i = 0; i < _n_qpoints; ++i)
154  {
155  _new_weights[i] = alpha * _new_weights[i] + (1.0 - alpha) * _physical_volfrac;
156  if (_new_weights[i] < 0.0) // We can end up with small (roundoff) negative weights
157  _new_weights[i] = 0.0;
158  }
159  }
160  break;
161  }
162  case Xfem::DIRECT: // remove q-points outside the partial element's physical domain
163  {
164  bool nonzero = false;
165  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
166  {
167  // q_points contains quadrature point locations in physical coordinates
168  if (isPointPhysical(q_points[qp]))
169  {
170  nonzero = true;
171  _new_weights[qp] = 1.0;
172  }
173  else
174  _new_weights[qp] = 0.0;
175  }
176  if (!nonzero)
177  {
178  // Set the weights to a small value (1e-3) to avoid having DOFs
179  // with zeros on the diagonal, which occurs for nodes that are
180  // connected to no physical material.
181  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
182  _new_weights[qp] = 1e-3;
183  }
184  break;
185  }
186  default:
187  mooseError("Undefined option for XFEM_QRULE");
188  }
189  _have_weights = true;
190 }
unsigned int _n_qpoints
Definition: XFEMCutElem.h:43
void mooseError(Args &&... args)
std::vector< Point > _qp_points
Definition: XFEMCutElem.h:46
Real getPhysicalVolumeFraction() const
Returns the volume fraction of the element fragment.
Definition: XFEMCutElem.C:51
Real _physical_volfrac
Definition: XFEMCutElem.h:50
std::vector< Real > _qp_weights
Definition: XFEMCutElem.h:47
std::vector< Real > _new_weights
quadrature weights from volume fraction and moment fitting
Definition: XFEMCutElem.h:55
bool _have_weights
Definition: XFEMCutElem.h:52
virtual void computePhysicalVolumeFraction()=0
Computes the volume fraction of the element fragment.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::string alpha
Definition: NS.h:126
virtual void computeMomentFittingWeights()=0
bool isPointPhysical(const Point &p) const
Definition: XFEMCutElem.C:193

◆ getCrackTipOriginAndDirection()

void XFEMCutElem3D::getCrackTipOriginAndDirection ( unsigned  tip_id,
Point &  origin,
Point &  direction 
) const
virtual

Implements XFEMCutElem.

Definition at line 232 of file XFEMCutElem3D.C.

235 {
236  // TODO: not implemented for 3D
237  mooseError("getCrackTipOriginAndDirection not yet implemented for XFEMCutElem3D");
238 }
void mooseError(Args &&... args)

◆ getCutPlaneNormal()

Point XFEMCutElem3D::getCutPlaneNormal ( unsigned int  plane_id,
MeshBase *  displaced_mesh = nullptr 
) const
virtual

Implements XFEMCutElem.

Definition at line 190 of file XFEMCutElem3D.C.

Referenced by getIntersectionInfo().

191 {
192  Point normal(0.0, 0.0, 0.0);
193  std::vector<std::vector<EFANode *>> cut_plane_nodes;
194  for (unsigned int i = 0; i < _efa_elem3d.getFragment(0)->numFaces(); ++i)
195  {
197  {
198  EFAFace * face = _efa_elem3d.getFragment(0)->getFace(i);
199  std::vector<EFANode *> node_line;
200  for (unsigned int j = 0; j < face->numNodes(); ++j)
201  node_line.push_back(face->getNode(j));
202  cut_plane_nodes.push_back(node_line);
203  }
204  }
205  if (cut_plane_nodes.size() == 0)
206  mooseError("no cut plane found in this element");
207  if (plane_id < cut_plane_nodes.size()) // valid plane_id
208  {
209  std::vector<Point> cut_plane_points;
210  for (unsigned int i = 0; i < cut_plane_nodes[plane_id].size(); ++i)
211  cut_plane_points.push_back(getNodeCoordinates(cut_plane_nodes[plane_id][i], displaced_mesh));
212 
213  Point center(0.0, 0.0, 0.0);
214  for (unsigned int i = 0; i < cut_plane_points.size(); ++i)
215  center += cut_plane_points[i];
216  center /= cut_plane_points.size();
217 
218  for (unsigned int i = 0; i < cut_plane_points.size(); ++i)
219  {
220  unsigned int iplus1 = i < cut_plane_points.size() - 1 ? i + 1 : 0;
221  Point ray1 = cut_plane_points[i] - center;
222  Point ray2 = cut_plane_points[iplus1] - center;
223  normal += ray1.cross(ray2);
224  }
225  normal /= cut_plane_points.size();
226  }
227  Xfem::normalizePoint(normal);
228  return normal;
229 }
EFAFragment3D * getFragment(unsigned int frag_id) const
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
void mooseError(Args &&... args)
EFANode * getNode(unsigned int node_id) const
Definition: EFAFace.C:99
bool isFaceInterior(unsigned int face_id) const
EFAFace * getFace(unsigned int face_id) const
unsigned int numNodes() const
Definition: EFAFace.C:87
virtual Point getNodeCoordinates(EFANode *node, MeshBase *displaced_mesh=nullptr) const
Definition: XFEMCutElem3D.C:36
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
void normalizePoint(Point &p)
Definition: XFEMFuncs.C:621
static const std::string center
Definition: NS.h:28
unsigned int numFaces() const

◆ getCutPlaneOrigin()

Point XFEMCutElem3D::getCutPlaneOrigin ( unsigned int  plane_id,
MeshBase *  displaced_mesh = nullptr 
) const
virtual

Implements XFEMCutElem.

Definition at line 159 of file XFEMCutElem3D.C.

160 {
161  Point orig(0.0, 0.0, 0.0);
162  std::vector<std::vector<EFANode *>> cut_plane_nodes;
163  for (unsigned int i = 0; i < _efa_elem3d.getFragment(0)->numFaces(); ++i)
164  {
166  {
167  EFAFace * face = _efa_elem3d.getFragment(0)->getFace(i);
168  std::vector<EFANode *> node_line;
169  for (unsigned int j = 0; j < face->numNodes(); ++j)
170  node_line.push_back(face->getNode(j));
171  cut_plane_nodes.push_back(node_line);
172  }
173  }
174  if (cut_plane_nodes.size() == 0)
175  mooseError("no cut plane found in this element");
176  if (plane_id < cut_plane_nodes.size()) // valid plane_id
177  {
178  std::vector<Point> cut_plane_points;
179  for (unsigned int i = 0; i < cut_plane_nodes[plane_id].size(); ++i)
180  cut_plane_points.push_back(getNodeCoordinates(cut_plane_nodes[plane_id][i], displaced_mesh));
181 
182  for (unsigned int i = 0; i < cut_plane_points.size(); ++i)
183  orig += cut_plane_points[i];
184  orig /= cut_plane_points.size();
185  }
186  return orig;
187 }
EFAFragment3D * getFragment(unsigned int frag_id) const
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
void mooseError(Args &&... args)
EFANode * getNode(unsigned int node_id) const
Definition: EFAFace.C:99
bool isFaceInterior(unsigned int face_id) const
EFAFace * getFace(unsigned int face_id) const
unsigned int numNodes() const
Definition: EFAFace.C:87
virtual Point getNodeCoordinates(EFANode *node, MeshBase *displaced_mesh=nullptr) const
Definition: XFEMCutElem3D.C:36
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
unsigned int numFaces() const

◆ getEFAElement()

const EFAElement * XFEMCutElem3D::getEFAElement ( ) const
virtual

Implements XFEMCutElem.

Definition at line 249 of file XFEMCutElem3D.C.

250 {
251  return &_efa_elem3d;
252 }
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41

◆ getFaceWeightMultipliers()

void XFEMCutElem::getFaceWeightMultipliers ( MooseArray< Real > &  face_weights,
QBase *  qrule,
Xfem::XFEM_QRULE  xfem_qrule,
const MooseArray< Point > &  q_points,
unsigned int  side 
)
inherited

Definition at line 77 of file XFEMCutElem.C.

Referenced by XFEM::getXFEMFaceWeights().

82 {
83  if (!_have_face_weights[side])
84  computeXFEMFaceWeights(qrule, xfem_qrule, q_points, side);
85 
86  face_weights.resize(_new_face_weights[side].size());
87  for (unsigned int qp = 0; qp < _new_face_weights[side].size(); ++qp)
88  face_weights[qp] = _new_face_weights[side][qp];
89 }
std::vector< std::vector< Real > > _new_face_weights
face quadrature weights from surface area fraction
Definition: XFEMCutElem.h:57
std::vector< bool > _have_face_weights
Definition: XFEMCutElem.h:53
void computeXFEMFaceWeights(QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points, unsigned int side)
Computes face integration weights for the cut element side.
Definition: XFEMCutElem.C:92
void resize(unsigned int size)

◆ getFragmentFaces()

void XFEMCutElem3D::getFragmentFaces ( std::vector< std::vector< Point >> &  frag_faces,
MeshBase *  displaced_mesh = nullptr 
) const
virtual

Implements XFEMCutElem.

Definition at line 241 of file XFEMCutElem3D.C.

243 {
244  // TODO: not implemented for 3D
245  mooseError("getFragmentFaces not yet implemented for XFEMCutElem3D");
246 }
void mooseError(Args &&... args)

◆ getIntersectionInfo()

void XFEMCutElem3D::getIntersectionInfo ( unsigned int  plane_id,
Point &  normal,
std::vector< Point > &  intersectionPoints,
MeshBase *  displaced_mesh = nullptr 
) const
virtual

Implements XFEMCutElem.

Definition at line 265 of file XFEMCutElem3D.C.

269 {
270  intersectionPoints.clear();
271  std::vector<std::vector<EFANode *>> cut_plane_nodes;
272  for (unsigned int i = 0; i < _efa_elem3d.getFragment(0)->numFaces(); ++i)
273  {
275  {
276  EFAFace * face = _efa_elem3d.getFragment(0)->getFace(i);
277  std::vector<EFANode *> node_line;
278  for (unsigned int j = 0; j < face->numNodes(); ++j)
279  node_line.push_back(face->getNode(j));
280  cut_plane_nodes.push_back(node_line);
281  }
282  }
283  if (cut_plane_nodes.size() == 0)
284  mooseError("No cut plane found in this element");
285 
286  if (plane_id < cut_plane_nodes.size()) // valid plane_id
287  {
288  intersectionPoints.resize(cut_plane_nodes[plane_id].size());
289  for (unsigned int i = 0; i < cut_plane_nodes[plane_id].size(); ++i)
290  intersectionPoints[i] = getNodeCoordinates(cut_plane_nodes[plane_id][i], displaced_mesh);
291  }
292 
293  normal = getCutPlaneNormal(plane_id, displaced_mesh);
294 }
EFAFragment3D * getFragment(unsigned int frag_id) const
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
void mooseError(Args &&... args)
virtual Point getCutPlaneNormal(unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const
EFANode * getNode(unsigned int node_id) const
Definition: EFAFace.C:99
bool isFaceInterior(unsigned int face_id) const
EFAFace * getFace(unsigned int face_id) const
unsigned int numNodes() const
Definition: EFAFace.C:87
virtual Point getNodeCoordinates(EFANode *node, MeshBase *displaced_mesh=nullptr) const
Definition: XFEMCutElem3D.C:36
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
unsigned int numFaces() const

◆ getMomentFittingWeight()

Real XFEMCutElem::getMomentFittingWeight ( unsigned int  i_qp) const
inherited

◆ getNodeCoordinates()

Point XFEMCutElem3D::getNodeCoordinates ( EFANode node,
MeshBase *  displaced_mesh = nullptr 
) const
privatevirtual

Implements XFEMCutElem.

Definition at line 36 of file XFEMCutElem3D.C.

Referenced by computePhysicalFaceAreaFraction(), computePhysicalVolumeFraction(), getCutPlaneNormal(), getCutPlaneOrigin(), and getIntersectionInfo().

37 {
38  Point node_coor(0.0, 0.0, 0.0);
39  std::vector<EFANode *> master_nodes;
40  std::vector<Point> master_points;
41  std::vector<double> master_weights;
42 
43  _efa_elem3d.getMasterInfo(CEMnode, master_nodes, master_weights);
44  for (unsigned int i = 0; i < master_nodes.size(); ++i)
45  {
46  if (master_nodes[i]->category() == EFANode::N_CATEGORY_LOCAL_INDEX)
47  {
48  Node * node = _nodes[master_nodes[i]->id()];
49  if (displaced_mesh)
50  node = displaced_mesh->node_ptr(node->id());
51  Point node_p((*node)(0), (*node)(1), (*node)(2));
52  master_points.push_back(node_p);
53  }
54  else
55  mooseError("master nodes must be local");
56  }
57  for (unsigned int i = 0; i < master_nodes.size(); ++i)
58  node_coor += master_weights[i] * master_points[i];
59  return node_coor;
60 }
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
void mooseError(Args &&... args)
virtual void getMasterInfo(EFANode *node, std::vector< EFANode *> &master_nodes, std::vector< double > &master_weights) const
Definition: EFAElement3D.C:373
std::vector< Node * > _nodes
Definition: XFEMCutElem.h:45

◆ getPhysicalFaceAreaFraction()

Real XFEMCutElem::getPhysicalFaceAreaFraction ( unsigned int  side) const
inherited

Returns the surface area fraction of the element side.

Parameters
sideThe side of the element

Definition at line 57 of file XFEMCutElem.C.

Referenced by XFEMCutElem::computeXFEMFaceWeights().

58 {
59  return _physical_areafrac[side];
60 }
std::vector< Real > _physical_areafrac
Definition: XFEMCutElem.h:51

◆ getPhysicalVolumeFraction()

Real XFEMCutElem::getPhysicalVolumeFraction ( ) const
inherited

Returns the volume fraction of the element fragment.

Definition at line 51 of file XFEMCutElem.C.

Referenced by XFEMCutElem::computeXFEMWeights(), and XFEM::getPhysicalVolumeFraction().

52 {
53  return _physical_volfrac;
54 }
Real _physical_volfrac
Definition: XFEMCutElem.h:50

◆ getWeightMultipliers()

void XFEMCutElem::getWeightMultipliers ( MooseArray< Real > &  weights,
QBase *  qrule,
Xfem::XFEM_QRULE  xfem_qrule,
const MooseArray< Point > &  q_points 
)
inherited

Definition at line 63 of file XFEMCutElem.C.

Referenced by XFEM::getXFEMWeights().

67 {
68  if (!_have_weights)
69  computeXFEMWeights(qrule, xfem_qrule, q_points);
70 
71  weights.resize(_new_weights.size());
72  for (unsigned int qp = 0; qp < _new_weights.size(); ++qp)
73  weights[qp] = _new_weights[qp];
74 }
std::vector< Real > _new_weights
quadrature weights from volume fraction and moment fitting
Definition: XFEMCutElem.h:55
bool _have_weights
Definition: XFEMCutElem.h:52
void computeXFEMWeights(QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
Computes integration weights for the cut element.
Definition: XFEMCutElem.C:109
void resize(unsigned int size)

◆ isPointPhysical()

bool XFEMCutElem::isPointPhysical ( const Point &  p) const
inherited

Definition at line 193 of file XFEMCutElem.C.

Referenced by XFEMCutElem::computeXFEMWeights(), XFEM::healMesh(), and XFEM::isPointInsidePhysicalDomain().

194 {
195  // determine whether the point is inside the physical domain of a partial element
196  bool physical_flag = true;
197  unsigned int n_cut_planes = numCutPlanes();
198  for (unsigned int plane_id = 0; plane_id < n_cut_planes; ++plane_id)
199  {
200  Point origin = getCutPlaneOrigin(plane_id);
201  Point normal = getCutPlaneNormal(plane_id);
202  Point origin2qp = p - origin;
203  if (origin2qp * normal > 0.0)
204  {
205  physical_flag = false; // Point outside pysical domain
206  break;
207  }
208  }
209  return physical_flag;
210 }
virtual unsigned int numCutPlanes() const =0
virtual Point getCutPlaneOrigin(unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const =0
virtual Point getCutPlaneNormal(unsigned int plane_id, MeshBase *displaced_mesh=nullptr) const =0

◆ numCutPlanes()

unsigned int XFEMCutElem3D::numCutPlanes ( ) const
virtual

Implements XFEMCutElem.

Definition at line 255 of file XFEMCutElem3D.C.

256 {
257  unsigned int counter = 0;
258  for (unsigned int i = 0; i < _efa_elem3d.getFragment(0)->numFaces(); ++i)
260  counter += 1;
261  return counter;
262 }
EFAFragment3D * getFragment(unsigned int frag_id) const
EFAElement3D _efa_elem3d
Definition: XFEMCutElem3D.h:41
bool isFaceInterior(unsigned int face_id) const
unsigned int numFaces() const

◆ setQuadraturePointsAndWeights()

void XFEMCutElem::setQuadraturePointsAndWeights ( const std::vector< Point > &  qp_points,
const std::vector< Real > &  qp_weights 
)
inherited

Member Data Documentation

◆ _efa_elem3d

EFAElement3D XFEMCutElem3D::_efa_elem3d
private

◆ _elem_side_area

std::vector<Real> XFEMCutElem::_elem_side_area
protectedinherited

◆ _elem_volume

Real XFEMCutElem::_elem_volume
protectedinherited

◆ _have_face_weights

std::vector<bool> XFEMCutElem::_have_face_weights
protectedinherited

◆ _have_weights

bool XFEMCutElem::_have_weights
protectedinherited

◆ _n_nodes

unsigned int XFEMCutElem::_n_nodes
protectedinherited

◆ _n_qpoints

unsigned int XFEMCutElem::_n_qpoints
protectedinherited

◆ _n_sides

unsigned int XFEMCutElem::_n_sides
protectedinherited

Definition at line 44 of file XFEMCutElem.h.

Referenced by XFEMCutElem::XFEMCutElem().

◆ _new_face_weights

std::vector<std::vector<Real> > XFEMCutElem::_new_face_weights
protectedinherited

face quadrature weights from surface area fraction

Definition at line 57 of file XFEMCutElem.h.

Referenced by XFEMCutElem::computeXFEMFaceWeights(), and XFEMCutElem::getFaceWeightMultipliers().

◆ _new_weights

std::vector<Real> XFEMCutElem::_new_weights
protectedinherited

quadrature weights from volume fraction and moment fitting

Definition at line 55 of file XFEMCutElem.h.

Referenced by XFEMCutElem2D::computeMomentFittingWeights(), computeMomentFittingWeights(), XFEMCutElem::computeXFEMWeights(), and XFEMCutElem::getWeightMultipliers().

◆ _nodes

std::vector<Node *> XFEMCutElem::_nodes
protectedinherited

◆ _physical_areafrac

std::vector<Real> XFEMCutElem::_physical_areafrac
protectedinherited

◆ _physical_volfrac

Real XFEMCutElem::_physical_volfrac
protectedinherited

◆ _qp_points

std::vector<Point> XFEMCutElem::_qp_points
protectedinherited

◆ _qp_weights

std::vector<Real> XFEMCutElem::_qp_weights
protectedinherited

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