www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
XFEMCutElem Class Referenceabstract

#include <XFEMCutElem.h>

Inheritance diagram for XFEMCutElem:
[legend]

Public Member Functions

 XFEMCutElem (Elem *elem, unsigned int n_qpoints, unsigned int n_sides)
 Constructor initializes XFEMCutElem object. More...
 
virtual ~XFEMCutElem ()
 
void setQuadraturePointsAndWeights (const std::vector< Point > &qp_points, const std::vector< Real > &qp_weights)
 
virtual void computePhysicalVolumeFraction ()=0
 Computes the volume fraction of the element fragment. More...
 
Real getPhysicalVolumeFraction () const
 Returns the volume fraction of the element fragment. More...
 
virtual void computePhysicalFaceAreaFraction (unsigned int side)=0
 Computes the surface area fraction of the element side. More...
 
Real getPhysicalFaceAreaFraction (unsigned int side) const
 Returns the surface area fraction of the element side. More...
 
virtual void computeMomentFittingWeights ()=0
 
Real getMomentFittingWeight (unsigned int i_qp) const
 
virtual Point getCutPlaneOrigin (unsigned int plane_id, MeshBase *displaced_mesh=NULL) const =0
 
virtual Point getCutPlaneNormal (unsigned int plane_id, MeshBase *displaced_mesh=NULL) const =0
 
virtual void getCrackTipOriginAndDirection (unsigned tip_id, Point &origin, Point &direction) const =0
 
virtual void getFragmentFaces (std::vector< std::vector< Point >> &frag_faces, MeshBase *displaced_mesh=NULL) const =0
 
virtual const EFAElementgetEFAElement () const =0
 
virtual unsigned int numCutPlanes () const =0
 
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
 
virtual void getIntersectionInfo (unsigned int plane_id, Point &normal, std::vector< Point > &intersectionPoints, MeshBase *displaced_mesh=NULL) const =0
 

Protected Member Functions

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

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

Detailed Description

Definition at line 30 of file XFEMCutElem.h.

Constructor & Destructor Documentation

◆ XFEMCutElem()

XFEMCutElem::XFEMCutElem ( Elem *  elem,
unsigned int  n_qpoints,
unsigned int  n_sides 
)

Constructor initializes XFEMCutElem object.

Parameters
elemThe element on which XFEMCutElem is built
n_qpointsThe number of quadrature points
n_sidesThe number of sides which the element has

Definition at line 22 of file XFEMCutElem.C.

23  : _n_nodes(elem->n_nodes()),
24  _n_qpoints(n_qpoints),
25  _n_sides(n_sides),
26  _nodes(_n_nodes, NULL),
27  _elem_side_area(n_sides),
28  _physical_areafrac(n_sides),
29  _have_weights(false),
30  _have_face_weights(n_sides),
31  _new_face_weights(n_sides)
32 {
33  for (unsigned int i = 0; i < _n_nodes; ++i)
34  _nodes[i] = elem->get_node(i);
35 
36  for (unsigned int i = 0; i < _n_sides; ++i)
37  {
38  _have_face_weights[i] = false;
39  _physical_areafrac[i] = 1.0;
40  _elem_side_area[i] = elem->side(i)->volume();
41  }
42 
43  _elem_volume = elem->volume();
44 }
std::vector< Real > _elem_side_area
Definition: XFEMCutElem.h:50
std::vector< std::vector< Real > > _new_face_weights
face quadrature weights from surface area fraction
Definition: XFEMCutElem.h:58
std::vector< bool > _have_face_weights
Definition: XFEMCutElem.h:54
unsigned int _n_sides
Definition: XFEMCutElem.h:45
unsigned int _n_qpoints
Definition: XFEMCutElem.h:44
bool _have_weights
Definition: XFEMCutElem.h:53
unsigned int _n_nodes
Definition: XFEMCutElem.h:43
std::vector< Node * > _nodes
Definition: XFEMCutElem.h:46
Real _elem_volume
Definition: XFEMCutElem.h:49
std::vector< Real > _physical_areafrac
Definition: XFEMCutElem.h:52

◆ ~XFEMCutElem()

XFEMCutElem::~XFEMCutElem ( )
virtual

Definition at line 46 of file XFEMCutElem.C.

46 {}

Member Function Documentation

◆ computeMomentFittingWeights()

virtual void XFEMCutElem::computeMomentFittingWeights ( )
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by computeXFEMWeights().

◆ computePhysicalFaceAreaFraction()

virtual void XFEMCutElem::computePhysicalFaceAreaFraction ( unsigned int  side)
pure virtual

Computes the surface area fraction of the element side.

Parameters
sideThe side of the element

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by computeXFEMFaceWeights().

◆ computePhysicalVolumeFraction()

virtual void XFEMCutElem::computePhysicalVolumeFraction ( )
pure virtual

Computes the volume fraction of the element fragment.

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

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

◆ computeXFEMFaceWeights()

void XFEMCutElem::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.

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 90 of file XFEMCutElem.C.

Referenced by getFaceWeightMultipliers().

94 {
95  _new_face_weights[side].clear();
96  _new_face_weights[side].resize(qrule->n_points(), 1.0);
97 
99  Real surffrac = getPhysicalFaceAreaFraction(side);
100  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
101  _new_face_weights[side][qp] = surffrac;
102 
103  _have_face_weights[side] = true;
104 }
std::vector< std::vector< Real > > _new_face_weights
face quadrature weights from surface area fraction
Definition: XFEMCutElem.h:58
std::vector< bool > _have_face_weights
Definition: XFEMCutElem.h:54
virtual void computePhysicalFaceAreaFraction(unsigned int side)=0
Computes the surface area fraction of the element side.
Real getPhysicalFaceAreaFraction(unsigned int side) const
Returns the surface area fraction of the element side.
Definition: XFEMCutElem.C:55

