20 #include "libmesh/quadrature_gauss.h" 21 #include "libmesh/enum_quadrature_type.h" 40 return std::make_unique<QGauss>(*this);
44 const unsigned int n_pts)
55 unsigned int offset=0;
58 for (
unsigned int p=0; p<n_pts; ++p)
62 libmesh_assert_not_equal_to (rule_data[p][0], static_cast<Real>(0.0));
65 libmesh_assert_not_equal_to (rule_data[p][3], static_cast<Real>(0.0));
71 unsigned int pointtype=1;
73 if (rule_data[p][1] != static_cast<Real>(0.0))
75 if (rule_data[p][2] != static_cast<Real>(0.0))
83 if (rule_data[p][2] != static_cast<Real>(0.0))
93 libmesh_assert_less (offset + 0,
_points.size());
95 const Real a = rule_data[p][0];
101 _weights[offset + 0] = rule_data[p][3];
110 libmesh_assert_less (offset + 3,
_points.size());
112 const Real a = rule_data[p][0];
113 const Real b = rule_data[p][1];
114 const Real wt = rule_data[p][3];
123 for (
unsigned int j=0; j<pointtype; ++j)
133 libmesh_assert_less (offset + 5,
_points.size());
135 const Real a = rule_data[p][0];
136 const Real b = rule_data[p][2];
137 const Real wt = rule_data[p][3];
147 for (
unsigned int j=0; j<pointtype; ++j)
158 libmesh_assert_less (offset + 11,
_points.size());
160 const Real a = rule_data[p][0];
161 const Real b = rule_data[p][1];
162 const Real c = rule_data[p][2];
163 const Real wt = rule_data[p][3];
173 for (
unsigned int j=0; j<pointtype; ++j)
181 libmesh_error_msg(
"Don't know what to do with this many permutation points!");
198 const unsigned int * permutation_ids,
204 unsigned int total_pts = 0;
205 for (
unsigned int p=0; p<n_wts; ++p)
206 total_pts += permutation_ids[p];
213 unsigned int offset=0;
215 for (
unsigned int p=0; p<n_wts; ++p)
217 switch (permutation_ids[p])
238 for (
unsigned int j=0; j<3; ++j)
255 for (
unsigned int j=0; j<6; ++j)
263 libmesh_error_msg(
"Unknown permutation id: " << permutation_ids[p] <<
"!");
271 const unsigned int n_pts)
280 unsigned int offset=0;
283 for (
unsigned int p=0; p<n_pts; ++p)
287 libmesh_assert_not_equal_to ( rule_data[p][0], static_cast<Real>(0.0) );
290 libmesh_assert_not_equal_to ( rule_data[p][3], static_cast<Real>(0.0) );
296 unsigned int pointtype=1;
298 if (rule_data[p][1] != static_cast<Real>(0.0))
300 if (rule_data[p][2] != static_cast<Real>(0.0))
311 libmesh_assert_less (offset + 0,
_points.size());
314 _points[offset + 0] =
Point(rule_data[p][0], rule_data[p][0]);
317 _weights[offset + 0] = rule_data[p][3];
326 libmesh_assert_less (offset + 2,
_points.size());
330 _points[offset + 0] =
Point(rule_data[p][0], rule_data[p][1]);
331 _points[offset + 1] =
Point(rule_data[p][1], rule_data[p][0]);
332 _points[offset + 2] =
Point(rule_data[p][1], rule_data[p][1]);
334 for (
unsigned int j=0; j<3; ++j)
335 _weights[offset + j] = rule_data[p][3];
344 libmesh_assert_less (offset + 5,
_points.size());
347 _points[offset + 0] =
Point(rule_data[p][0], rule_data[p][1]);
348 _points[offset + 1] =
Point(rule_data[p][0], rule_data[p][2]);
349 _points[offset + 2] =
Point(rule_data[p][1], rule_data[p][0]);
350 _points[offset + 3] =
Point(rule_data[p][1], rule_data[p][2]);
351 _points[offset + 4] =
Point(rule_data[p][2], rule_data[p][0]);
352 _points[offset + 5] =
Point(rule_data[p][2], rule_data[p][1]);
354 for (
unsigned int j=0; j<6; ++j)
355 _weights[offset + j] = rule_data[p][3];
362 libmesh_error_msg(
"Don't know what to do with this many permutation points!");
virtual QuadratureType type() const override
QuadratureType
Defines an enum for currently available quadrature rules.
The libMesh namespace provides an interface to certain functionality in the library.
std::vector< Point > _points
The locations of the quadrature points in reference element space.
std::vector< Real > _weights
The quadrature weights.
void dunavant_rule2(const Real *wts, const Real *a, const Real *b, const unsigned int *permutation_ids, const unsigned int n_wts)
virtual std::unique_ptr< QBase > clone() const override
void dunavant_rule(const Real rule_data[][4], const unsigned int n_pts)
The Dunavant rules are for triangles.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void keast_rule(const Real rule_data[][4], const unsigned int n_pts)
The Keast rules are for tets.
A Point defines a location in LIBMESH_DIM dimensional Real space.