libMesh
Public Member Functions | Private Member Functions | Private Attributes | List of all members
libMesh::QComposite< QSubCell > Class Template Referencefinal

This class implements generic composite quadrature rules. More...

#include <quadrature_composite.h>

Inheritance diagram for libMesh::QComposite< QSubCell >:
[legend]

Public Member Functions

 QComposite (unsigned int dim, Order order=INVALID_ORDER)
 Constructor. More...
 
 QComposite (const QComposite &)=delete
 This class contains a unique_ptr member, so it can't be default copy constructed or assigned. More...
 
QCompositeoperator= (const QComposite &)=delete
 
 QComposite (QComposite &&)=default
 Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this simple class. More...
 
QCompositeoperator= (QComposite &&)=default
 
virtual ~QComposite ()=default
 
virtual QuadratureType type () const override
 
virtual void init (const Elem &elem, const std::vector< Real > &vertex_distance_func, unsigned int p_level=0) override
 Overrides the base class init() function, and uses the ElemCutter to subdivide the element into "inside" and "outside" subelements. More...
 

Private Member Functions

void add_subelem_values (const std::vector< Elem const *> &subelem)
 Helper function called from init() to collect all the points and weights of the subelement quadrature rules. More...
 

Private Attributes

QSubCell _q_subcell
 Subcell quadrature object. More...
 
ElemCutter _elem_cutter
 ElemCutter object. More...
 
std::unique_ptr< FEBase_lagrange_fe
 Lagrange FE to use for subcell mapping. More...
 

Detailed Description

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

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 A quadrature rule for subdivided elements.

Definition at line 51 of file quadrature_composite.h.

Constructor & Destructor Documentation

◆ QComposite() [1/3]

template<class QSubCell >
libMesh::QComposite< QSubCell >::QComposite ( unsigned int  dim,
Order  order = INVALID_ORDER 
)

Constructor.

Declares the order of the quadrature rule.

Definition at line 45 of file quadrature_composite.C.

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

46  :
47  QSubCell(d,o), // explicitly call base class constructor
48  _q_subcell(d,o),
50 {
51  // explicitly call the init function in 1D since the
52  // other tensor-product rules require this one.
53  // note that EDGE will not be used internally, however
54  // if we called the function with INVALID_ELEM it would try to
55  // be smart and return, thinking it had already done the work.
56  if (_dim == 1)
58 
59  libmesh_assert (_lagrange_fe.get() != nullptr);
60 
61  _lagrange_fe->attach_quadrature_rule (&_q_subcell);
62 }
QSubCell _q_subcell
Subcell quadrature object.
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
libmesh_assert(ctx)
std::unique_ptr< FEBase > _lagrange_fe
Lagrange FE to use for subcell mapping.

◆ QComposite() [2/3]

template<class QSubCell>
libMesh::QComposite< QSubCell >::QComposite ( const QComposite< QSubCell > &  )
delete

This class contains a unique_ptr member, so it can't be default copy constructed or assigned.

◆ QComposite() [3/3]

template<class QSubCell>
libMesh::QComposite< QSubCell >::QComposite ( QComposite< QSubCell > &&  )
default

Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this simple class.

◆ ~QComposite()

template<class QSubCell>
virtual libMesh::QComposite< QSubCell >::~QComposite ( )
virtualdefault

Member Function Documentation

◆ add_subelem_values()

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 127 of file quadrature_composite.C.

References libMesh::err.

129 {
130  const std::vector<Real> & subelem_weights = _lagrange_fe->get_JxW();
131  const std::vector<Point> & subelem_points = _lagrange_fe->get_xyz();
132 
133  for (const auto & elem : subelem)
134  {
135  // tetgen seems to create 0-volume cells on occasion, but we *should*
136  // be catching that appropriately now inside the ElemCutter class.
137  // Just in case trap here, describe the error, and abort.
138  libmesh_try
139  {
140  _lagrange_fe->reinit(elem);
141  _weights.insert(_weights.end(),
142  subelem_weights.begin(), subelem_weights.end());
143 
144  _points.insert(_points.end(),
145  subelem_points.begin(), subelem_points.end());
146  }
147  libmesh_catch (...)
148  {
149  libMesh::err << "ERROR: found a bad cut cell!\n";
150 
151  for (auto n : elem->node_index_range())
152  libMesh::err << elem->point(n) << std::endl;
153 
154  libmesh_error_msg("Tetgen may have created a 0-volume cell during Cutcell integration.");
155  }
156  }
157 }
OStreamProxy err
std::unique_ptr< FEBase > _lagrange_fe
Lagrange FE to use for subcell mapping.

