6#ifndef DUNE_GEOMETRY_AXISALIGNED_CUBE_GEOMETRY_HH
7#define DUNE_GEOMETRY_AXISALIGNED_CUBE_GEOMETRY_HH
15#include <dune/common/fvector.hh>
16#include <dune/common/fmatrix.hh>
17#include <dune/common/diagonalmatrix.hh>
48 template <
class CoordType,
unsigned int dim,
unsigned int coorddim>
79 typedef typename std::conditional<dim==coorddim,
80 DiagonalMatrix<ctype,dim>,
89 typedef typename std::conditional<dim==coorddim,
90 DiagonalMatrix<ctype,dim>,
100 using Jacobian = std::conditional_t<dim==coorddim, DiagonalMatrix<ctype,dim>, FieldMatrix<ctype,coorddim,dim> >;
109 using JacobianInverse = std::conditional_t<dim==coorddim, DiagonalMatrix<ctype,dim>, FieldMatrix<ctype,dim,coorddim> >;
116 const Dune::FieldVector<ctype,coorddim> upper)
121 static_assert(dim==coorddim,
"Use this constructor only if dim==coorddim!");
123 axes_ = (1<<coorddim)-1;
134 const Dune::FieldVector<ctype,coorddim> upper,
135 const std::bitset<coorddim>& axes)
140 assert(axes.count()==dim);
141 for (
size_t i=0; i<coorddim; i++)
143 upper_[i] = lower_[i];
164 if (dim == coorddim) {
165 for (
size_t i=0; i<coorddim; i++)
166 result[i] = lower_[i] +
local[i]*(upper_[i] - lower_[i]);
167 }
else if (dim == 0) {
171 for (
size_t i=0; i<coorddim; i++)
172 result[i] = (axes_[i])
173 ? lower_[i] +
local[lc++]*(upper_[i] - lower_[i])
183 if (dim == coorddim) {
184 for (
size_t i=0; i<dim; i++)
185 result[i] = (
global[i] - lower_[i]) / (upper_[i] - lower_[i]);
186 }
else if (dim != 0) {
188 for (
size_t i=0; i<coorddim; i++)
190 result[lc++] = (
global[i] - lower_[i]) / (upper_[i] - lower_[i]);
247 for (
size_t i=0; i<coorddim; i++)
248 result[i] = CoordType(0.5) * (lower_[i] + upper_[i]);
263 if (dim == coorddim) {
264 for (
size_t i=0; i<coorddim; i++)
265 result[i] = (k & (1<<i)) ? upper_[i] : lower_[i];
266 }
else if (dim == 0) {
269 unsigned int mask = 1;
271 for (
size_t i=0; i<coorddim; i++) {
273 result[i] = lower_[i];
275 result[i] = (k & mask) ? upper_[i] : lower_[i];
289 if (dim == coorddim) {
290 for (
size_t i=0; i<dim; i++)
291 vol *= upper_[i] - lower_[i];
293 }
else if (dim != 0) {
294 for (
size_t i=0; i<coorddim; i++)
296 vol *= upper_[i] - lower_[i];
316 for (
size_t i=0; i<dim; i++)
327 for (
size_t i=0; i<coorddim; i++)
335 for (
size_t i=0; i<dim; i++)
346 for (
size_t i=0; i<coorddim; i++)
351 Dune::FieldVector<ctype,coorddim> lower_;
353 Dune::FieldVector<ctype,coorddim> upper_;
355 std::bitset<coorddim> axes_;
A unique label for each type of element that can occur in a grid.
constexpr GeometryType cube(unsigned int dim)
Returns a GeometryType representing a hypercube of dimension dim.
Definition type.hh:473
unspecified-type ReferenceElement
Returns the type of reference element for the argument types T...
Definition referenceelements.hh:417
Definition affinegeometry.hh:21
static const ReferenceElement & cube()
get hypercube reference elements
Definition referenceelements.hh:210
A geometry implementation for axis-aligned hypercubes.
Definition axisalignedcubegeometry.hh:50
Volume volume() const
Return the element volume.
Definition axisalignedcubegeometry.hh:286
AxisAlignedCubeGeometry(const Dune::FieldVector< ctype, coorddim > lower, const Dune::FieldVector< ctype, coorddim > upper, const std::bitset< coorddim > &axes)
Constructor from a lower left and an upper right corner.
Definition axisalignedcubegeometry.hh:133
AxisAlignedCubeGeometry(const Dune::FieldVector< ctype, coorddim > lower, const Dune::FieldVector< ctype, coorddim > upper)
Constructor from a lower left and an upper right corner.
Definition axisalignedcubegeometry.hh:115
JacobianInverse jacobianInverse(const LocalCoordinate &local) const
Inverse Jacobian of the transformation from local to global coordinates.
Definition axisalignedcubegeometry.hh:226
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
Inverse Jacobian transposed of the transformation from local to global coordinates.
Definition axisalignedcubegeometry.hh:208
std::conditional< dim==coorddim, DiagonalMatrix< ctype, dim >, FieldMatrix< ctype, dim, coorddim > >::type JacobianTransposed
Return type of jacobianTransposed.
Definition axisalignedcubegeometry.hh:81
AxisAlignedCubeGeometry(const Dune::FieldVector< ctype, coorddim > lower)
Constructor from a single point only.
Definition axisalignedcubegeometry.hh:150
static constexpr int mydimension
Dimension of the cube element.
Definition axisalignedcubegeometry.hh:56
static constexpr int coorddimension
Dimension of the world space that the cube element is embedded in.
Definition axisalignedcubegeometry.hh:59
GlobalCoordinate corner(int k) const
Return world coordinates of the k-th corner of the element.
Definition axisalignedcubegeometry.hh:260
ctype Volume
Type used for volume.
Definition axisalignedcubegeometry.hh:71
FieldVector< ctype, dim > LocalCoordinate
Type used for a vector of element coordinates.
Definition axisalignedcubegeometry.hh:65
friend Dune::Transitional::ReferenceElement< ctype, Dim< dim > > referenceElement(const AxisAlignedCubeGeometry &)
Definition axisalignedcubegeometry.hh:307
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
Jacobian transposed of the transformation from local to global coordinates.
Definition axisalignedcubegeometry.hh:196
FieldVector< ctype, coorddim > GlobalCoordinate
Type used for a vector of world coordinates.
Definition axisalignedcubegeometry.hh:68
LocalCoordinate local(const GlobalCoordinate &global) const
Map a point in global (world) coordinates to element coordinates.
Definition axisalignedcubegeometry.hh:180
CoordType ctype
Type used for single coordinate coefficients.
Definition axisalignedcubegeometry.hh:62
std::conditional_t< dim==coorddim, DiagonalMatrix< ctype, dim >, FieldMatrix< ctype, coorddim, dim > > Jacobian
Return type of jacobian.
Definition axisalignedcubegeometry.hh:100
GeometryType type() const
Type of the cube. Here: a hypercube of the correct dimension.
Definition axisalignedcubegeometry.hh:155
int corners() const
Return the number of corners of the element.
Definition axisalignedcubegeometry.hh:254
Jacobian jacobian(const LocalCoordinate &local) const
Jacobian of the transformation from local to global coordinates.
Definition axisalignedcubegeometry.hh:220
std::conditional< dim==coorddim, DiagonalMatrix< ctype, dim >, FieldMatrix< ctype, coorddim, dim > >::type JacobianInverseTransposed
Return type of jacobianInverseTransposed.
Definition axisalignedcubegeometry.hh:91
Volume integrationElement(const LocalCoordinate &local) const
Return the integration element, i.e., the determinant term in the integral transformation formula.
Definition axisalignedcubegeometry.hh:234
GlobalCoordinate center() const
Return center of mass of the element.
Definition axisalignedcubegeometry.hh:240
bool affine() const
Return if the element is affine. Here: yes.
Definition axisalignedcubegeometry.hh:302
std::conditional_t< dim==coorddim, DiagonalMatrix< ctype, dim >, FieldMatrix< ctype, dim, coorddim > > JacobianInverse
Return type of jacobianInverse.
Definition axisalignedcubegeometry.hh:109
GlobalCoordinate global(const LocalCoordinate &local) const
Map a point in local (element) coordinates to world coordinates.
Definition axisalignedcubegeometry.hh:161
Unique label for each type of entities that can occur in DUNE grids.
Definition type.hh:126