3 #ifndef DUNE_POLYHEDRALGRID_GEOMETRY_HH
4 #define DUNE_POLYHEDRALGRID_GEOMETRY_HH
8 #include <dune/common/version.hh>
9 #include <dune/common/fmatrix.hh>
10 #include <dune/grid/common/geometry.hh>
12 #include <dune/geometry/referenceelements.hh>
13 #include <dune/geometry/type.hh>
14 #include <dune/geometry/multilineargeometry.hh>
23 template<
int,
int,
class >
class PolyhedralGridGeometry;
24 template<
int,
int,
class >
class PolyhedralGridLocalGeometry;
29 template<
int mydim,
int cdim,
class Gr
id >
32 static const int dimension = Grid::dimension;
33 static const int mydimension = mydim;
34 static const int codimension = dimension - mydimension;
36 static const int dimensionworld = Grid::dimensionworld;
37 static const int coorddimension = dimensionworld;
39 typedef typename Grid::ctype ctype;
40 typedef Dune::FieldVector< ctype, coorddimension > GlobalCoordinate;
41 typedef Dune::FieldVector< ctype, mydimension > LocalCoordinate;
43 typedef typename Grid::Traits::ExtraData ExtraData;
44 typedef typename Grid::Traits::template Codim<codimension>::EntitySeed EntitySeed;
48 :
public Dune::MultiLinearGeometryTraits< ct >
53 :
public Dune::ForwardIteratorFacade< Iterator, GlobalCoordinate, GlobalCoordinate >
57 explicit Iterator(
const Storage* ptr,
int count ) : data_( ptr ), count_( count ) {}
59 GlobalCoordinate dereference()
const {
return data_->corner( count_ ); }
60 void increment() { ++count_; }
62 bool equals(
const Iterator& other )
const {
return count_ == other.count_; }
69 Storage( ExtraData data, EntitySeed seed )
70 : data_( data ), seed_( seed )
74 : data_( data ), seed_()
77 ExtraData data()
const {
return data_; }
78 bool isValid ()
const {
return seed_.isValid(); }
80 GlobalCoordinate operator [] (
const int i)
const {
return corner( i ); }
82 Iterator begin()
const {
return Iterator(
this, 0); }
83 Iterator end ()
const {
return Iterator(
this, corners()); }
85 int corners ()
const {
return data()->corners( seed_ ); }
86 GlobalCoordinate corner (
const int i )
const {
return data()->corner( seed_, i ); }
87 GlobalCoordinate center ()
const {
return data()->centroids( seed_ ); }
89 ctype volume()
const {
return data()->volumes( seed_ ); }
91 const EntitySeed& seed ()
const {
return seed_; }
94 template <
int mdim,
int cordim>
101 typedef Dune::MultiLinearGeometry< ctype, mydimension, coorddimension, PolyhedralMultiLinearGeometryTraits<ctype> >
102 MultiLinearGeometryType;
105 CornerStorage<mydimension, coorddimension >::Type CornerStorageType;
114 typedef Dune::Impl::FieldMatrixHelper< ctype > MatrixHelperType;
121 : storage_( data, seed )
123 GeometryType myType = type();
124 if( ! myType.isNone() && storage_.isValid() )
126 geometryImpl_.reset(
new MultiLinearGeometryType(myType, storage_) );
131 GeometryType type ()
const {
return data()->geometryType( storage_.seed() ); }
132 bool affine ()
const {
return (geometryImpl_) ? geometryImpl_->affine() :
false; }
134 int corners ()
const {
return storage_.corners(); }
135 GlobalCoordinate corner (
const int i )
const {
return storage_.corner( i ); }
136 GlobalCoordinate center ()
const
138 if( type().isNone() )
140 return storage_.center();
144 #if DUNE_VERSION_NEWER(DUNE_GRID,2,7)
145 const auto refElem = Dune::referenceElement< ctype, mydim > ( type() );
147 const auto& refElem = Dune::ReferenceElements< ctype, mydim >::general( type() );
149 return global( refElem.position(0,0) );
153 GlobalCoordinate global(
const LocalCoordinate&
local)
const
157 return geometryImpl_->global(
local );
165 LocalCoordinate
local(
const GlobalCoordinate& global)
const
169 return geometryImpl_->local( global );
173 return LocalCoordinate( 1 );
176 ctype integrationElement (
const LocalCoordinate &
local )
const
180 return geometryImpl_->integrationElement(
local );
185 ctype volume ()
const
189 return geometryImpl_->volume();
191 return storage_.volume();
198 return geometryImpl_->jacobianTransposed(
local );
201 DUNE_THROW(NotImplemented,
"jacobianTransposed not implemented");
209 return geometryImpl_->jacobianInverseTransposed(
local );
212 DUNE_THROW(NotImplemented,
"jacobianInverseTransposed not implemented");
216 ExtraData data()
const {
return storage_.data(); }
219 CornerStorageType storage_;
220 std::shared_ptr< MultiLinearGeometryType > geometryImpl_;
227 template<
int mydim,
int cdim,
class Gr
id >
234 typedef typename Base::ExtraData ExtraData;
235 typedef typename Base::EntitySeed EntitySeed;
246 template<
int mydim,
int cdim,
class Gr
id >
253 typedef typename Base::ExtraData ExtraData;
Definition: geometry.hh:230
Definition: geometry.hh:249
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
Definition: geometry.hh:96
Definition: geometry.hh:54
Definition: geometry.hh:51
Definition: geometry.hh:49
Definition: geometry.hh:31
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: geometry.hh:108
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: geometry.hh:111
LocalCoordinate local(const GlobalCoordinate &global) const
Mapping from the cell to the reference domain.
Definition: geometry.hh:165