libMesh::QComposite< QSubCell > Class Template Reference

A quadrature rule for subdivided elements. More...

#include <quadrature_composite.h>

Inheritance diagram for libMesh::QComposite< QSubCell >:

Public Member Functions

 QComposite (const unsigned int _dim, const Order _order=INVALID_ORDER)
 
 ~QComposite ()
 
virtual QuadratureType type () const libmesh_override
 
virtual void init (const Elem &elem, const std::vector< Real > &vertex_distance_func, unsigned int p_level=0) libmesh_override
 

Private Member Functions

void add_subelem_values (const std::vector< Elem const * > &subelem)
 

Private Attributes

QSubCell _q_subcell
 
ElemCutter _elem_cutter
 
UniquePtr< FEBase_lagrange_fe
 

Detailed Description

template<class QSubCell>
class libMesh::QComposite< QSubCell >

A quadrature rule for subdivided elements.

This class implements generic composite quadrature rules. Composite quadrature rules are constructed from any of the supported rules by breaking an element into subelements and applying the base rule on each subelement. This class uses the ElemCutter, which is only available if libmesh is configured with –disable-strict-lgpl.

Author
Benjamin Kirk
Date
2013

Definition at line 49 of file quadrature_composite.h.

Constructor & Destructor Documentation

template<class QSubCell >
libMesh::QComposite< QSubCell >::QComposite ( const unsigned int  _dim,
const Order  _order = INVALID_ORDER 
)

Constructor. Declares the order of the quadrature rule.

Definition at line 36 of file quadrature_composite.C.

References libMesh::QComposite< QSubCell >::_lagrange_fe, libMesh::QComposite< QSubCell >::_q_subcell, libMesh::EDGE2, libMesh::TriangleWrapper::init(), libMesh::libmesh_assert(), and libmesh_nullptr.

37  :
38  QSubCell(d,o), // explicitly call base class constructor
39  _q_subcell(d,o),
41 {
42  // explicitly call the init function in 1D since the
43  // other tensor-product rules require this one.
44  // note that EDGE will not be used internally, however
45  // if we called the function with INVALID_ELEM it would try to
46  // be smart and return, thinking it had already done the work.
47  if (_dim == 1)
49 
51 
52  _lagrange_fe->attach_quadrature_rule (&_q_subcell);
53 }
const class libmesh_nullptr_t libmesh_nullptr
libmesh_assert(j)
UniquePtr< FEBase > _lagrange_fe
void init(triangulateio &t)
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
template<class QSubCell >
libMesh::QComposite< QSubCell >::~QComposite ( )

Destructor.

Definition at line 58 of file quadrature_composite.C.

59 {}

Member Function Documentation

template<class QSubCell >
void libMesh::QComposite< QSubCell >::add_subelem_values ( const std::vector< Elem const * > &  subelem)
private

Helper function called from init() to collect all the points and weights of the subelement quadrature rules.

Definition at line 120 of file quadrature_composite.C.

References libMesh::QComposite< QSubCell >::_lagrange_fe, and libMesh::err.

Referenced by libMesh::QComposite< QSubCell >::init(), and libMesh::QComposite< QSubCell >::type().

122 {
123  const std::vector<Real> & subelem_weights = _lagrange_fe->get_JxW();
124  const std::vector<Point> & subelem_points = _lagrange_fe->get_xyz();
125 
126  for (std::vector<Elem const *>::const_iterator it = subelem.begin();
127  it!=subelem.end(); ++it)
128  {
129  // tetgen seems to create 0-volume cells on occasion, but we *should*
130  // be catching that appropriately now inside the ElemCutter class.
131  // Just in case trap here, describe the error, and abort.
132  libmesh_try
133  {
134  _lagrange_fe->reinit(*it);
135  _weights.insert(_weights.end(),
136  subelem_weights.begin(), subelem_weights.end());
137 
138  _points.insert(_points.end(),
139  subelem_points.begin(), subelem_points.end());
140  }
141  libmesh_catch (...)
142  {
143  libMesh::err << "ERROR: found a bad cut cell!\n";
144 
145  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
146  libMesh::err << (*it)->point(n) << std::endl;
147 
148  libmesh_error_msg("Tetgen may have created a 0-volume cell during Cutcell integration.");
149  }
150  }
151 }
UniquePtr< FEBase > _lagrange_fe
OStreamProxy err(std::cerr)
template<class QSubCell >
void libMesh::QComposite< QSubCell >::init ( const Elem elem,
const std::vector< Real > &  vertex_distance_func,
unsigned int  p_level = 0 
)
virtual

