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>
38template <
class Scalar,
unsigned dim>
42 enum { numCorners = (1 << dim) };
44 using LocalPosition = Dune::FieldVector<Scalar, dim>;
45 using GlobalPosition = Dune::FieldVector<Scalar, dim>;
47 Dune::GeometryType type()
const
48 {
return Dune::GeometryType((1 << dim) - 1, dim); }
50 template <
class CornerContainer>
55 for (
unsigned j = 0; j < dim; ++j)
75 GlobalPosition
global(
const LocalPosition& localPos)
const
77 GlobalPosition globalPos(0.0);
91 const LocalPosition& localPos)
const
95 for (
unsigned k = 0;
k < dim; ++
k) {
97 for (
unsigned j = 0; j < dim; ++j) {
119 Dune::FieldMatrix<Scalar, dim, dim>
jac;
121 return jac.determinant();
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_;
Quadrature geometry for quadrilaterals.
Definition quadraturegeometries.hh:40
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
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
const GlobalPosition & corner(unsigned cornerIdx) const
Return the position of the corner with a given index.
Definition quadraturegeometries.hh:127
GlobalPosition global(const LocalPosition &localPos) const
Convert a local coordinate into a global one.
Definition quadraturegeometries.hh:75
const GlobalPosition & center() const
Returns the center of weight of the polyhedron.
Definition quadraturegeometries.hh:69
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
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242