dune-grid  2.6-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  return GeometryTypes::simplex( mydimension );
151  }
152 
154  bool affine () const { return true; }
155 
157  int corners () const
158  {
159  return numCorners;
160  }
161 
163  GlobalCoordinate corner ( const int i ) const
164  {
165  assert( (i >= 0) && (i < corners()) );
166  return coord_[ i ];
167  }
168 
170  GlobalCoordinate center () const
171  {
172  return centroid_;
173  }
174 
176  GlobalCoordinate global ( const LocalCoordinate &local ) const;
177 
179  LocalCoordinate local ( const GlobalCoordinate &global ) const;
180 
186  ctype integrationElement () const
187  {
188  assert( calcedDet_ );
189  return elDet_;
190  }
191 
193  ctype integrationElement ( const LocalCoordinate &local ) const
194  {
195  return integrationElement();
196  }
197 
199  ctype volume () const
200  {
201  return integrationElement() / ctype( Factorial< mydimension >::factorial );
202  }
203 
209  const JacobianTransposed &jacobianTransposed () const;
210 
212  const JacobianTransposed &
213  jacobianTransposed ( const LocalCoordinate &local ) const
214  {
215  return jacobianTransposed();
216  }
217 
223  const JacobianInverseTransposed &jacobianInverseTransposed () const;
224 
226  const JacobianInverseTransposed &
227  jacobianInverseTransposed ( const LocalCoordinate &local ) const
228  {
229  return jacobianInverseTransposed();
230  }
231 
232  //***********************************************************************
233  // Methods that not belong to the Interface, but have to be public
234  //***********************************************************************
235 
237  void invalidate ()
238  {
239  builtJT_ = false;
240  builtJTInv_ = false;
241  calcedDet_ = false;
242  }
243 
245  template< class CoordReader >
246  void build ( const CoordReader &coordReader );
247 
248  void print ( std::ostream &out ) const;
249 
250  private:
251  // calculates the volume of the element
252  ctype elDeterminant () const
253  {
254  return std::abs( Alberta::determinant( jacobianTransposed() ) );
255  }
256 
258  CoordMatrix coord_;
259 
261  GlobalCoordinate centroid_;
262 
263  // storage for the transposed of the jacobian
264  mutable JacobianTransposed jT_;
265 
266  // storage for the transposed inverse of the jacboian
267  mutable JacobianInverseTransposed jTInv_;
268 
269  // has jT_ been computed, yet?
270  mutable bool builtJT_;
271  // has jTInv_ been computed, yet?
272  mutable bool builtJTInv_;
273 
274  mutable bool calcedDet_;
275  mutable ctype elDet_;
276  };
277 
278 
279 
280  // AlbertaGridGlobalGeometry
281  // -------------------------
282 
283  template< int mydim, int cdim, class GridImp >
285  : public AlbertaGridGeometry< mydim, cdim, GridImp >
286  {
289 
290  public:
292  : Base()
293  {}
294 
295  template< class CoordReader >
296  AlbertaGridGlobalGeometry ( const CoordReader &coordReader )
297  : Base( coordReader )
298  {}
299  };
300 
301 
302 #if !DUNE_ALBERTA_CACHE_COORDINATES
303  template< int dim, int cdim >
304  class AlbertaGridGlobalGeometry< dim, cdim, const AlbertaGrid< dim, cdim > >
305  {
307 
308  // remember type of the grid
310 
311  // dimension of barycentric coordinates
312  static const int dimbary = dim + 1;
313 
315 
316  public:
319 
320  static const int dimension = Grid::dimension;
321  static const int mydimension = dimension;
322  static const int codimension = dimension - mydimension;
323  static const int coorddimension = cdim;
324 
325  typedef FieldVector< ctype, mydimension > LocalCoordinate;
326  typedef FieldVector< ctype, coorddimension > GlobalCoordinate;
327 
328  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
329  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
330 
331  private:
332  static const int numCorners = mydimension + 1;
333 
334  typedef FieldMatrix< ctype, numCorners, coorddimension > CoordMatrix;
335 
336  public:
338  {
339  invalidate();
340  }
341 
342  template< class CoordReader >
343  AlbertaGridGlobalGeometry ( const CoordReader &coordReader )
344  {
345  build( coordReader );
346  }
347 
350  {
351  return GeometryTypes::simplex( mydimension );
352  }
353 
355  int corners () const
356  {
357  return numCorners;
358  }
359 
361  GlobalCoordinate corner ( const int i ) const
362  {
363  assert( (i >= 0) && (i < corners()) );
364  const Alberta::GlobalCoordinate &x = elementInfo_.coordinate( i );
365  GlobalCoordinate y;
366  for( int j = 0; j < coorddimension; ++j )
367  y[ j ] = x[ j ];
368  return y;
369  }
370 
372  GlobalCoordinate center () const
373  {
374  GlobalCoordinate centroid_ = corner( 0 );
375  for( int i = 1; i < numCorners; ++i )
376  centroid_ += corner( i );
377  centroid_ *= ctype( 1 ) / ctype( numCorners );
378  return centroid_;
379  }
380 
382  GlobalCoordinate global ( const LocalCoordinate &local ) const;
383 
385  LocalCoordinate local ( const GlobalCoordinate &global ) const;
386 
392  ctype integrationElement () const
393  {
394  return elementInfo_.geometryCache().integrationElement();
395  }
396 
398  ctype integrationElement ( const LocalCoordinate &local ) const
399  {
400  return integrationElement();
401  }
402 
404  ctype volume () const
405  {
406  return integrationElement() / ctype( Factorial< mydimension >::factorial );
407  }
408 
414  const JacobianTransposed &jacobianTransposed () const
415  {
416  return elementInfo_.geometryCache().jacobianTransposed();
417  }
418 
420  const JacobianTransposed &
421  jacobianTransposed ( const LocalCoordinate &local ) const
422  {
423  return jacobianTransposed();
424  }
425 
431  const JacobianInverseTransposed &jacobianInverseTransposed () const
432  {
433  return elementInfo_.geometryCache().jacobianInverseTransposed();
434  }
435 
437  const JacobianInverseTransposed &
438  jacobianInverseTransposed ( const LocalCoordinate &local ) const
439  {
440  return jacobianInverseTransposed();
441  }
442 
443  //***********************************************************************
444  // Methods that not belong to the Interface, but have to be public
445  //***********************************************************************
446 
448  void invalidate ()
449  {
450  elementInfo_ = ElementInfo();
451  }
452 
454  template< class CoordReader >
455  void build ( const CoordReader &coordReader )
456  {
457  elementInfo_ = coordReader.elementInfo();
458  }
459 
460  private:
461  ElementInfo elementInfo_;
462  };
463 #endif // #if !DUNE_ALBERTA_CACHE_COORDINATES
464 
465 
466 
467  // AlbertaGridLocalGeometryProvider
468  // --------------------------------
469 
470  template< class Grid >
472  {
474 
475  public:
476  typedef typename Grid::ctype ctype;
477 
478  static const int dimension = Grid::dimension;
479 
480  template< int codim >
481  struct Codim
482  {
483  typedef AlbertaGridGeometry< dimension-codim, dimension, Grid > LocalGeometry;
484  };
485 
488 
489  static const int numChildren = 2;
490  static const int numFaces = dimension + 1;
491 
492  static const int minFaceTwist = Alberta::Twist< dimension, dimension-1 >::minTwist;
493  static const int maxFaceTwist = Alberta::Twist< dimension, dimension-1 >::maxTwist;
494  static const int numFaceTwists = maxFaceTwist - minFaceTwist + 1;
495 
496  private:
497  struct GeoInFatherCoordReader;
498  struct FaceCoordReader;
499 
501  {
502  buildGeometryInFather();
503  buildFaceGeometry();
504  }
505 
507  {
508  for( int child = 0; child < numChildren; ++child )
509  {
510  delete geometryInFather_[ child ][ 0 ];
511  delete geometryInFather_[ child ][ 1 ];
512  }
513 
514  for( int i = 0; i < numFaces; ++i )
515  {
516  for( int j = 0; j < numFaceTwists; ++j )
517  delete faceGeometry_[ i ][ j ];
518  }
519  }
520 
521  void buildGeometryInFather();
522  void buildFaceGeometry();
523 
524  public:
525  const LocalElementGeometry &
526  geometryInFather ( int child, const int orientation = 1 ) const
527  {
528  assert( (child >= 0) && (child < numChildren) );
529  assert( (orientation == 1) || (orientation == -1) );
530  return *geometryInFather_[ child ][ (orientation + 1) / 2 ];
531  }
532 
533  const LocalFaceGeometry &
534  faceGeometry ( int face, int twist = 0 ) const
535  {
536  assert( (face >= 0) && (face < numFaces) );
537  assert( (twist >= minFaceTwist) && (twist <= maxFaceTwist) );
538  return *faceGeometry_[ face ][ twist - minFaceTwist ];
539  }
540 
541  static const This &instance ()
542  {
543  static This theInstance;
544  return theInstance;
545  }
546 
547  private:
548  template< int codim >
549  static int mapVertices ( int subEntity, int i )
550  {
552  }
553 
554  const LocalElementGeometry *geometryInFather_[ numChildren ][ 2 ];
555  const LocalFaceGeometry *faceGeometry_[ numFaces ][ numFaceTwists ];
556  };
557 
558 } // namespace Dune
559 
560 #endif // #if HAVE_ALBERTA
561 
562 #endif // #ifndef DUNE_ALBERTA_GEOMETRY_HH
Include standard header files.
Definition: agrid.hh:58
AlbertaGridGeometry< dimension-codim, dimension, Grid > LocalGeometry
Definition: albertagrid/geometry.hh:483
AlbertaGridGlobalGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:296
void coordinate(int i, Coordinate &x) const
Definition: albertagrid/geometry.hh:54
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:193
GlobalCoordinate center() const
return center of geometry
Definition: albertagrid/geometry.hh:170
Definition: albertagrid/geometry.hh:471
const ElementInfo & elementInfo() const
Definition: albertagrid/geometry.hh:49
static const int codimension
Definition: albertagrid/geometry.hh:32
ct ctype
Define type used for coordinates in grid module.
Definition: common/grid.hh:522
static K determinant(const FieldMatrix< K, 0, m > &matrix)
Definition: algebra.hh:28
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:438
geometry implementation for AlbertaGrid
Definition: albertagrid/geometry.hh:105
int corners() const
number of corner the geometry
Definition: albertagrid/geometry.hh:355
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
transposed inverse of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:227
bool affine() const
returns always true since we only have simplices
Definition: albertagrid/geometry.hh:154
Definition: misc.hh:441
FieldVector< ctype, coorddimension > Coordinate
Definition: albertagrid/geometry.hh:39
static const int mydimension
Definition: albertagrid/geometry.hh:33
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Definition: albertagrid/geometry.hh:329
const JacobianTransposed & jacobianTransposed() const
transposed of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:414
void invalidate()
invalidate the geometry
Definition: albertagrid/geometry.hh:237
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
Definition: albertagrid/geometry.hh:27
[ provides Dune::Grid ]
Definition: agrid.hh:136
ALBERTA REAL Real
Definition: misc.hh:46
ctype integrationElement() const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:186
ctype integrationElement() const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:392
Alberta::Real ctype
Definition: albertagrid/geometry.hh:36
ctype volume() const
volume of geometry
Definition: albertagrid/geometry.hh:404
AlbertaGridCoordinateReader(const GridImp &grid, const ElementInfo &elementInfo, int subEntity)
Definition: albertagrid/geometry.hh:41
GlobalCoordinate center() const
return center of geometry
Definition: albertagrid/geometry.hh:372
std::remove_const< GridImp >::type Grid
Definition: albertagrid/geometry.hh:29
provides a wrapper for ALBERTA&#39;s el_info structure
Alberta::Real ctype
type of coordinates
Definition: albertagrid/geometry.hh:117
Definition: albertagrid/geometry.hh:481
const LocalElementGeometry & geometryInFather(int child, const int orientation=1) const
Definition: albertagrid/geometry.hh:526
void abs(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:326
FieldVector< ctype, coorddimension > GlobalCoordinate
Definition: albertagrid/geometry.hh:125
AlbertaGridGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:142
GeometryType
Type representing VTK&#39;s entity geometry types.
Definition: common.hh:178
AlbertaGridGeometry()
Definition: albertagrid/geometry.hh:136
Definition: misc.hh:529
GeometryType type() const
obtain the type of reference element
Definition: albertagrid/geometry.hh:349
AlbertaGridGlobalGeometry(const CoordReader &coordReader)
Definition: albertagrid/geometry.hh:343
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Definition: albertagrid/geometry.hh:328
ALBERTA REAL_D GlobalVector
Definition: misc.hh:48
FieldVector< ctype, mydimension > LocalCoordinate
Definition: albertagrid/geometry.hh:124
int corners() const
number of corner the geometry
Definition: albertagrid/geometry.hh:157
FieldVector< ctype, coorddimension > GlobalCoordinate
Definition: albertagrid/geometry.hh:326
The dimension of the world the grid lives in.
Definition: common/grid.hh:393
Grid::ctype ctype
Definition: albertagrid/geometry.hh:476
Codim< 1 >::LocalGeometry LocalFaceGeometry
Definition: albertagrid/geometry.hh:487
ctype integrationElement(const LocalCoordinate &local) const
integration element of the geometry mapping
Definition: albertagrid/geometry.hh:398
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Definition: albertagrid/geometry.hh:128
const JacobianInverseTransposed & jacobianInverseTransposed() const
transposed inverse of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:431
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Definition: albertagrid/geometry.hh:127
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: albertagrid/geometry.hh:163
static const int dimension
Definition: albertagrid/geometry.hh:31
Alberta::ElementInfo< dimension > ElementInfo
Definition: albertagrid/geometry.hh:38
GlobalCoordinate corner(const int i) const
obtain the i-th corner of this geometry
Definition: albertagrid/geometry.hh:361
static const int coorddimension
Definition: albertagrid/geometry.hh:34
GeometryType type() const
obtain the type of reference element
Definition: albertagrid/geometry.hh:148
ctype volume() const
volume of geometry
Definition: albertagrid/geometry.hh:199
The dimension of the grid.
Definition: common/grid.hh:387
Definition: albertagrid/geometry.hh:284
FieldVector< ctype, mydimension > LocalCoordinate
Definition: albertagrid/geometry.hh:325
bool hasDeterminant() const
Definition: albertagrid/geometry.hh:65
void invalidate()
invalidate the geometry
Definition: albertagrid/geometry.hh:448
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:421
AlbertaGridGlobalGeometry()
Definition: albertagrid/geometry.hh:291
const LocalFaceGeometry & faceGeometry(int face, int twist=0) const
Definition: albertagrid/geometry.hh:534
static const This & instance()
Definition: albertagrid/geometry.hh:541
void build(const CoordReader &coordReader)
build the geometry from a coordinate reader
Definition: albertagrid/geometry.hh:455
ctype determinant() const
Definition: albertagrid/geometry.hh:70
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
transposed of the geometry mapping&#39;s Jacobian
Definition: albertagrid/geometry.hh:213
Codim< 0 >::LocalGeometry LocalElementGeometry
Definition: albertagrid/geometry.hh:486
Alberta::Real ctype
type of coordinates
Definition: albertagrid/geometry.hh:318
Wrapper and interface classes for element geometries.