dune-grid  2.4
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 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 GenericGeometry::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  {
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
310  typedef AlbertaGrid< dim, cdim > 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 GenericGeometry::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
remove_const< GridImp >::type Grid
Definition: albertagrid/geometry.hh:29
GlobalCoordinate global(const LocalCoordinate &local) const
map a point from the refence element to the geometry
Definition: geometry.cc:32
ALBERTA REAL_D GlobalVector
Definition: misc.hh:47
geometry implementation for AlbertaGrid
Definition: albertagrid/geometry.hh:105
void coordinate(int i, Coordinate &x) const
Definition: albertagrid/geometry.hh:54
const LocalElementGeometry & geometryInFather(int child, const int orientation=1) const
Definition: albertagrid/geometry.hh:528
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: albertagrid/geometry.hh:164
FieldVector< ctype, mydimension > LocalCoordinate
Definition: albertagrid/geometry.hh:124
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:194
AlbertaGridCoordinateReader(const GridImp &grid, const ElementInfo &elementInfo, int subEntity)
Definition: albertagrid/geometry.hh:41
void build(const CoordReader &coordReader)
build the geometry from a coordinate reader
Definition: albertagrid/geometry.hh:457
Alberta::Real ctype
Definition: albertagrid/geometry.hh:36
void invalidate()
invalidate the geometry
Definition: albertagrid/geometry.hh:238
static const int dimension
Definition: albertagrid/geometry.hh:119
GlobalCoordinate center() const
return center of geometry
Definition: albertagrid/geometry.hh:171
static const int minFaceTwist
Definition: albertagrid/geometry.hh:494
static const int maxFaceTwist
Definition: albertagrid/geometry.hh:495
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
ctype integrationElement() const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:187
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Definition: albertagrid/geometry.hh:330
int corners() const
number of corner the geometry
Definition: albertagrid/geometry.hh:158
Definition: albertagrid/geometry.hh:285
Wrapper and interface classes for element geometries.
static const int numFaces
Definition: albertagrid/geometry.hh:492
Definition: misc.hh:440
static const int codimension
Definition: albertagrid/geometry.hh:121
Include standard header files.
Definition: agrid.hh:59
bool hasDeterminant() const
Definition: albertagrid/geometry.hh:65
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping's Jacobian
Definition: albertagrid/geometry.hh:440
GlobalCoordinate center() const
return center of geometry
Definition: albertagrid/geometry.hh:374
void print(std::ostream &out) const
Definition: geometry.cc:18
const JacobianTransposed & jacobianTransposed() const
transposed of the geometry mapping's Jacobian
Definition: geometry.cc:53
static const int dimension
Definition: albertagrid/geometry.hh:31
bool affine() const
returns always true since we only have simplices
Definition: albertagrid/geometry.hh:155
static const This & instance()
Definition: albertagrid/geometry.hh:543
static const int numChildren
Definition: albertagrid/geometry.hh:491
static const int dimension
Definition: albertagrid/geometry.hh:480
Alberta::ElementInfo< dimension > ElementInfo
Definition: albertagrid/geometry.hh:38
static const int mydimension
Definition: albertagrid/geometry.hh:33
void invalidate()
invalidate the geometry
Definition: albertagrid/geometry.hh:450
static const int mydimension
Definition: albertagrid/geometry.hh:120
const LocalFaceGeometry & faceGeometry(int face, int twist=0) const
Definition: albertagrid/geometry.hh:536
ctype integrationElement() const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:394
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Definition: albertagrid/geometry.hh:127
AlbertaGridGlobalGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:344
Definition: misc.hh:528
The dimension of the world the grid lives in.
Definition: common/grid.hh:408
FieldVector< ctype, coorddimension > GlobalCoordinate
Definition: albertagrid/geometry.hh:327
Definition: albertagrid/geometry.hh:27
static const int codimension
Definition: albertagrid/geometry.hh:32
GeometryType type() const
obtain the type of reference element
Definition: albertagrid/geometry.hh:350
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: albertagrid/geometry.hh:363
ctype volume() const
volume of geometry
Definition: albertagrid/geometry.hh:200
LocalCoordinate local(const GlobalCoordinate &global) const
map a point from the geometry to the reference element
Definition: geometry.cc:43
ct ctype
Define type used for coordinates in grid module.
Definition: common/grid.hh:548
AlbertaGridGeometry< dimension-codim, dimension, Grid > LocalGeometry
Definition: albertagrid/geometry.hh:485
ctype determinant() const
Definition: albertagrid/geometry.hh:70
The dimension of the grid.
Definition: common/grid.hh:402
static K determinant(const FieldMatrix< K, 0, m > &matrix)
Definition: algebra.hh:28
AlbertaGridGlobalGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:297
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping's Jacobian
Definition: albertagrid/geometry.hh:228
const ElementInfo & elementInfo() const
Definition: albertagrid/geometry.hh:49
FieldVector< ctype, mydimension > LocalCoordinate
Definition: albertagrid/geometry.hh:326
int corners() const
number of corner the geometry
Definition: albertagrid/geometry.hh:357
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Definition: albertagrid/geometry.hh:329
Definition: albertagrid/geometry.hh:473
Codim< 1 >::LocalGeometry LocalFaceGeometry
Definition: albertagrid/geometry.hh:489
ctype volume() const
volume of geometry
Definition: albertagrid/geometry.hh:406
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:400
[ provides Dune::Grid ]
Definition: agrid.hh:137
AlbertaGridGlobalGeometry()
Definition: albertagrid/geometry.hh:292
static const int coorddimension
Definition: albertagrid/geometry.hh:34
const JacobianTransposed & jacobianTransposed() const
transposed of the geometry mapping's Jacobian
Definition: albertagrid/geometry.hh:416
static const int numFaceTwists
Definition: albertagrid/geometry.hh:496
void abs(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:326
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping's Jacobian
Definition: albertagrid/geometry.hh:214
GeometryType type() const
obtain the type of reference element
Definition: albertagrid/geometry.hh:148
FieldVector< ctype, coorddimension > GlobalCoordinate
Definition: albertagrid/geometry.hh:125
Alberta::Real ctype
type of coordinates
Definition: albertagrid/geometry.hh:117
const JacobianInverseTransposed & jacobianInverseTransposed() const
transposed inverse of the geometry mapping's Jacobian
Definition: geometry.cc:71
Definition: albertagrid/geometry.hh:483
provides a wrapper for ALBERTA's el_info structure
Codim< 0 >::LocalGeometry LocalElementGeometry
Definition: albertagrid/geometry.hh:488
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping's Jacobian
Definition: albertagrid/geometry.hh:423
Grid::ctype ctype
Definition: albertagrid/geometry.hh:478
const JacobianInverseTransposed & jacobianInverseTransposed() const
transposed inverse of the geometry mapping's Jacobian
Definition: albertagrid/geometry.hh:433
AlbertaGridGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:142
static const int coorddimension
Definition: albertagrid/geometry.hh:122
ALBERTA REAL Real
Definition: misc.hh:45
void build(const CoordReader &coordReader)
build the geometry from a coordinate reader
Definition: geometry.cc:88
Alberta::Real ctype
type of coordinates
Definition: albertagrid/geometry.hh:319
FieldVector< ctype, coorddimension > Coordinate
Definition: albertagrid/geometry.hh:39
AlbertaGridGeometry()
Definition: albertagrid/geometry.hh:136
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Definition: albertagrid/geometry.hh:128