Overrides the base class init() function, and uses the ElemCutter to subdivide the element into "inside" and "outside" subelements.

Definition at line 64 of file quadrature_composite.C.

References libMesh::QComposite< QSubCell >::_elem_cutter, libMesh::QComposite< QSubCell >::_q_subcell, libMesh::QComposite< QSubCell >::add_subelem_values(), libMesh::Elem::dim(), libMesh::ElemCutter::inside_elements(), libMesh::ElemCutter::is_cut(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::Elem::n_vertices(), libMesh::ElemCutter::outside_elements(), libMesh::Elem::reference_elem(), and libMesh::Elem::type().

Referenced by libMesh::QComposite< QSubCell >::type().

67 {
68  libmesh_assert_equal_to (vertex_distance_func.size(), elem.n_vertices());
69  libmesh_assert_equal_to (_dim, elem.dim());
70 
71  // if we are not cut, revert to simple base class init() method.
72  if (!_elem_cutter.is_cut (elem, vertex_distance_func))
73  {
74  _q_subcell.init (elem.type(), p_level);
75  _points = _q_subcell.get_points();
76  _weights = _q_subcell.get_weights();
77 
78  //this->print_info();
79  return;
80  }
81 
82  // Get a pointer to the element's reference element. We want to
83  // perform cutting on the reference element such that the quadrature
84  // point locations of the subelements live in the reference
85  // coordinate system, thereby eliminating the need for inverse
86  // mapping.
87  const Elem * reference_elem = elem.reference_elem();
88 
89  libmesh_assert (reference_elem != libmesh_nullptr);
90 
91  _elem_cutter(*reference_elem, vertex_distance_func);
92  //_elem_cutter(elem, vertex_distance_func);
93 
94  // clear our state & accumulate points from subelements
95  _points.clear();
96  _weights.clear();
97 
98  // inside subelem
99  {
100  const std::vector<Elem const *> & inside_elem (_elem_cutter.inside_elements());
101  std::cout << inside_elem.size() << " elements inside\n";
102 
103  this->add_subelem_values(inside_elem);
104  }
105 
106  // outside subelem
107  {
108  const std::vector<Elem const *> & outside_elem (_elem_cutter.outside_elements());
109  std::cout << outside_elem.size() << " elements outside\n";
110 
111  this->add_subelem_values(outside_elem);
112  }
113 
114  this->print_info();
115 }
const class libmesh_nullptr_t libmesh_nullptr
const std::vector< Elem const * > & inside_elements() const
Definition: elem_cutter.h:115
bool is_cut(const Elem &elem, const std::vector< Real > &vertex_distance_func) const
Definition: elem_cutter.C:89
libmesh_assert(j)
void add_subelem_values(const std::vector< Elem const * > &subelem)
const std::vector< Elem const * > & outside_elements() const
Definition: elem_cutter.h:123
template<class QSubCell >
virtual QuadratureType libMesh::QComposite< QSubCell >::type ( ) const
inlinevirtual

Member Data Documentation

template<class QSubCell >
ElemCutter libMesh::QComposite< QSubCell >::_elem_cutter
private

ElemCutter object.

Definition at line 101 of file quadrature_composite.h.

Referenced by libMesh::QComposite< QSubCell >::init().

template<class QSubCell >
UniquePtr<FEBase> libMesh::QComposite< QSubCell >::_lagrange_fe
private

Lagrange FE to use for subcell mapping.

Definition at line 106 of file quadrature_composite.h.

Referenced by libMesh::QComposite< QSubCell >::add_subelem_values(), and libMesh::QComposite< QSubCell >::QComposite().

template<class QSubCell >
QSubCell libMesh::QComposite< QSubCell >::_q_subcell
private

Subcell quadrature object.

Definition at line 96 of file quadrature_composite.h.

Referenced by libMesh::QComposite< QSubCell >::init(), and libMesh::QComposite< QSubCell >::QComposite().


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