1#ifndef DUNE_FEM_FEMQUADRATURES_INLINE_HH
2#define DUNE_FEM_FEMQUADRATURES_INLINE_HH
12#define SimplexPointsAdapter ParDGSimplexQuadrature::Quadrature
15 template <
class ct,
int dim>
18 order_((order <= 0) ? 1 : order)
21 assert(
order_ <= SimplexMaxOrder :: maxOrder( dim ) );
23 SimplexPointsAdapter<dim> points(
order);
27 for (
int i = 0; i < points.numPoints(); ++i) {
32 template <
class ct,
int dim>
35 order_((order <= 0) ? 1 : order)
52 if (m==0) DUNE_THROW(NotImplemented,
"order " <<
order_ <<
" not implemented");
56 int n = gp.
power(m,dim);
57 for (
int i = 0; i < n; i++) {
61 for (
int k=0; k<dim; ++k) {
69 for (
int k = 0; k < dim; k++)
71 local[k] = gp.
point(m,x[k]);
84 order_((order <= 0) ? 1 : order)
104 order_((order <= 0) ? 1 : order)
106 SimplexPointsAdapter<2> simplexPoints(
order);
107 int simplexOrder = simplexPoints.order();
114 if (gp.
order(i)>=order_) {
119 if (m==0) DUNE_THROW(NotImplemented,
"order not implemented");
121 int gaussOrder = gp.
order(m);
122 int minOrder = ((simplexOrder < gaussOrder) ? simplexOrder : gaussOrder);
125 int numSimplexPoints = simplexPoints.numPoints();
126 int numGaussPoints = gp.
power(m,1);
128 FieldVector<ct, 3> local;
129 double weight, simplexWeight;
131 for (
int i = 0; i < numSimplexPoints; ++i) {
132 local[0] = simplexPoints.point(i)[0];
133 local[1] = simplexPoints.point(i)[1];
134 simplexWeight = simplexPoints.weight(i);
135 for (
int j = 0; j < numGaussPoints; ++j) {
136 local[2] = gp.
point(m,j);
147 order_((order <= 0) ? 1 : order)
153 if (points.
order(i) >= order_) {
159 if (m==0) DUNE_THROW(NotImplemented,
"order not implemented");
160 order_ = points.
order(m);
163 for (
int i = 0; i < points.
numPoints(m); ++i) {
169 template <
class ct,
int dim>
172 geometryType_( geomType ),
173 order_((order <= 0) ? 1 : order)
Definition: bindguard.hh:11
int order_
Definition: femquadratures.hh:66
virtual int order() const
obtain order of the integration point list
Definition: femquadratures.hh:86
SimplexQuadrature(const GeometryType &geometry, int order, size_t id)
constructor filling the list of points and weights
Definition: femquadratures_inline.hh:16
CubeQuadrature(const GeometryType &geometry, int order, size_t id)
constructor filling the list of points and weights
Definition: femquadratures_inline.hh:33
int order_
Definition: femquadratures.hh:132
BaseType::CoordinateType CoordinateType
type of local coordinates
Definition: femquadratures.hh:129
PrismQuadrature(const GeometryType &geometry, int order, size_t id)
constructor filling the list of points and weights
Definition: femquadratures_inline.hh:102
virtual int order() const
obtain order of the integration point list
Definition: femquadratures.hh:208
PyramidQuadrature(const GeometryType &geometry, int order, size_t id)
constructor filling the list of points and weights
Definition: femquadratures_inline.hh:145
PolyhedronQuadrature(const GeometryType &geometry, int order, size_t id)
constructor filling the list of points and weights
Definition: femquadratures_inline.hh:170
one-dimensional Gauss points and their weights
Definition: gausspoints.hh:26
double point(int m, int i) const
obtain the i-th point of the m-th quadratre
Definition: gausspoints.hh:51
@ highestOrder
Definition: gausspoints.hh:32
int order(int m) const
obtain the order of the m-th quadratre
Definition: gausspoints.hh:76
@ MAXP
Definition: gausspoints.hh:29
int power(int y, int d) const
a simple power method
Definition: gausspoints.hh:90
double weight(int m, int i) const
obtain the i-th weight of the m-th quadratre
Definition: gausspoints.hh:64
Definition: pyramidpoints.hh:17
const FieldVector< double, 3 > & point(int m, int i) const
Access to the ith point of quadrature rule m.
Definition: pyramidpoints.hh:24
int numPoints(int m) const
Number of points in the quadrature rule m.
Definition: pyramidpoints.hh:42
int order(int m) const
Actual order of quadrature rule m.
Definition: pyramidpoints.hh:36
double weight(int m, int i) const
Access to the ith weight of quadrature rule m.
Definition: pyramidpoints.hh:30
@ numQuads
Definition: pyramidpoints.hh:19
Generic implementation of a Dune quadrature.
Definition: quadratureimp.hh:196
void addQuadraturePoint(const CoordinateType &point, const FieldType weight)
Adds a point-weight pair to the quadrature.
Definition: quadratureimp.hh:270
const FieldType & weight(size_t i) const
obtain weight of i-th integration point
Definition: quadratureimp.hh:251