◆ computeXFEMWeights()

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

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 107 of file XFEMCutElem.C.

Referenced by getWeightMultipliers().

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

◆ getCrackTipOriginAndDirection()

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

◆ getCutPlaneNormal()

virtual Point XFEMCutElem::getCutPlaneNormal ( unsigned int  plane_id,
MeshBase *  displaced_mesh = NULL 
) const
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by XFEM::getCutPlane(), and isPointPhysical().

◆ getCutPlaneOrigin()

virtual Point XFEMCutElem::getCutPlaneOrigin ( unsigned int  plane_id,
MeshBase *  displaced_mesh = NULL 
) const
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by XFEM::getCutPlane(), and isPointPhysical().

◆ getEFAElement()

virtual const EFAElement* XFEMCutElem::getEFAElement ( ) const
pure virtual

◆ getFaceWeightMultipliers()

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

Definition at line 75 of file XFEMCutElem.C.

Referenced by XFEM::getXFEMFaceWeights().

80 {
81  if (!_have_face_weights[side])
82  computeXFEMFaceWeights(qrule, xfem_qrule, q_points, side);
83 
84  face_weights.resize(_new_face_weights[side].size());
85  for (unsigned int qp = 0; qp < _new_face_weights[side].size(); ++qp)
86  face_weights[qp] = _new_face_weights[side][qp];
87 }
std::vector< std::vector< Real > > _new_face_weights
face quadrature weights from surface area fraction
Definition: XFEMCutElem.h:58
std::vector< bool > _have_face_weights
Definition: XFEMCutElem.h:54
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:90

◆ getFragmentFaces()

virtual void XFEMCutElem::getFragmentFaces ( std::vector< std::vector< Point >> &  frag_faces,
MeshBase *  displaced_mesh = NULL 
) const
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by XFEM::getFragmentFaces().

◆ getIntersectionInfo()

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

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by XFEM::getXFEMIntersectionInfo().

◆ getMomentFittingWeight()

Real XFEMCutElem::getMomentFittingWeight ( unsigned int  i_qp) const

◆ getNodeCoordinates()

virtual Point XFEMCutElem::getNodeCoordinates ( EFANode node,
MeshBase *  displaced_mesh = NULL 
) const
protectedpure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

◆ getPhysicalFaceAreaFraction()

Real XFEMCutElem::getPhysicalFaceAreaFraction ( unsigned int  side) const

Returns the surface area fraction of the element side.

Parameters
sideThe side of the element

Definition at line 55 of file XFEMCutElem.C.

Referenced by computeXFEMFaceWeights().

56 {
57  return _physical_areafrac[side];
58 }
std::vector< Real > _physical_areafrac
Definition: XFEMCutElem.h:52

◆ getPhysicalVolumeFraction()

Real XFEMCutElem::getPhysicalVolumeFraction ( ) const

Returns the volume fraction of the element fragment.

Definition at line 49 of file XFEMCutElem.C.

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

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

◆ getWeightMultipliers()

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

Definition at line 61 of file XFEMCutElem.C.

Referenced by XFEM::getXFEMWeights().

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

◆ isPointPhysical()

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

Definition at line 191 of file XFEMCutElem.C.

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

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

◆ numCutPlanes()

virtual unsigned int XFEMCutElem::numCutPlanes ( ) const
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by isPointPhysical().

◆ setQuadraturePointsAndWeights()

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

Member Data Documentation

◆ _elem_side_area

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

◆ _elem_volume

Real XFEMCutElem::_elem_volume
protected

◆ _have_face_weights

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

Definition at line 54 of file XFEMCutElem.h.

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

◆ _have_weights

bool XFEMCutElem::_have_weights
protected

Definition at line 53 of file XFEMCutElem.h.

Referenced by computeXFEMWeights(), and getWeightMultipliers().

◆ _n_nodes

unsigned int XFEMCutElem::_n_nodes
protected

Definition at line 43 of file XFEMCutElem.h.

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

◆ _n_qpoints

unsigned int XFEMCutElem::_n_qpoints
protected

◆ _n_sides

unsigned int XFEMCutElem::_n_sides
protected

Definition at line 45 of file XFEMCutElem.h.

Referenced by XFEMCutElem().

◆ _new_face_weights

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

face quadrature weights from surface area fraction

Definition at line 58 of file XFEMCutElem.h.

Referenced by computeXFEMFaceWeights(), and getFaceWeightMultipliers().

◆ _new_weights

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

quadrature weights from volume fraction and moment fitting

Definition at line 56 of file XFEMCutElem.h.

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

◆ _nodes

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

◆ _physical_areafrac

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

◆ _physical_volfrac

Real XFEMCutElem::_physical_volfrac
protected

◆ _qp_points

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

Definition at line 47 of file XFEMCutElem.h.

Referenced by computeXFEMWeights(), and XFEMCutElem2D::solveMomentFitting().

◆ _qp_weights

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

Definition at line 48 of file XFEMCutElem.h.

Referenced by computeXFEMWeights(), and XFEMCutElem2D::solveMomentFitting().


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