3 #ifndef DUNE_POLYHEDRALGRID_GEOMETRY_HH 4 #define DUNE_POLYHEDRALGRID_GEOMETRY_HH 6 #include <dune/common/version.hh> 7 #include <dune/common/fmatrix.hh> 8 #include <dune/grid/common/geometry.hh> 10 #if DUNE_VERSION_NEWER(DUNE_GEOMETRY, 2, 5 ) 11 #include <dune/geometry/type.hh> 13 #include <dune/geometry/genericgeometry/geometrytraits.hh> 14 #include <dune/geometry/genericgeometry/matrixhelper.hh> 30 template<
int mydim,
int cdim,
class Gr
id >
33 static const int dimension = Grid::dimension;
34 static const int mydimension = mydim;
35 static const int codimension = dimension - mydimension;
37 static const int dimensionworld = Grid::dimensionworld;
38 static const int coorddimension = dimensionworld;
40 typedef typename Grid::ctype ctype;
41 typedef Dune::FieldVector< ctype, coorddimension > GlobalCoordinate;
42 typedef Dune::FieldVector< ctype, mydimension > LocalCoordinate;
51 #if DUNE_VERSION_NEWER(DUNE_GRID,2,5) 52 typedef Dune::Impl::FieldMatrixHelper< double > MatrixHelperType;
54 typedef Dune::GenericGeometry::MatrixHelper< Dune::GenericGeometry::DuneCoordTraits<double> > MatrixHelperType;
57 typedef typename Grid::Traits::ExtraData ExtraData;
58 typedef typename Grid::Traits::template Codim<codimension>::EntitySeed EntitySeed;
70 GeometryType type ()
const {
return data()->geomTypes(codimension)[0]; }
71 bool affine ()
const {
return false; }
73 int corners ()
const {
return data()->corners( seed_ ); }
74 GlobalCoordinate corner (
const int i )
const {
return data()->corner( seed_, i ); }
75 GlobalCoordinate center ()
const {
return data()->centroids( seed_ ); }
77 GlobalCoordinate global(
const LocalCoordinate&
local)
const 79 const GeometryType geomType = type();
80 if( geomType.isCube() )
82 assert( mydimension == 3 );
83 assert( coorddimension == 3 );
86 LocalCoordinate uvw[2] = { LocalCoordinate(1.0),
local };
93 const int pat[8][3] = { { 0, 0, 0 },
102 const int nCorners = corners();
105 GlobalCoordinate xyz(0.0);
106 for (
int i = 0; i < nCorners ; ++i)
108 GlobalCoordinate cornerContrib = corner(i);
111 for (
int j = 0; j < mydimension; ++j)
114 factor *= uvw[ pat[ i ][ j ] ][ j ];
116 cornerContrib *= factor;
117 xyz += cornerContrib;
121 else if ( geomType.isNone() )
128 DUNE_THROW(NotImplemented,
"global for geometry type " << geomType <<
" is not supported!");
129 return GlobalCoordinate( 0 );
135 LocalCoordinate
local(
const GlobalCoordinate& y)
const 137 const GeometryType geomType = type();
138 if( geomType.isCube() )
142 const ctype epsilon = 1e-12;
143 #if DUNE_VERSION_NEWER(DUNE_GRID,2,3) 144 const ReferenceElement< ctype , mydimension > & refElement =
145 ReferenceElements< ctype, mydimension >::general(type());
147 const GenericReferenceElement< ctype , mydimension > & refElement =
148 GenericReferenceElements< ctype, mydimension >::general(type());
151 LocalCoordinate x = refElement.position(0,0);
156 GlobalCoordinate z = global(x);
158 MatrixHelperType::template xTRightInvA<mydimension,coorddimension>(JT, z, dx );
160 }
while (dx.two_norm2() > epsilon*epsilon);
163 else if ( geomType.isNone() )
166 return LocalCoordinate( 1 );
170 DUNE_THROW(NotImplemented,
"local for geometry type " << geomType <<
" is not supported!");
171 return LocalCoordinate( 0 );
175 ctype integrationElement (
const LocalCoordinate & )
const {
return volume(); }
176 ctype volume ()
const {
return data()->volumes( seed_ ); }
178 #if DUNE_VERSION_NEWER(DUNE_GRID,2,4) 181 DUNE_THROW(NotImplemented,
"jacobianTransposed not implemented");
187 DUNE_THROW(NotImplemented,
"jacobianInverseTransposed not implemented");
193 DUNE_THROW(NotImplemented,
"jacobianTransposed not implemented");
200 DUNE_THROW(NotImplemented,
"jacobianInverseTransposed not implemented");
206 ExtraData data()
const {
return data_; }
218 template<
int mydim,
int cdim,
class Gr
id >
219 class PolyhedralGridGeometry
220 :
public PolyhedralGridBasicGeometry< mydim, cdim, Grid >
222 typedef PolyhedralGridBasicGeometry< mydim, cdim, Grid > Base;
225 typedef typename Base::ExtraData ExtraData;
226 typedef typename Base::EntitySeed EntitySeed;
228 explicit PolyhedralGridGeometry ( ExtraData data )
232 PolyhedralGridGeometry ( ExtraData data,
const EntitySeed& seed )
237 template<
int mydim,
int cdim,
class Gr
id >
238 class PolyhedralGridLocalGeometry
239 :
public PolyhedralGridBasicGeometry< mydim, cdim, Grid >
241 typedef PolyhedralGridBasicGeometry< mydim, cdim, Grid > Base;
244 typedef typename Base::ExtraData ExtraData;
246 explicit PolyhedralGridLocalGeometry ( ExtraData data )
252 #if ! DUNE_VERSION_NEWER(DUNE_GRID,2,4) 253 namespace FacadeOptions
258 template<
int mydim,
int cdim,
class Gr
idImp >
262 static const bool v =
false;
265 template<
int mydim,
int cdim,
class Gr
idImp >
269 static const bool v =
false;
278 #endif // #ifndef DUNE_POLYHEDRALGRID_GEOMETRY_HH Definition: geometry.hh:24
LocalCoordinate local(const GlobalCoordinate &y) const
Mapping from the cell to the reference domain.
Definition: geometry.hh:135
Holds the implementation of the CpGrid as a pimple.
Definition: OpmParserIncludes.hpp:42
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: geometry.hh:45
Definition: geometry.hh:31
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: geometry.hh:48
Definition: geometry.hh:25