dune-grid  2.5-git
albertagrid/geometry.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_ALBERTA_GEOMETRY_HH
4 #define DUNE_ALBERTA_GEOMETRY_HH
5 
9 
10 #if HAVE_ALBERTA
11 
12 namespace Dune
13 {
14 
15  // Forward Declarations
16  // --------------------
17 
18  template< int dim, int dimworld >
19  class AlbertaGrid;
20 
21 
22 
23  // AlbertaGridCoordinateReader
24  // ---------------------------
25 
26  template< int codim, class GridImp >
28  {
29  typedef typename std::remove_const< GridImp >::type Grid;
30 
31  static const int dimension = Grid::dimension;
32  static const int codimension = codim;
33  static const int mydimension = dimension - codimension;
35 
37 
39  typedef FieldVector< ctype, coorddimension > Coordinate;
40 
41  AlbertaGridCoordinateReader ( const GridImp &grid,
42  const ElementInfo &elementInfo,
43  int subEntity )
44  : grid_( grid ),
45  elementInfo_( elementInfo ),
46  subEntity_( subEntity )
47  {}
48 
49  const ElementInfo &elementInfo () const
50  {
51  return elementInfo_;
52  }
53 
54  void coordinate ( int i, Coordinate &x ) const
55  {
56  assert( !elementInfo_ == false );
57  assert( (i >= 0) && (i <= mydimension) );
58 
59  const int k = mapVertices( subEntity_, i );
60  const Alberta::GlobalVector &coord = grid_.getCoord( elementInfo_, k );
61  for( int j = 0; j < coorddimension; ++j )
62  x[ j ] = coord[ j ];
63  }
64 
65  bool hasDeterminant () const
66  {
67  return false;
68  }
69 
70  ctype determinant () const
71  {
72  assert( hasDeterminant() );
73  return ctype( 0 );
74  }
75 
76  private:
77  static int mapVertices ( int subEntity, int i )
78  {
80  }
81 
82  const Grid &grid_;
83  const ElementInfo &elementInfo_;
84  const int subEntity_;
85  };
86 
87 
88 
89  // AlbertaGridGeometry
90  // -------------------
91 
104  template< int mydim, int cdim, class GridImp >
106  {
108 
109  // remember type of the grid
110  typedef GridImp Grid;
111 
112  // dimension of barycentric coordinates
113  static const int dimbary = mydim + 1;
114 
115  public:
118 
119  static const int dimension = Grid :: dimension;
120  static const int mydimension = mydim;
121  static const int codimension = dimension - mydimension;
122  static const int coorddimension = cdim;
123 
124  typedef FieldVector< ctype, mydimension > LocalCoordinate;
125  typedef FieldVector< ctype, coorddimension > GlobalCoordinate;
126 
127  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
128  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
129 
130  private:
131  static const int numCorners = mydimension + 1;
132 
133  typedef FieldMatrix< ctype, numCorners, coorddimension > CoordMatrix;
134 
135  public:
137  {
138  invalidate();
139  }
140 
141  template< class CoordReader >
142  AlbertaGridGeometry ( const CoordReader &coordReader )
143  {
144  build( coordReader );
145  }
146 
149  {
150  typedef typename Impl::SimplexTopology< mydimension >::type Topology;
151  return GeometryType( Topology() );
152  }
153 
155  bool affine () const { return true; }
156 
158  int corners () const
159  {
160  return numCorners;
161  }
162 
164  GlobalCoordinate corner ( const int i ) const
165  {
166  assert( (i >= 0) && (i < corners()) );
167  return coord_[ i ];
168  }
169 
171  GlobalCoordinate center () const
172  {
173  return centroid_;
174  }
175 
177  GlobalCoordinate global ( const LocalCoordinate &local ) const;
178 
180  LocalCoordinate local ( const GlobalCoordinate &global ) const;
181 
187  ctype integrationElement () const
188  {
189  assert( calcedDet_ );
190  return elDet_;
191  }
192 
194  ctype integrationElement ( const LocalCoordinate &local ) const
195  {
196  return integrationElement();
197  }
198 
200  ctype volume () const
201  {
202  return integrationElement() / ctype( Factorial< mydimension >::factorial );
203  }
204 
210  const JacobianTransposed &jacobianTransposed () const;
211 
213  const JacobianTransposed &
214  jacobianTransposed ( const LocalCoordinate &local ) const
215  {
216  return jacobianTransposed();
217  }
218 
224  const JacobianInverseTransposed &jacobianInverseTransposed () const;
225 
227  const JacobianInverseTransposed &
228  jacobianInverseTransposed ( const LocalCoordinate &local ) const
229  {
230  return jacobianInverseTransposed();
231  }
232 
233  //***********************************************************************
234  // Methods that not belong to the Interface, but have to be public
235  //***********************************************************************
236 
238  void invalidate ()
239  {
240  builtJT_ = false;
241  builtJTInv_ = false;
242  calcedDet_ = false;
243  }
244 
246  template< class CoordReader >
247  void build ( const CoordReader &coordReader );
248 
249  void print ( std::ostream &out ) const;
250 
251  private:
252  // calculates the volume of the element
253  ctype elDeterminant () const
254  {
255  return std::abs( Alberta::determinant( jacobianTransposed() ) );
256  }
257 
259  CoordMatrix coord_;
260 
262  GlobalCoordinate centroid_;
263 
264  // storage for the transposed of the jacobian
265  mutable JacobianTransposed jT_;
266 
267  // storage for the transposed inverse of the jacboian
268  mutable JacobianInverseTransposed jTInv_;
269 
270  // has jT_ been computed, yet?
271  mutable bool builtJT_;
272  // has jTInv_ been computed, yet?
273  mutable bool builtJTInv_;
274 
275  mutable bool calcedDet_;
276  mutable ctype elDet_;
277  };
278 
279 
280 
281  // AlbertaGridGlobalGeometry
282  // -------------------------
283 
284  template< int mydim, int cdim, class GridImp >
286  : public AlbertaGridGeometry< mydim, cdim, GridImp >
287  {
290 
291  public:
293  : Base()
294  {}
295 
296  template< class CoordReader >
297  AlbertaGridGlobalGeometry ( const CoordReader &coordReader )
298  : Base( coordReader )
299  {}
300  };
301 
302 
303 #if !DUNE_ALBERTA_CACHE_COORDINATES
304  template< int dim, int cdim >
305  class AlbertaGridGlobalGeometry< dim, cdim, const AlbertaGrid< dim, cdim > >
306  {
308 
309  // remember type of the grid
311 
312  // dimension of barycentric coordinates
313  static const int dimbary = dim + 1;
314 
316 
317  public:
320 
321  static const int dimension = Grid::dimension;
322  static const int mydimension = dimension;
323  static const int codimension = dimension - mydimension;
324  static const int coorddimension = cdim;
325 
326  typedef FieldVector< ctype, mydimension > LocalCoordinate;
327  typedef FieldVector< ctype, coorddimension > GlobalCoordinate;
328 
329  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
330  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
331 
332  private:
333  static const int numCorners = mydimension + 1;
334 
335  typedef FieldMatrix< ctype, numCorners, coorddimension > CoordMatrix;
336 
337  public:
339  {
340  invalidate();
341  }
342 
343  template< class CoordReader >
344  AlbertaGridGlobalGeometry ( const CoordReader &coordReader )
345  {
346  build( coordReader );
347  }
348 
351  {
352  typedef typename Impl::SimplexTopology< mydimension >::type Topology;
353  return GeometryType( Topology() );
354  }
355 
357  int corners () const
358  {
359  return numCorners;
360  }
361 
363  GlobalCoordinate corner ( const int i ) const
364  {
365  assert( (i >= 0) && (i < corners()) );
366  const Alberta::GlobalCoordinate &x = elementInfo_.coordinate( i );
367  GlobalCoordinate y;
368  for( int j = 0; j < coorddimension; ++j )
369  y[ j ] = x[ j ];
370  return y;
371  }
372 
374  GlobalCoordinate center () const
375  {
376  GlobalCoordinate centroid_ = corner( 0 );
377  for( int i = 1; i < numCorners; ++i )
378  centroid_ += corner( i );
379  centroid_ *= ctype( 1 ) / ctype( numCorners );
380  return centroid_;
381  }
382 
384  GlobalCoordinate global ( const LocalCoordinate &local ) const;
385 
387  LocalCoordinate local ( const GlobalCoordinate &global ) const;
388 
394  ctype integrationElement () const
395  {
396  return elementInfo_.geometryCache().integrationElement();
397  }
398 
400  ctype integrationElement ( const LocalCoordinate &local ) const
401  {
402  return integrationElement();
403  }
404 
406  ctype volume () const
407  {
408  return integrationElement() / ctype( Factorial< mydimension >::factorial );
409  }
410 
416  const JacobianTransposed &jacobianTransposed () const
417  {
418  return elementInfo_.geometryCache().jacobianTransposed();
419  }
420 
422  const JacobianTransposed &
423  jacobianTransposed ( const LocalCoordinate &local ) const
424  {
425  return jacobianTransposed();
426  }
427 
433  const JacobianInverseTransposed &jacobianInverseTransposed () const
434  {
435  return elementInfo_.geometryCache().jacobianInverseTransposed();
436  }
437 
439  const JacobianInverseTransposed &
440  jacobianInverseTransposed ( const LocalCoordinate &local ) const
441  {
442  return jacobianInverseTransposed();
443  }
444 
445  //***********************************************************************
446  // Methods that not belong to the Interface, but have to be public
447  //***********************************************************************
448 
450  void invalidate ()
451  {
452  elementInfo_ = ElementInfo();
453  }
454 
456  template< class CoordReader >
457  void build ( const CoordReader &coordReader )
458  {
459  elementInfo_ = coordReader.elementInfo();
460  }
461 
462  private:
463  ElementInfo elementInfo_;
464  };
465 #endif // #if !DUNE_ALBERTA_CACHE_COORDINATES
466 
467 
468 
469  // AlbertaGridLocalGeometryProvider
470  // --------------------------------
471 
472  template< class Grid >
474  {
476 
477  public:
478  typedef typename Grid::ctype ctype;
479 
480  static const int dimension = Grid::dimension;
481 
482  template< int codim >
483  struct Codim
484  {
485  typedef AlbertaGridGeometry< dimension-codim, dimension, Grid > LocalGeometry;
486  };
487 
490 
491  static const int numChildren = 2;
492  static const int numFaces = dimension + 1;
493 
494  static const int minFaceTwist = Alberta::Twist< dimension, dimension-1 >::minTwist;
495  static const int maxFaceTwist = Alberta::Twist< dimension, dimension-1 >::maxTwist;
496  static const int numFaceTwists = maxFaceTwist - minFaceTwist + 1;
497 
498  private:
499  struct GeoInFatherCoordReader;
500  struct FaceCoordReader;
501 
503  {
504  buildGeometryInFather();
505  buildFaceGeometry();
506  }
507 
509  {
510  for( int child = 0; child < numChildren; ++child )
511  {
512  delete geometryInFather_[ child ][ 0 ];
513  delete geometryInFather_[ child ][ 1 ];
514  }
515 
516  for( int i = 0; i < numFaces; ++i )
517  {
518  for( int j = 0; j < numFaceTwists; ++j )
519  delete faceGeometry_[ i ][ j ];
520  }
521  }
522 
523  void buildGeometryInFather();
524  void buildFaceGeometry();
525 
526  public:
527  const LocalElementGeometry &
528  geometryInFather ( int child, const int orientation = 1 ) const
529  {
530  assert( (child >= 0) && (child < numChildren) );
531  assert( (orientation == 1) || (orientation == -1) );
532  return *geometryInFather_[ child ][ (orientation + 1) / 2 ];
533  }
534 
535  const LocalFaceGeometry &
536  faceGeometry ( int face, int twist = 0 ) const
537  {
538  assert( (face >= 0) && (face < numFaces) );
539  assert( (twist >= minFaceTwist) && (twist <= maxFaceTwist) );
540  return *faceGeometry_[ face ][ twist - minFaceTwist ];
541  }
542 
543  static const This &instance ()
544  {
545  static This theInstance;
546  return theInstance;
547  }
548 
549  private:
550  template< int codim >
551  static int mapVertices ( int subEntity, int i )
552  {
554  }
555 
556  const LocalElementGeometry *geometryInFather_[ numChildren ][ 2 ];
557  const LocalFaceGeometry *faceGeometry_[ numFaces ][ numFaceTwists ];
558  };
559 
560 } // namespace Dune
561 
562 #endif // #if HAVE_ALBERTA
563 
564 #endif // #ifndef DUNE_ALBERTA_GEOMETRY_HH
Definition: misc.hh:440
int corners() const
number of corner the geometry
Definition: albertagrid/geometry.hh:158
provides a wrapper for ALBERTA&#39;s el_info structure
ct ctype
Define type used for coordinates in grid module.
Definition: common/grid.hh:522
ctype integrationElement() const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:187
static const int dimension
Definition: albertagrid/geometry.hh:31
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:228
GlobalCoordinate center() const
return center of geometry
Definition: albertagrid/geometry.hh:374
const JacobianInverseTransposed & jacobianInverseTransposed() const
transposed inverse of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:433
std::remove_const< GridImp >::type Grid
Definition: albertagrid/geometry.hh:29
ctype determinant() const
Definition: albertagrid/geometry.hh:70
AlbertaGridGlobalGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:344
void invalidate()
invalidate the geometry
Definition: albertagrid/geometry.hh:450
AlbertaGridGeometry()
Definition: albertagrid/geometry.hh:136
ctype integrationElement() const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:394
GeometryType
Type representing VTK&#39;s entity geometry types.
Definition: common.hh:178
FieldVector< ctype, coorddimension > GlobalCoordinate
Definition: albertagrid/geometry.hh:125
FieldVector< ctype, coorddimension > GlobalCoordinate
Definition: albertagrid/geometry.hh:327
static const int codimension
Definition: albertagrid/geometry.hh:32
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:423
Codim< 1 >::LocalGeometry LocalFaceGeometry
Definition: albertagrid/geometry.hh:489
FieldVector< ctype, coorddimension > Coordinate
Definition: albertagrid/geometry.hh:39
const LocalElementGeometry & geometryInFather(int child, const int orientation=1) const
Definition: albertagrid/geometry.hh:528
ctype volume() const
volume of geometry
Definition: albertagrid/geometry.hh:406
AlbertaGridGlobalGeometry()
Definition: albertagrid/geometry.hh:292
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: albertagrid/geometry.hh:164
ctype volume() const
volume of geometry
Definition: albertagrid/geometry.hh:200
Alberta::Real ctype
type of coordinates
Definition: albertagrid/geometry.hh:319
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Definition: albertagrid/geometry.hh:127
void abs(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:326
GeometryType type() const
obtain the type of reference element
Definition: albertagrid/geometry.hh:148
FieldVector< ctype, mydimension > LocalCoordinate
Definition: albertagrid/geometry.hh:326
Definition: albertagrid/geometry.hh:27
[ provides Dune::Grid ]
Definition: agrid.hh:137
FieldVector< ctype, mydimension > LocalCoordinate
Definition: albertagrid/geometry.hh:124
Grid::ctype ctype
Definition: albertagrid/geometry.hh:478
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:214
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Definition: albertagrid/geometry.hh:128
GlobalCoordinate center() const
return center of geometry
Definition: albertagrid/geometry.hh:171
AlbertaGridGeometry< dimension-codim, dimension, Grid > LocalGeometry
Definition: albertagrid/geometry.hh:485
Definition: albertagrid/geometry.hh:483
Codim< 0 >::LocalGeometry LocalElementGeometry
Definition: albertagrid/geometry.hh:488
Alberta::Real ctype
Definition: albertagrid/geometry.hh:36
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Definition: albertagrid/geometry.hh:329
void coordinate(int i, Coordinate &x) const
Definition: albertagrid/geometry.hh:54
Alberta::Real ctype
type of coordinates
Definition: albertagrid/geometry.hh:117
void build(const CoordReader &coordReader)
build the geometry from a coordinate reader
Definition: albertagrid/geometry.hh:457
const LocalFaceGeometry & faceGeometry(int face, int twist=0) const
Definition: albertagrid/geometry.hh:536
Definition: albertagrid/geometry.hh:285
bool hasDeterminant() const
Definition: albertagrid/geometry.hh:65
int corners() const
number of corner the geometry
Definition: albertagrid/geometry.hh:357
Include standard header files.
Definition: agrid.hh:59
static const int coorddimension
Definition: albertagrid/geometry.hh:34
void invalidate()
invalidate the geometry
Definition: albertagrid/geometry.hh:238
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: albertagrid/geometry.hh:363
AlbertaGridCoordinateReader(const GridImp &grid, const ElementInfo &elementInfo, int subEntity)
Definition: albertagrid/geometry.hh:41
const JacobianTransposed & jacobianTransposed() const
transposed of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:416
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Definition: albertagrid/geometry.hh:330
bool affine() const
returns always true since we only have simplices
Definition: albertagrid/geometry.hh:155
static K determinant(const FieldMatrix< K, 0, m > &matrix)
Definition: algebra.hh:28
Wrapper and interface classes for element geometries.
Definition: albertagrid/geometry.hh:473
geometry implementation for AlbertaGrid
Definition: albertagrid/geometry.hh:105
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:400
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:194
static const int mydimension
Definition: albertagrid/geometry.hh:33
GeometryType type() const
obtain the type of reference element
Definition: albertagrid/geometry.hh:350
The dimension of the grid.
Definition: common/grid.hh:387
AlbertaGridGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:142
ALBERTA REAL Real
Definition: misc.hh:45
ALBERTA REAL_D GlobalVector
Definition: misc.hh:47
Grid abstract base classThis class is the base class for all grid implementations. Although no virtual functions are used we call it abstract since its methods do not contain an implementation but forward to the methods of the derived class via the Barton-Nackman trick.
Definition: common/grid.hh:373
const ElementInfo & elementInfo() const
Definition: albertagrid/geometry.hh:49
AlbertaGridGlobalGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:297
Alberta::ElementInfo< dimension > ElementInfo
Definition: albertagrid/geometry.hh:38
Definition: misc.hh:528
static const This & instance()
Definition: albertagrid/geometry.hh:543
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:440
The dimension of the world the grid lives in.
Definition: common/grid.hh:393