dune-grid 2.8.0
Loading...
Searching...
No Matches
identitygridgeometry.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_IDENTITYGRIDGEOMETRY_HH
4#define DUNE_IDENTITYGRIDGEOMETRY_HH
5
10#include <dune/common/fmatrix.hh>
11#include <dune/common/typetraits.hh>
13
14namespace Dune {
15
16 template<int mydim, int coorddim, class GridImp>
18 public GeometryDefaultImplementation <mydim, coorddim, GridImp, IdentityGridGeometry>
19 {
20 private:
21
22 typedef typename GridImp::ctype ctype;
23
24
25 public:
26
27 // The codimension of this entitypointer wrt the host grid
28 enum {CodimInHostGrid = GridImp::HostGridType::dimension - mydim};
29 enum {DimensionWorld = GridImp::HostGridType::dimensionworld};
30
31 // select appropriate hostgrid geometry via typeswitch
32 typedef typename GridImp::HostGridType::Traits::template Codim<CodimInHostGrid>::Geometry HostGridGeometryType;
33 typedef typename GridImp::HostGridType::Traits::template Codim<CodimInHostGrid>::Geometry HostGridLocalGeometryType;
34
35 typedef typename std::conditional<coorddim==DimensionWorld, HostGridGeometryType, HostGridLocalGeometryType>::type HostGridGeometry;
36
38 typedef typename HostGridGeometryType::JacobianInverseTransposed JacobianInverseTransposed;
39 typedef typename HostGridGeometryType::JacobianTransposed JacobianTransposed;
40
41
45 : hostGeometry_(hostGeometry)
46 {}
47
48
51 GeometryType type () const {
52 return hostGeometry_.type();
53 }
54
55 // return wether we have an affine mapping
56 bool affine() const {
57 return hostGeometry_.affine();
58 }
59
61 int corners () const {
62 return hostGeometry_.corners();
63 }
64
65
67 const FieldVector<ctype, coorddim> corner (int i) const {
68 return hostGeometry_.corner(i);
69 }
70
71
74 FieldVector<ctype, coorddim> global (const FieldVector<ctype, mydim>& local) const {
75 return hostGeometry_.global(local);
76 }
77
81 jacobianTransposed ( const FieldVector<ctype, mydim>& local ) const {
82 return hostGeometry_.jacobianTransposed(local);
83 }
84
87 FieldVector<ctype, mydim> local (const FieldVector<ctype, coorddim>& global) const {
88 return hostGeometry_.local(global);
89 }
90
91
93 bool checkInside(const FieldVector<ctype, mydim> &local) const {
94 return hostGeometry_.checkInside(local);
95 }
96
97
100 ctype integrationElement (const FieldVector<ctype, mydim>& local) const {
101 return hostGeometry_.integrationElement(local);
102 }
103
104
106 JacobianInverseTransposed jacobianInverseTransposed (const FieldVector<ctype, mydim>& local) const {
107 return hostGeometry_.jacobianInverseTransposed(local);
108 }
109
110
112
113 };
114
115} // namespace Dune
116
117#endif
Include standard header files.
Definition: agrid.hh:58
Default implementation for class Geometry.
Definition: common/geometry.hh:301
Definition: identitygridgeometry.hh:19
bool checkInside(const FieldVector< ctype, mydim > &local) const
Returns true if the point is in the current element.
Definition: identitygridgeometry.hh:93
JacobianInverseTransposed jacobianInverseTransposed(const FieldVector< ctype, mydim > &local) const
The Jacobian matrix of the mapping from the reference element to this element.
Definition: identitygridgeometry.hh:106
JacobianTransposed jacobianTransposed(const FieldVector< ctype, mydim > &local) const
Return the transposed of the Jacobian.
Definition: identitygridgeometry.hh:81
std::conditional< coorddim==DimensionWorld, HostGridGeometryType, HostGridLocalGeometryType >::type HostGridGeometry
Definition: identitygridgeometry.hh:35
bool affine() const
Definition: identitygridgeometry.hh:56
FieldVector< ctype, coorddim > global(const FieldVector< ctype, mydim > &local) const
Maps a local coordinate within reference element to global coordinate in element
Definition: identitygridgeometry.hh:74
IdentityGridGeometry(const HostGridGeometry &hostGeometry)
Definition: identitygridgeometry.hh:44
GeometryType type() const
Return the element type identifier.
Definition: identitygridgeometry.hh:51
@ DimensionWorld
Definition: identitygridgeometry.hh:29
const FieldVector< ctype, coorddim > corner(int i) const
access to coordinates of corners. Index is the number of the corner
Definition: identitygridgeometry.hh:67
HostGridGeometryType::JacobianTransposed JacobianTransposed
Definition: identitygridgeometry.hh:39
ctype integrationElement(const FieldVector< ctype, mydim > &local) const
Definition: identitygridgeometry.hh:100
HostGridGeometry hostGeometry_
Definition: identitygridgeometry.hh:111
GridImp::HostGridType::Traits::template Codim< CodimInHostGrid >::Geometry HostGridLocalGeometryType
Definition: identitygridgeometry.hh:33
GridImp::HostGridType::Traits::template Codim< CodimInHostGrid >::Geometry HostGridGeometryType
Definition: identitygridgeometry.hh:32
@ CodimInHostGrid
Definition: identitygridgeometry.hh:28
HostGridGeometryType::JacobianInverseTransposed JacobianInverseTransposed
type of jacobian transposed
Definition: identitygridgeometry.hh:38
int corners() const
return the number of corners of this element. Corners are numbered 0...n-1
Definition: identitygridgeometry.hh:61
Wrapper and interface classes for element geometries.