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)
 
virtual ~XFEMCutElem ()
 
void setQuadraturePointsAndWeights (const std::vector< Point > &qp_points, const std::vector< Real > &qp_weights)
 
virtual void computePhysicalVolumeFraction ()=0
 
Real getPhysicalVolumeFraction () const
 
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 computeXFEMWeights (QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
 
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
 
std::vector< Node * > _nodes
 
std::vector< Point > _qp_points
 
std::vector< Real > _qp_weights
 
Real _elem_volume
 
Real _physical_volfrac
 
bool _have_weights
 
std::vector< Real > _new_weights
 

Detailed Description

Definition at line 28 of file XFEMCutElem.h.

Constructor & Destructor Documentation

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

Definition at line 20 of file XFEMCutElem.C.

21  : _n_nodes(elem->n_nodes()), _n_qpoints(n_qpoints), _nodes(_n_nodes, NULL), _have_weights(false)
22 {
23  for (unsigned int i = 0; i < _n_nodes; ++i)
24  _nodes[i] = elem->get_node(i);
25  _elem_volume = elem->volume();
26 }
unsigned int _n_qpoints
Definition: XFEMCutElem.h:36
bool _have_weights
Definition: XFEMCutElem.h:42
unsigned int _n_nodes
Definition: XFEMCutElem.h:35
std::vector< Node * > _nodes
Definition: XFEMCutElem.h:37
Real _elem_volume
Definition: XFEMCutElem.h:40
XFEMCutElem::~XFEMCutElem ( )
virtual

Definition at line 28 of file XFEMCutElem.C.

28 {}

Member Function Documentation

virtual void XFEMCutElem::computeMomentFittingWeights ( )
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by computeXFEMWeights().

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

Definition at line 51 of file XFEMCutElem.C.

Referenced by getWeightMultipliers().

54 {
55  _new_weights.clear();
56 
57  _new_weights.resize(qrule->n_points(), 1.0);
58 
59  switch (xfem_qrule)
60  {
61  case Xfem::VOLFRAC:
62  {
64  Real volfrac = getPhysicalVolumeFraction();
65  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
66  _new_weights[qp] = volfrac;
67  break;
68  }
70  {
71  // These are the coordinates in parametric coordinates
72  _qp_points = qrule->get_points();
73  _qp_weights = qrule->get_weights();
74 
76 
77  // Blend weights from moment fitting and volume fraction to avoid negative weights
78  Real alpha = 1.0;
79  if (*std::min_element(_new_weights.begin(), _new_weights.end()) < 0.0)
80  {
81  // One or more of the weights computed by moment fitting is negative.
82  // Blend moment and volume fraction weights to keep them nonnegative.
83  // Find the largest value of alpha that will keep all of the weights nonnegative.
84  for (unsigned int i = 0; i < _n_qpoints; ++i)
85  {
86  const Real denominator = _physical_volfrac - _new_weights[i];
87  if (denominator > 0.0) // Negative values would give a negative value for alpha, which
88  // must be between 0 and 1.
89  {
90  const Real alpha_i = _physical_volfrac / denominator;
91  if (alpha_i < alpha)
92  alpha = alpha_i;
93  }
94  }
95  for (unsigned int i = 0; i < _n_qpoints; ++i)
96  {
97  _new_weights[i] = alpha * _new_weights[i] + (1.0 - alpha) * _physical_volfrac;
98  if (_new_weights[i] < 0.0) // We can end up with small (roundoff) negative weights
99  _new_weights[i] = 0.0;
100  }
101  }
102  break;
103  }
104  case Xfem::DIRECT: // remove q-points outside the partial element's physical domain
105  {
106  bool nonzero = false;
107  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
108  {
109  // q_points contains quadrature point locations in physical coordinates
110  if (isPointPhysical(q_points[qp]))
111  {
112  nonzero = true;
113  _new_weights[qp] = 1.0;
114  }
115  else
116  _new_weights[qp] = 0.0;
117  }
118  if (!nonzero)
119  {
120  // Set the weights to a small value (1e-3) to avoid having DOFs
121  // with zeros on the diagonal, which occurs for nodes that are
122  // connected to no physical material.
123  for (unsigned qp = 0; qp < qrule->n_points(); ++qp)
124  _new_weights[qp] = 1e-3;
125  }
126  break;
127  }
128  default:
129  mooseError("Undefined option for XFEM_QRULE");
130  }
131  _have_weights = true;
132 }
unsigned int _n_qpoints
Definition: XFEMCutElem.h:36
std::vector< Point > _qp_points
Definition: XFEMCutElem.h:38
Real _physical_volfrac
Definition: XFEMCutElem.h:41
Real getPhysicalVolumeFraction() const
Definition: XFEMCutElem.C:31
bool isPointPhysical(const Point &p) const
Definition: XFEMCutElem.C:135
std::vector< Real > _qp_weights
Definition: XFEMCutElem.h:39
std::vector< Real > _new_weights
Definition: XFEMCutElem.h:43
bool _have_weights
Definition: XFEMCutElem.h:42
virtual void computePhysicalVolumeFraction()=0
virtual void computeMomentFittingWeights()=0
virtual void XFEMCutElem::getCrackTipOriginAndDirection ( unsigned  tip_id,
Point &  origin,
Point &  direction 
) const
pure virtual
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().

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

virtual const EFAElement* XFEMCutElem::getEFAElement ( ) const
pure virtual
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().

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

Real XFEMCutElem::getMomentFittingWeight ( unsigned int  i_qp) const
virtual Point XFEMCutElem::getNodeCoordinates ( EFANode node,
MeshBase *  displaced_mesh = NULL 
) const
protectedpure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Real XFEMCutElem::getPhysicalVolumeFraction ( ) const

Definition at line 31 of file XFEMCutElem.C.

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

32 {
33  return _physical_volfrac;
34 }
Real _physical_volfrac
Definition: XFEMCutElem.h:41
void XFEMCutElem::getWeightMultipliers ( MooseArray< Real > &  weights,
QBase *  qrule,
Xfem::XFEM_QRULE  xfem_qrule,
const MooseArray< Point > &  q_points 
)

Definition at line 37 of file XFEMCutElem.C.

Referenced by XFEM::getXFEMWeights().

41 {
42  if (!_have_weights)
43  computeXFEMWeights(qrule, xfem_qrule, q_points);
44 
45  weights.resize(_new_weights.size());
46  for (unsigned int qp = 0; qp < _new_weights.size(); ++qp)
47  weights[qp] = _new_weights[qp];
48 }
std::vector< Real > _new_weights
Definition: XFEMCutElem.h:43
bool _have_weights
Definition: XFEMCutElem.h:42
void computeXFEMWeights(QBase *qrule, Xfem::XFEM_QRULE xfem_qrule, const MooseArray< Point > &q_points)
Definition: XFEMCutElem.C:51
bool XFEMCutElem::isPointPhysical ( const Point &  p) const

Definition at line 135 of file XFEMCutElem.C.

Referenced by computeXFEMWeights().

136 {
137  // determine whether the point is inside the physical domain of a partial element
138  bool physical_flag = true;
139  unsigned int n_cut_planes = numCutPlanes();
140  for (unsigned int plane_id = 0; plane_id < n_cut_planes; ++plane_id)
141  {
142  Point origin = getCutPlaneOrigin(plane_id);
143  Point normal = getCutPlaneNormal(plane_id);
144  Point origin2qp = p - origin;
145  if (origin2qp * normal > 0.0)
146  {
147  physical_flag = false; // Point outside pysical domain
148  break;
149  }
150  }
151  return physical_flag;
152 }
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 unsigned int numCutPlanes() const =0
virtual unsigned int XFEMCutElem::numCutPlanes ( ) const
pure virtual

Implemented in XFEMCutElem2D, and XFEMCutElem3D.

Referenced by isPointPhysical().

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

Member Data Documentation

Real XFEMCutElem::_elem_volume
protected
bool XFEMCutElem::_have_weights
protected

Definition at line 42 of file XFEMCutElem.h.

Referenced by computeXFEMWeights(), and getWeightMultipliers().

unsigned int XFEMCutElem::_n_nodes
protected

Definition at line 35 of file XFEMCutElem.h.

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

unsigned int XFEMCutElem::_n_qpoints
protected
std::vector<Real> XFEMCutElem::_new_weights
protected
std::vector<Node *> XFEMCutElem::_nodes
protected
Real XFEMCutElem::_physical_volfrac
protected
std::vector<Point> XFEMCutElem::_qp_points
protected

Definition at line 38 of file XFEMCutElem.h.

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

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

Definition at line 39 of file XFEMCutElem.h.

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


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