◆ init()

template<class QSubCell >
void libMesh::QComposite< QSubCell >::init ( const Elem elem,
const std::vector< Real > &  vertex_distance_func,
unsigned int  p_level = 0 
)
overridevirtual

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

Definition at line 67 of file quadrature_composite.C.

References libMesh::Elem::dim(), libMesh::LAGRANGE_MAP, libMesh::libmesh_assert(), libMesh::Elem::mapping_type(), libMesh::Elem::n_vertices(), libMesh::Elem::reference_elem(), and libMesh::Elem::type().

Referenced by integrate_function().

70 {
71  libmesh_assert_equal_to (vertex_distance_func.size(), elem.n_vertices());
72  libmesh_assert_equal_to (_dim, elem.dim());
73 
74  // We already initialized a Lagrange map; we're not supporting
75  // others yet
76  libmesh_assert_equal_to (elem.mapping_type(), LAGRANGE_MAP);
77 
78  // if we are not cut, revert to simple base class init() method.
79  if (!_elem_cutter.is_cut (elem, vertex_distance_func))
80  {
81  _q_subcell.init (elem.type(), p_level);
82  _points = _q_subcell.get_points();
83  _weights = _q_subcell.get_weights();
84 
85  //this->print_info();
86  return;
87  }
88 
89  // Get a pointer to the element's reference element. We want to
90  // perform cutting on the reference element such that the quadrature
91  // point locations of the subelements live in the reference
92  // coordinate system, thereby eliminating the need for inverse
93  // mapping.
94  const Elem * reference_elem = elem.reference_elem();
95 
96  libmesh_assert (reference_elem != nullptr);
97 
98  _elem_cutter(*reference_elem, vertex_distance_func);
99  //_elem_cutter(elem, vertex_distance_func);
100 
101  // clear our state & accumulate points from subelements
102  _points.clear();
103  _weights.clear();
104 
105  // inside subelem
106  {
107  const std::vector<Elem const *> & inside_elem (_elem_cutter.inside_elements());
108  // std::cout << inside_elem.size() << " elements inside\n";
109 
110  this->add_subelem_values(inside_elem);
111  }
112 
113  // outside subelem
114  {
115  const std::vector<Elem const *> & outside_elem (_elem_cutter.outside_elements());
116  // std::cout << outside_elem.size() << " elements outside\n";
117 
118  this->add_subelem_values(outside_elem);
119  }
120 
121  // this->print_info();
122 }
QSubCell _q_subcell
Subcell quadrature object.
bool is_cut(const Elem &elem, const std::vector< Real > &vertex_distance_func) const
Definition: elem_cutter.C:89
const std::vector< Elem const * > & outside_elements() const
Definition: elem_cutter.h:125
libmesh_assert(ctx)
const std::vector< Elem const * > & inside_elements() const
Definition: elem_cutter.h:117
ElemCutter _elem_cutter
ElemCutter object.
void add_subelem_values(const std::vector< Elem const *> &subelem)
Helper function called from init() to collect all the points and weights of the subelement quadrature...

◆ operator=() [1/2]

template<class QSubCell>
QComposite& libMesh::QComposite< QSubCell >::operator= ( const QComposite< QSubCell > &  )
delete

◆ operator=() [2/2]

template<class QSubCell>
QComposite& libMesh::QComposite< QSubCell >::operator= ( QComposite< QSubCell > &&  )
default

◆ type()

template<class QSubCell >
QuadratureType libMesh::QComposite< QSubCell >::type ( ) const
overridevirtual
Returns
QCOMPOSITE.

Definition at line 37 of file quadrature_composite.C.

References libMesh::QCOMPOSITE.

38 {
39  return QCOMPOSITE;
40 }

Member Data Documentation

◆ _elem_cutter

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

ElemCutter object.

Definition at line 113 of file quadrature_composite.h.

◆ _lagrange_fe

template<class QSubCell>
std::unique_ptr<FEBase> libMesh::QComposite< QSubCell >::_lagrange_fe
private

Lagrange FE to use for subcell mapping.

Definition at line 118 of file quadrature_composite.h.

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

◆ _q_subcell

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

Subcell quadrature object.

Definition at line 108 of file quadrature_composite.h.

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


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