27 #ifndef EWOMS_QUADRATURE_GEOMETRIES_HH
28 #define EWOMS_QUADRATURE_GEOMETRIES_HH
30 #include <dune/common/fvector.hh>
31 #include <dune/common/fmatrix.hh>
32 #include <dune/geometry/type.hh>
38 template <
class Scalar,
unsigned dim>
42 enum { numCorners = (1 << dim) };
44 typedef Dune::FieldVector<Scalar, dim> LocalPosition;
45 typedef Dune::FieldVector<Scalar, dim> GlobalPosition;
47 Dune::GeometryType type()
const
48 {
return Dune::GeometryType(Dune::GeometryType::cube, dim); }
50 template <
class CornerContainer>
51 void setCorners(
const CornerContainer& corners,
unsigned numCorners)
54 for (cornerIdx = 0; cornerIdx < numCorners; ++cornerIdx) {
55 for (
unsigned j = 0; j < dim; ++j)
56 corners_[cornerIdx][j] = corners[cornerIdx][j];
58 assert(cornerIdx == numCorners);
61 for (cornerIdx = 0; cornerIdx < numCorners; ++cornerIdx)
62 center_ += corners_[cornerIdx];
63 center_ /= numCorners;
75 GlobalPosition
global(
const LocalPosition& localPos)
const
77 GlobalPosition globalPos(0.0);
79 for (
unsigned cornerIdx = 0; cornerIdx < numCorners; ++cornerIdx)
90 void jacobian(Dune::FieldMatrix<Scalar, dim, dim>& jac,
91 const LocalPosition& localPos)
const
94 for (
unsigned cornerIdx = 0; cornerIdx < numCorners; ++cornerIdx) {
95 for (
unsigned k = 0; k < dim; ++k) {
96 Scalar dWeight_dk = (cornerIdx& (1 << k)) ? 1 : -1;
97 for (
unsigned j = 0; j < dim; ++j) {
99 if (cornerIdx& (1 << j))
100 dWeight_dk *= localPos[j];
102 dWeight_dk *= 1 - localPos[j];
107 jac[k].axpy(dWeight_dk, corners_[cornerIdx]);
119 Dune::FieldMatrix<Scalar, dim, dim> jac;
121 return jac.determinant();
127 const GlobalPosition&
corner(
unsigned cornerIdx)
const
128 {
return corners_[cornerIdx]; }
134 Scalar
cornerWeight(
const LocalPosition& localPos,
unsigned cornerIdx)
const
136 GlobalPosition globalPos(0.0);
141 for (
unsigned j = 0; j < dim; ++j)
142 weight *= (cornerIdx& (1 << j)) ? localPos[j] : (1 - localPos[j]);
148 GlobalPosition corners_[numCorners];
149 GlobalPosition center_;
154 #endif // EWOMS_QUADRATURE_GEOMETRY_HH
Quadrature geometry for quadrilaterals.
Definition: quadraturegeometries.hh:39
void jacobian(Dune::FieldMatrix< Scalar, dim, dim > &jac, const LocalPosition &localPos) const
Returns the Jacobian matrix of the local to global mapping at a given local position.
Definition: quadraturegeometries.hh:90
Scalar cornerWeight(const LocalPosition &localPos, unsigned cornerIdx) const
Return the weight of an individual corner for the local to global mapping.
Definition: quadraturegeometries.hh:134
GlobalPosition global(const LocalPosition &localPos) const
Convert a local coordinate into a global one.
Definition: quadraturegeometries.hh:75
Scalar integrationElement(const LocalPosition &localPos) const
Return the determinant of the Jacobian of the mapping from local to global coordinates at a given loc...
Definition: quadraturegeometries.hh:117
const GlobalPosition & center() const
Returns the center of weight of the polyhedron.
Definition: quadraturegeometries.hh:69
const GlobalPosition & corner(unsigned cornerIdx) const
Return the position of the corner with a given index.
Definition: quadraturegeometries.hh:127