5#ifndef DUNE_GEOMETRY_MAPPEDGEOMETRY_HH
6#define DUNE_GEOMETRY_MAPPEDGEOMETRY_HH
14#include <dune/common/copyableoptional.hh>
15#include <dune/common/exceptions.hh>
16#include <dune/common/fmatrix.hh>
17#include <dune/common/fvector.hh>
18#include <dune/common/math.hh>
19#include <dune/common/transpose.hh>
63template <
class Map,
class Geo>
71 using GlobalCoordinate = std::remove_reference_t<decltype(std::declval<Map>()(std::declval<typename Geo::GlobalCoordinate>()))>;
74 using ctype =
typename Geo::ctype;
83 using Volume = std::remove_reference_t<decltype(Dune::power(std::declval<ctype>(),
mydimension))>;
86 using Jacobian = FieldMatrix<ctype, coorddimension, mydimension>;
122 template <
class Geo_,
class Map_,
123 std::enable_if_t<Dune::IsInteroperable<Map, Map_>::value,
int> = 0,
124 std::enable_if_t<Dune::IsInteroperable<Geo, Geo_>::value,
int> = 0>
126 : mapping_(
std::forward<Map_>(mapping))
127 , dMapping_(derivative(*mapping_))
128 , geometry_(
std::forward<Geo_>(geometry))
145 return geometry_.type();
151 return geometry_.corners();
157 assert( (i >= 0) && (i <
corners()) );
158 return mapping()(geometry_.corner(i));
164 return mapping()(geometry_.center());
178 return mapping()(geometry_.global(
local));
200 Impl::GaussNewtonErrorCode err = Impl::gaussNewton(
206 if (err != Impl::GaussNewtonErrorCode::OK)
207 DUNE_THROW(Dune::Exception,
208 "Local coordinate can not be recovered from global coordinate, error code = " <<
int(err) <<
"\n"
209 <<
" (global(x) - y).two_norm() = " << (
global(x) - y).two_norm()
210 <<
" > tol = " << opts.absTol);
248 for (
int p = 2; p < opts.maxIt; ++p) {
250 if (abs(vol1 - vol0) < opts.absTol)
262 for (
const auto& qp : quadRule)
274 auto&& jLocal = geometry_.jacobian(
local);
275 auto&& jMapping = (*dMapping_)(geometry_.global(
local));
276 return jMapping * jLocal;
330 const Mapping& mapping ()
const
343 CopyableOptional<Mapping> mapping_;
346 CopyableOptional<DerivativeMapping> dMapping_;
356template <
class Map,
class Geo>
360template <
class Map,
class Geo>
A unique label for each type of element that can occur in a grid.
An implementation of the Geometry interface for affine geometries.
Definition affinegeometry.hh:21
Class providing access to the singletons of the reference elements.
Definition referenceelements.hh:128
typename Container::ReferenceElement ReferenceElement
The reference element type.
Definition referenceelements.hh:146
Geometry parametrized by a LocalFunction and a LocalGeometry.
Definition mappedgeometry.hh:65
std::remove_reference_t< decltype(std::declval< Map >()(std::declval< typename Geo::GlobalCoordinate >()))> GlobalCoordinate
type of global coordinates
Definition mappedgeometry.hh:71
GlobalCoordinate center() const
Map the center of the wrapped geometry.
Definition mappedgeometry.hh:162
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
Obtain the transposed of the Jacobian's inverse.
Definition mappedgeometry.hh:310
JacobianInverse jacobianInverse(const LocalCoordinate &local) const
Obtain the Jacobian's inverse.
Definition mappedgeometry.hh:296
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
Obtain the transposed of the Jacobian.
Definition mappedgeometry.hh:284
MappedGeometry(Map_ &&mapping, Geo_ &&geometry, bool affine=false)
Constructor from mapping to parametrize the geometry.
Definition mappedgeometry.hh:125
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
type of jacobian transposed
Definition mappedgeometry.hh:89
static constexpr int mydimension
geometry dimension
Definition mappedgeometry.hh:77
GeometryType type() const
Obtain the geometry type from the reference element.
Definition mappedgeometry.hh:143
typename Geo::LocalCoordinate LocalCoordinate
type of local coordinates
Definition mappedgeometry.hh:68
typename Geo::ctype ctype
coordinate type
Definition mappedgeometry.hh:74
LocalCoordinate local(const GlobalCoordinate &y, Impl::GaussNewtonOptions< ctype > opts={}) const
Evaluate the inverse coordinate mapping.
Definition mappedgeometry.hh:197
Geo Geometry
Definition mappedgeometry.hh:108
GlobalCoordinate corner(int i) const
Obtain coordinates of the i-th corner.
Definition mappedgeometry.hh:155
static constexpr int coorddimension
coordinate dimension
Definition mappedgeometry.hh:80
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
type of jacobian inverse transposed
Definition mappedgeometry.hh:95
Map Mapping
Definition mappedgeometry.hh:105
GlobalCoordinate global(const LocalCoordinate &local) const
Evaluate the coordinate mapping.
Definition mappedgeometry.hh:176
ReferenceElement refElement() const
Definition mappedgeometry.hh:323
Volume volume(const QuadratureRule< ctype, mydimension > &quadRule) const
Obtain the volume of the mapping's image by given quadrature rules.
Definition mappedgeometry.hh:259
int corners() const
Obtain number of corners of the corresponding reference element.
Definition mappedgeometry.hh:149
FieldMatrix< ctype, mydimension, coorddimension > JacobianInverse
type of jacobian inverse
Definition mappedgeometry.hh:92
friend ReferenceElement referenceElement(const MappedGeometry &geometry)
Obtain the reference-element related to this geometry.
Definition mappedgeometry.hh:316
Volume volume(Impl::ConvergenceOptions< ctype > opts={}) const
Obtain the volume of the mapping's image.
Definition mappedgeometry.hh:241
std::remove_reference_t< decltype(Dune::power(std::declval< ctype >(), mydimension))> Volume
type of volume
Definition mappedgeometry.hh:83
bool affine() const
Is this mapping affine? Not in general, since we don't know anything about the mapping....
Definition mappedgeometry.hh:137
std::remove_reference_t< decltype(derivative(std::declval< Map >()))> DerivativeMapping
Definition mappedgeometry.hh:111
Jacobian jacobian(const LocalCoordinate &local) const
Obtain the Jacobian.
Definition mappedgeometry.hh:272
Impl::FieldMatrixHelper< ctype > MatrixHelper
Definition mappedgeometry.hh:102
FieldMatrix< ctype, coorddimension, mydimension > Jacobian
type of jacobian
Definition mappedgeometry.hh:86
ctype integrationElement(const LocalCoordinate &local) const
Obtain the integration element.
Definition mappedgeometry.hh:225
Abstract base class for quadrature rules.
Definition quadraturerules.hh:214
A container for all quadrature rules of dimension dim
Definition quadraturerules.hh:260
Unique label for each type of entities that can occur in DUNE grids.
Definition type.hh:114