4 #ifndef DUNE_GEOMETRY_AXISALIGNED_CUBE_GEOMETRY_HH
5 #define DUNE_GEOMETRY_AXISALIGNED_CUBE_GEOMETRY_HH
13 #include <dune/common/fvector.hh>
14 #include <dune/common/fmatrix.hh>
15 #include <dune/common/diagonalmatrix.hh>
16 #include <dune/common/unused.hh>
47 template <
class CoordType,
unsigned int dim,
unsigned int coorddim>
78 typedef typename std::conditional<dim==coorddim,
79 DiagonalMatrix<ctype,dim>,
88 typedef typename std::conditional<dim==coorddim,
89 DiagonalMatrix<ctype,dim>,
97 const Dune::FieldVector<ctype,coorddim> upper)
103 axes_ = (1<<coorddim)-1;
114 const Dune::FieldVector<ctype,coorddim> upper,
115 const std::bitset<coorddim>& axes)
120 assert(axes.count()==dim);
121 for (
size_t i=0; i<coorddim; i++)
123 upper_[i] = lower_[i];
137 lower_ = other.lower_;
138 upper_ = other.upper_;
153 if (dim == coorddim) {
154 for (
size_t i=0; i<coorddim; i++)
155 result[i] = lower_[i] +
local[i]*(upper_[i] - lower_[i]);
156 }
else if (dim == 0) {
160 for (
size_t i=0; i<coorddim; i++)
161 result[i] = (axes_[i])
162 ? lower_[i] +
local[lc++]*(upper_[i] - lower_[i])
172 if (dim == coorddim) {
173 for (
size_t i=0; i<dim; i++)
174 result[i] = (
global[i] - lower_[i]) / (upper_[i] - lower_[i]);
175 }
else if (dim != 0) {
177 for (
size_t i=0; i<coorddim; i++)
179 result[lc++] = (
global[i] - lower_[i]) / (upper_[i] - lower_[i]);
224 for (
size_t i=0; i<coorddim; i++)
225 result[i] = CoordType(0.5) * (lower_[i] + upper_[i]);
240 if (dim == coorddim) {
241 for (
size_t i=0; i<coorddim; i++)
242 result[i] = (k & (1<<i)) ? upper_[i] : lower_[i];
243 }
else if (dim == 0) {
246 unsigned int mask = 1;
248 for (
size_t i=0; i<coorddim; i++) {
250 result[i] = lower_[i];
252 result[i] = (k & mask) ? upper_[i] : lower_[i];
266 if (dim == coorddim) {
267 for (
size_t i=0; i<dim; i++)
268 vol *= upper_[i] - lower_[i];
270 }
else if (dim != 0) {
271 for (
size_t i=0; i<coorddim; i++)
273 vol *= upper_[i] - lower_[i];
293 for (
size_t i=0; i<dim; i++)
304 for (
size_t i=0; i<coorddim; i++)
312 for (
size_t i=0; i<dim; i++)
323 for (
size_t i=0; i<coorddim; i++)
328 Dune::FieldVector<ctype,coorddim> lower_;
330 Dune::FieldVector<ctype,coorddim> upper_;
332 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:775
unspecified-type ReferenceElement
Returns the type of reference element for the argument types T...
Definition: referenceelements.hh:415
Definition: affinegeometry.hh:19
static const ReferenceElement & cube()
get hypercube reference elements
Definition: referenceelements.hh:208
A geometry implementation for axis-aligned hypercubes.
Definition: axisalignedcubegeometry.hh:49
Volume volume() const
Return the element volume.
Definition: axisalignedcubegeometry.hh:263
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:113
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:96
AxisAlignedCubeGeometry(const Dune::FieldVector< ctype, coorddim > lower)
Constructor from a single point only.
Definition: axisalignedcubegeometry.hh:130
GlobalCoordinate corner(int k) const
Return world coordinates of the k-th corner of the element.
Definition: axisalignedcubegeometry.hh:237
JacobianTransposed jacobianTransposed(DUNE_UNUSED const LocalCoordinate &local) const
Jacobian transposed of the transformation from local to global coordinates.
Definition: axisalignedcubegeometry.hh:185
ctype Volume
Type used for volume.
Definition: axisalignedcubegeometry.hh:70
@ mydimension
Definition: axisalignedcubegeometry.hh:55
FieldVector< ctype, dim > LocalCoordinate
Type used for a vector of element coordinates.
Definition: axisalignedcubegeometry.hh:64
AxisAlignedCubeGeometry & operator=(const AxisAlignedCubeGeometry &other)
Assignment operator.
Definition: axisalignedcubegeometry.hh:135
friend Dune::Transitional::ReferenceElement< ctype, Dim< dim > > referenceElement(const AxisAlignedCubeGeometry &geometry)
Definition: axisalignedcubegeometry.hh:284
FieldVector< ctype, coorddim > GlobalCoordinate
Type used for a vector of world coordinates.
Definition: axisalignedcubegeometry.hh:67
JacobianInverseTransposed jacobianInverseTransposed(DUNE_UNUSED const LocalCoordinate &local) const
Jacobian transposed of the transformation from local to global coordinates.
Definition: axisalignedcubegeometry.hh:197
LocalCoordinate local(const GlobalCoordinate &global) const
Map a point in global (world) coordinates to element coordinates.
Definition: axisalignedcubegeometry.hh:169
CoordType ctype
Type used for single coordinate coefficients.
Definition: axisalignedcubegeometry.hh:61
@ coorddimension
Definition: axisalignedcubegeometry.hh:58
std::conditional< dim==coorddim, DiagonalMatrix< ctype, dim >, FieldMatrix< ctype, dim, coorddim > >::type JacobianTransposed
Return type of jacobianTransposed.
Definition: axisalignedcubegeometry.hh:80
GeometryType type() const
Type of the cube. Here: a hypercube of the correct dimension.
Definition: axisalignedcubegeometry.hh:144
int corners() const
Return the number of corners of the element.
Definition: axisalignedcubegeometry.hh:231
ctype integrationElement(DUNE_UNUSED const LocalCoordinate &local) const
Return the integration element, i.e., the determinant term in the integral transformation formula.
Definition: axisalignedcubegeometry.hh:211
std::conditional< dim==coorddim, DiagonalMatrix< ctype, dim >, FieldMatrix< ctype, coorddim, dim > >::type JacobianInverseTransposed
Return type of jacobianInverseTransposed.
Definition: axisalignedcubegeometry.hh:90
GlobalCoordinate center() const
Return center of mass of the element.
Definition: axisalignedcubegeometry.hh:217
bool affine() const
Return if the element is affine. Here: yes.
Definition: axisalignedcubegeometry.hh:279
GlobalCoordinate global(const LocalCoordinate &local) const
Map a point in local (element) coordinates to world coordinates.
Definition: axisalignedcubegeometry.hh:150
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:279