dune-geometry  2.8.0
affinegeometry.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_GEOMETRY_AFFINEGEOMETRY_HH
4 #define DUNE_GEOMETRY_AFFINEGEOMETRY_HH
5 
11 #include <cmath>
12 
13 #include <dune/common/fmatrix.hh>
14 #include <dune/common/fvector.hh>
15 
16 #include <dune/geometry/type.hh>
17 
18 namespace Dune
19 {
20 
21  // External Forward Declarations
22  // -----------------------------
23 
24  namespace Geo
25  {
26 
27  template< typename Implementation >
28  class ReferenceElement;
29 
30  template< class ctype, int dim >
32 
33  template< class ctype, int dim >
34  struct ReferenceElements;
35 
36  }
37 
38 
39  namespace Impl
40  {
41 
42  // FieldMatrixHelper
43  // -----------------
44 
45  template< class ct >
46  struct FieldMatrixHelper
47  {
48  typedef ct ctype;
49 
50  template< int m, int n >
51  static void Ax ( const FieldMatrix< ctype, m, n > &A, const FieldVector< ctype, n > &x, FieldVector< ctype, m > &ret )
52  {
53  for( int i = 0; i < m; ++i )
54  {
55  ret[ i ] = ctype( 0 );
56  for( int j = 0; j < n; ++j )
57  ret[ i ] += A[ i ][ j ] * x[ j ];
58  }
59  }
60 
61  template< int m, int n >
62  static void ATx ( const FieldMatrix< ctype, m, n > &A, const FieldVector< ctype, m > &x, FieldVector< ctype, n > &ret )
63  {
64  for( int i = 0; i < n; ++i )
65  {
66  ret[ i ] = ctype( 0 );
67  for( int j = 0; j < m; ++j )
68  ret[ i ] += A[ j ][ i ] * x[ j ];
69  }
70  }
71 
72  template< int m, int n, int p >
73  static void AB ( const FieldMatrix< ctype, m, n > &A, const FieldMatrix< ctype, n, p > &B, FieldMatrix< ctype, m, p > &ret )
74  {
75  for( int i = 0; i < m; ++i )
76  {
77  for( int j = 0; j < p; ++j )
78  {
79  ret[ i ][ j ] = ctype( 0 );
80  for( int k = 0; k < n; ++k )
81  ret[ i ][ j ] += A[ i ][ k ] * B[ k ][ j ];
82  }
83  }
84  }
85 
86  template< int m, int n, int p >
87  static void ATBT ( const FieldMatrix< ctype, m, n > &A, const FieldMatrix< ctype, p, m > &B, FieldMatrix< ctype, n, p > &ret )
88  {
89  for( int i = 0; i < n; ++i )
90  {
91  for( int j = 0; j < p; ++j )
92  {
93  ret[ i ][ j ] = ctype( 0 );
94  for( int k = 0; k < m; ++k )
95  ret[ i ][ j ] += A[ k ][ i ] * B[ j ][ k ];
96  }
97  }
98  }
99 
100  template< int m, int n >
101  static void ATA_L ( const FieldMatrix< ctype, m, n > &A, FieldMatrix< ctype, n, n > &ret )
102  {
103  for( int i = 0; i < n; ++i )
104  {
105  for( int j = 0; j <= i; ++j )
106  {
107  ret[ i ][ j ] = ctype( 0 );
108  for( int k = 0; k < m; ++k )
109  ret[ i ][ j ] += A[ k ][ i ] * A[ k ][ j ];
110  }
111  }
112  }
113 
114  template< int m, int n >
115  static void ATA ( const FieldMatrix< ctype, m, n > &A, FieldMatrix< ctype, n, n > &ret )
116  {
117  for( int i = 0; i < n; ++i )
118  {
119  for( int j = 0; j <= i; ++j )
120  {
121  ret[ i ][ j ] = ctype( 0 );
122  for( int k = 0; k < m; ++k )
123  ret[ i ][ j ] += A[ k ][ i ] * A[ k ][ j ];
124  ret[ j ][ i ] = ret[ i ][ j ];
125  }
126 
127  ret[ i ][ i ] = ctype( 0 );
128  for( int k = 0; k < m; ++k )
129  ret[ i ][ i ] += A[ k ][ i ] * A[ k ][ i ];
130  }
131  }
132 
133  template< int m, int n >
134  static void AAT_L ( const FieldMatrix< ctype, m, n > &A, FieldMatrix< ctype, m, m > &ret )
135  {
136  /*
137  if (m==2) {
138  ret[0][0] = A[0]*A[0];
139  ret[1][1] = A[1]*A[1];
140  ret[1][0] = A[0]*A[1];
141  }
142  else
143  */
144  for( int i = 0; i < m; ++i )
145  {
146  for( int j = 0; j <= i; ++j )
147  {
148  ctype &retij = ret[ i ][ j ];
149  retij = A[ i ][ 0 ] * A[ j ][ 0 ];
150  for( int k = 1; k < n; ++k )
151  retij += A[ i ][ k ] * A[ j ][ k ];
152  }
153  }
154  }
155 
156  template< int m, int n >
157  static void AAT ( const FieldMatrix< ctype, m, n > &A, FieldMatrix< ctype, m, m > &ret )
158  {
159  for( int i = 0; i < m; ++i )
160  {
161  for( int j = 0; j < i; ++j )
162  {
163  ret[ i ][ j ] = ctype( 0 );
164  for( int k = 0; k < n; ++k )
165  ret[ i ][ j ] += A[ i ][ k ] * A[ j ][ k ];
166  ret[ j ][ i ] = ret[ i ][ j ];
167  }
168  ret[ i ][ i ] = ctype( 0 );
169  for( int k = 0; k < n; ++k )
170  ret[ i ][ i ] += A[ i ][ k ] * A[ i ][ k ];
171  }
172  }
173 
174  template< int n >
175  static void Lx ( const FieldMatrix< ctype, n, n > &L, const FieldVector< ctype, n > &x, FieldVector< ctype, n > &ret )
176  {
177  for( int i = 0; i < n; ++i )
178  {
179  ret[ i ] = ctype( 0 );
180  for( int j = 0; j <= i; ++j )
181  ret[ i ] += L[ i ][ j ] * x[ j ];
182  }
183  }
184 
185  template< int n >
186  static void LTx ( const FieldMatrix< ctype, n, n > &L, const FieldVector< ctype, n > &x, FieldVector< ctype, n > &ret )
187  {
188  for( int i = 0; i < n; ++i )
189  {
190  ret[ i ] = ctype( 0 );
191  for( int j = i; j < n; ++j )
192  ret[ i ] += L[ j ][ i ] * x[ j ];
193  }
194  }
195 
196  template< int n >
197  static void LTL ( const FieldMatrix< ctype, n, n > &L, FieldMatrix< ctype, n, n > &ret )
198  {
199  for( int i = 0; i < n; ++i )
200  {
201  for( int j = 0; j < i; ++j )
202  {
203  ret[ i ][ j ] = ctype( 0 );
204  for( int k = i; k < n; ++k )
205  ret[ i ][ j ] += L[ k ][ i ] * L[ k ][ j ];
206  ret[ j ][ i ] = ret[ i ][ j ];
207  }
208  ret[ i ][ i ] = ctype( 0 );
209  for( int k = i; k < n; ++k )
210  ret[ i ][ i ] += L[ k ][ i ] * L[ k ][ i ];
211  }
212  }
213 
214  template< int n >
215  static void LLT ( const FieldMatrix< ctype, n, n > &L, FieldMatrix< ctype, n, n > &ret )
216  {
217  for( int i = 0; i < n; ++i )
218  {
219  for( int j = 0; j < i; ++j )
220  {
221  ret[ i ][ j ] = ctype( 0 );
222  for( int k = 0; k <= j; ++k )
223  ret[ i ][ j ] += L[ i ][ k ] * L[ j ][ k ];
224  ret[ j ][ i ] = ret[ i ][ j ];
225  }
226  ret[ i ][ i ] = ctype( 0 );
227  for( int k = 0; k <= i; ++k )
228  ret[ i ][ i ] += L[ i ][ k ] * L[ i ][ k ];
229  }
230  }
231 
232  template< int n >
233  static bool cholesky_L ( const FieldMatrix< ctype, n, n > &A, FieldMatrix< ctype, n, n > &ret, const bool checkSingular = false )
234  {
235  using std::sqrt;
236  for( int i = 0; i < n; ++i )
237  {
238  ctype &rii = ret[ i ][ i ];
239 
240  ctype xDiag = A[ i ][ i ];
241  for( int j = 0; j < i; ++j )
242  xDiag -= ret[ i ][ j ] * ret[ i ][ j ];
243 
244  // in some cases A can be singular, e.g. when checking local for
245  // outside points during checkInside
246  if( checkSingular && ! ( xDiag > ctype( 0 )) )
247  return false ;
248 
249  // otherwise this should be true always
250  assert( xDiag > ctype( 0 ) );
251  rii = sqrt( xDiag );
252 
253  ctype invrii = ctype( 1 ) / rii;
254  for( int k = i+1; k < n; ++k )
255  {
256  ctype x = A[ k ][ i ];
257  for( int j = 0; j < i; ++j )
258  x -= ret[ i ][ j ] * ret[ k ][ j ];
259  ret[ k ][ i ] = invrii * x;
260  }
261  }
262 
263  // return true for meaning A is non-singular
264  return true;
265  }
266 
267  template< int n >
268  static ctype detL ( const FieldMatrix< ctype, n, n > &L )
269  {
270  ctype det( 1 );
271  for( int i = 0; i < n; ++i )
272  det *= L[ i ][ i ];
273  return det;
274  }
275 
276  template< int n >
277  static ctype invL ( FieldMatrix< ctype, n, n > &L )
278  {
279  ctype det( 1 );
280  for( int i = 0; i < n; ++i )
281  {
282  ctype &lii = L[ i ][ i ];
283  det *= lii;
284  lii = ctype( 1 ) / lii;
285  for( int j = 0; j < i; ++j )
286  {
287  ctype &lij = L[ i ][ j ];
288  ctype x = lij * L[ j ][ j ];
289  for( int k = j+1; k < i; ++k )
290  x += L[ i ][ k ] * L[ k ][ j ];
291  lij = (-lii) * x;
292  }
293  }
294  return det;
295  }
296 
297  // calculates x := L^{-1} x
298  template< int n >
299  static void invLx ( FieldMatrix< ctype, n, n > &L, FieldVector< ctype, n > &x )
300  {
301  for( int i = 0; i < n; ++i )
302  {
303  for( int j = 0; j < i; ++j )
304  x[ i ] -= L[ i ][ j ] * x[ j ];
305  x[ i ] /= L[ i ][ i ];
306  }
307  }
308 
309  // calculates x := L^{-T} x
310  template< int n >
311  static void invLTx ( FieldMatrix< ctype, n, n > &L, FieldVector< ctype, n > &x )
312  {
313  for( int i = n; i > 0; --i )
314  {
315  for( int j = i; j < n; ++j )
316  x[ i-1 ] -= L[ j ][ i-1 ] * x[ j ];
317  x[ i-1 ] /= L[ i-1 ][ i-1 ];
318  }
319  }
320 
321  template< int n >
322  static ctype spdDetA ( const FieldMatrix< ctype, n, n > &A )
323  {
324  // return A[0][0]*A[1][1]-A[1][0]*A[1][0];
325  FieldMatrix< ctype, n, n > L;
326  cholesky_L( A, L );
327  return detL( L );
328  }
329 
330  template< int n >
331  static ctype spdInvA ( FieldMatrix< ctype, n, n > &A )
332  {
333  FieldMatrix< ctype, n, n > L;
334  cholesky_L( A, L );
335  const ctype det = invL( L );
336  LTL( L, A );
337  return det;
338  }
339 
340  // calculate x := A^{-1} x
341  template< int n >
342  static bool spdInvAx ( FieldMatrix< ctype, n, n > &A, FieldVector< ctype, n > &x, const bool checkSingular = false )
343  {
344  FieldMatrix< ctype, n, n > L;
345  const bool invertible = cholesky_L( A, L, checkSingular );
346  if( ! invertible ) return invertible ;
347  invLx( L, x );
348  invLTx( L, x );
349  return invertible;
350  }
351 
352  template< int m, int n >
353  static ctype detATA ( const FieldMatrix< ctype, m, n > &A )
354  {
355  if( m >= n )
356  {
357  FieldMatrix< ctype, n, n > ata;
358  ATA_L( A, ata );
359  return spdDetA( ata );
360  }
361  else
362  return ctype( 0 );
363  }
364 
370  template< int m, int n >
371  static ctype sqrtDetAAT ( const FieldMatrix< ctype, m, n > &A )
372  {
373  using std::abs;
374  using std::sqrt;
375  // These special cases are here not only for speed reasons:
376  // The general implementation aborts if the matrix is almost singular,
377  // and the special implementation provide a stable way to handle that case.
378  if( (n == 2) && (m == 2) )
379  {
380  // Special implementation for 2x2 matrices: faster and more stable
381  return abs( A[ 0 ][ 0 ]*A[ 1 ][ 1 ] - A[ 1 ][ 0 ]*A[ 0 ][ 1 ] );
382  }
383  else if( (n == 3) && (m == 3) )
384  {
385  // Special implementation for 3x3 matrices
386  const ctype v0 = A[ 0 ][ 1 ] * A[ 1 ][ 2 ] - A[ 1 ][ 1 ] * A[ 0 ][ 2 ];
387  const ctype v1 = A[ 0 ][ 2 ] * A[ 1 ][ 0 ] - A[ 1 ][ 2 ] * A[ 0 ][ 0 ];
388  const ctype v2 = A[ 0 ][ 0 ] * A[ 1 ][ 1 ] - A[ 1 ][ 0 ] * A[ 0 ][ 1 ];
389  return abs( v0 * A[ 2 ][ 0 ] + v1 * A[ 2 ][ 1 ] + v2 * A[ 2 ][ 2 ] );
390  }
391  else if ( (n == 3) && (m == 2) )
392  {
393  // Special implementation for 2x3 matrices
394  const ctype v0 = A[ 0 ][ 0 ] * A[ 1 ][ 1 ] - A[ 0 ][ 1 ] * A[ 1 ][ 0 ];
395  const ctype v1 = A[ 0 ][ 0 ] * A[ 1 ][ 2 ] - A[ 1 ][ 0 ] * A[ 0 ][ 2 ];
396  const ctype v2 = A[ 0 ][ 1 ] * A[ 1 ][ 2 ] - A[ 0 ][ 2 ] * A[ 1 ][ 1 ];
397  return sqrt( v0*v0 + v1*v1 + v2*v2);
398  }
399  else if( n >= m )
400  {
401  // General case
402  FieldMatrix< ctype, m, m > aat;
403  AAT_L( A, aat );
404  return spdDetA( aat );
405  }
406  else
407  return ctype( 0 );
408  }
409 
410  // A^{-1}_L = (A^T A)^{-1} A^T
411  // => A^{-1}_L A = I
412  template< int m, int n >
413  static ctype leftInvA ( const FieldMatrix< ctype, m, n > &A, FieldMatrix< ctype, n, m > &ret )
414  {
415  static_assert((m >= n), "Matrix has no left inverse.");
416  FieldMatrix< ctype, n, n > ata;
417  ATA_L( A, ata );
418  const ctype det = spdInvA( ata );
419  ATBT( ata, A, ret );
420  return det;
421  }
422 
423  template< int m, int n >
424  static void leftInvAx ( const FieldMatrix< ctype, m, n > &A, const FieldVector< ctype, m > &x, FieldVector< ctype, n > &y )
425  {
426  static_assert((m >= n), "Matrix has no left inverse.");
427  FieldMatrix< ctype, n, n > ata;
428  ATx( A, x, y );
429  ATA_L( A, ata );
430  spdInvAx( ata, y );
431  }
432 
434  template< int m, int n >
435  static ctype rightInvA ( const FieldMatrix< ctype, m, n > &A, FieldMatrix< ctype, n, m > &ret )
436  {
437  static_assert((n >= m), "Matrix has no right inverse.");
438  using std::abs;
439  if( (n == 2) && (m == 2) )
440  {
441  const ctype det = (A[ 0 ][ 0 ]*A[ 1 ][ 1 ] - A[ 1 ][ 0 ]*A[ 0 ][ 1 ]);
442  const ctype detInv = ctype( 1 ) / det;
443  ret[ 0 ][ 0 ] = A[ 1 ][ 1 ] * detInv;
444  ret[ 1 ][ 1 ] = A[ 0 ][ 0 ] * detInv;
445  ret[ 1 ][ 0 ] = -A[ 1 ][ 0 ] * detInv;
446  ret[ 0 ][ 1 ] = -A[ 0 ][ 1 ] * detInv;
447  return abs( det );
448  }
449  else
450  {
451  FieldMatrix< ctype, m , m > aat;
452  AAT_L( A, aat );
453  const ctype det = spdInvA( aat );
454  ATBT( A , aat , ret );
455  return det;
456  }
457  }
458 
459  template< int m, int n >
460  static bool xTRightInvA ( const FieldMatrix< ctype, m, n > &A, const FieldVector< ctype, n > &x, FieldVector< ctype, m > &y )
461  {
462  static_assert((n >= m), "Matrix has no right inverse.");
463  FieldMatrix< ctype, m, m > aat;
464  Ax( A, x, y );
465  AAT_L( A, aat );
466  // check whether aat is singular and return true if non-singular
467  return spdInvAx( aat, y, true );
468  }
469  };
470 
471  } // namespace Impl
472 
473 
474 
480  template< class ct, int mydim, int cdim>
482  {
483  public:
484 
486  typedef ct ctype;
487 
489  static const int mydimension= mydim;
490 
492  static const int coorddimension = cdim;
493 
495  typedef FieldVector< ctype, mydimension > LocalCoordinate;
496 
498  typedef FieldVector< ctype, coorddimension > GlobalCoordinate;
499 
501  typedef ctype Volume;
502 
504  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
505 
507  typedef FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed;
508 
509  private:
512 
514 
515  // Helper class to compute a matrix pseudo inverse
516  typedef Impl::FieldMatrixHelper< ct > MatrixHelper;
517 
518  public:
520  AffineGeometry ( const ReferenceElement &refElement, const GlobalCoordinate &origin,
521  const JacobianTransposed &jt )
522  : refElement_(refElement), origin_(origin), jacobianTransposed_(jt)
523  {
524  integrationElement_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jacobianTransposed_, jacobianInverseTransposed_ );
525  }
526 
529  const JacobianTransposed &jt )
530  : AffineGeometry(ReferenceElements::general( gt ), origin, jt)
531  { }
532 
534  template< class CoordVector >
535  AffineGeometry ( const ReferenceElement &refElement, const CoordVector &coordVector )
536  : refElement_(refElement), origin_(coordVector[0])
537  {
538  for( int i = 0; i < mydimension; ++i )
539  jacobianTransposed_[ i ] = coordVector[ i+1 ] - origin_;
540  integrationElement_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jacobianTransposed_, jacobianInverseTransposed_ );
541  }
542 
544  template< class CoordVector >
545  AffineGeometry ( Dune::GeometryType gt, const CoordVector &coordVector )
546  : AffineGeometry(ReferenceElements::general( gt ), coordVector)
547  { }
548 
550  bool affine () const { return true; }
551 
553  Dune::GeometryType type () const { return refElement_.type(); }
554 
556  int corners () const { return refElement_.size( mydimension ); }
557 
559  GlobalCoordinate corner ( int i ) const
560  {
561  return global( refElement_.position( i, mydimension ) );
562  }
563 
565  GlobalCoordinate center () const { return global( refElement_.position( 0, 0 ) ); }
566 
573  GlobalCoordinate global ( const LocalCoordinate &local ) const
574  {
575  GlobalCoordinate global( origin_ );
576  jacobianTransposed_.umtv( local, global );
577  return global;
578  }
579 
593  LocalCoordinate local ( const GlobalCoordinate &global ) const
594  {
595  LocalCoordinate local;
596  jacobianInverseTransposed_.mtv( global - origin_, local );
597  return local;
598  }
599 
610  ctype integrationElement ([[maybe_unused]] const LocalCoordinate &local) const
611  {
612  return integrationElement_;
613  }
614 
616  Volume volume () const
617  {
618  return integrationElement_ * refElement_.volume();
619  }
620 
627  const JacobianTransposed &jacobianTransposed ([[maybe_unused]] const LocalCoordinate &local) const
628  {
629  return jacobianTransposed_;
630  }
631 
638  const JacobianInverseTransposed &jacobianInverseTransposed ([[maybe_unused]] const LocalCoordinate &local) const
639  {
640  return jacobianInverseTransposed_;
641  }
642 
644  {
645  return geometry.refElement_;
646  }
647 
648  private:
649  ReferenceElement refElement_;
650  GlobalCoordinate origin_;
651  JacobianTransposed jacobianTransposed_;
652  JacobianInverseTransposed jacobianInverseTransposed_;
653  ctype integrationElement_;
654  };
655 
656 } // namespace Dune
657 
658 #endif // #ifndef DUNE_GEOMETRY_AFFINEGEOMETRY_HH
A unique label for each type of element that can occur in a grid.
unspecified-type ReferenceElement
Returns the type of reference element for the argument type T.
Definition: referenceelements.hh:495
Definition: affinegeometry.hh:19
This class provides access to geometric and topological properties of a reference element.
Definition: referenceelement.hh:25
CoordinateField volume() const
obtain the volume of the reference element
Definition: referenceelement.hh:239
decltype(auto) type(int i, int c) const
obtain the type of subentity (i,c)
Definition: referenceelement.hh:169
int size(int c) const
number of subentities of codimension c
Definition: referenceelement.hh:92
decltype(auto) position(int i, int c) const
position of the barycenter of entity (i,c)
Definition: referenceelement.hh:201
Definition: affinegeometry.hh:31
Class providing access to the singletons of the reference elements.
Definition: referenceelements.hh:168
Implementation of the Geometry interface for affine geometries.
Definition: affinegeometry.hh:482
AffineGeometry(const ReferenceElement &refElement, const CoordVector &coordVector)
Create affine geometry from reference element and a vector of vertex coordinates.
Definition: affinegeometry.hh:535
AffineGeometry(Dune::GeometryType gt, const GlobalCoordinate &origin, const JacobianTransposed &jt)
Create affine geometry from GeometryType, one vertex, and the Jacobian matrix.
Definition: affinegeometry.hh:528
FieldVector< ctype, mydimension > LocalCoordinate
Type for local coordinate vector.
Definition: affinegeometry.hh:495
Dune::GeometryType type() const
Obtain the type of the reference element.
Definition: affinegeometry.hh:553
static const int mydimension
Dimension of the geometry.
Definition: affinegeometry.hh:489
AffineGeometry(const ReferenceElement &refElement, const GlobalCoordinate &origin, const JacobianTransposed &jt)
Create affine geometry from reference element, one vertex, and the Jacobian matrix.
Definition: affinegeometry.hh:520
ctype Volume
Type used for volume.
Definition: affinegeometry.hh:501
friend ReferenceElement referenceElement(const AffineGeometry &geometry)
Definition: affinegeometry.hh:643
AffineGeometry(Dune::GeometryType gt, const CoordVector &coordVector)
Create affine geometry from GeometryType and a vector of vertex coordinates.
Definition: affinegeometry.hh:545
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
Type for the transposed Jacobian matrix.
Definition: affinegeometry.hh:504
GlobalCoordinate corner(int i) const
Obtain coordinates of the i-th corner.
Definition: affinegeometry.hh:559
int corners() const
Obtain number of corners of the corresponding reference element.
Definition: affinegeometry.hh:556
FieldMatrix< ctype, coorddimension, mydimension > JacobianInverseTransposed
Type for the transposed inverse Jacobian matrix.
Definition: affinegeometry.hh:507
static const int coorddimension
Dimension of the world space.
Definition: affinegeometry.hh:492
GlobalCoordinate global(const LocalCoordinate &local) const
Evaluate the mapping.
Definition: affinegeometry.hh:573
GlobalCoordinate center() const
Obtain the centroid of the mapping's image.
Definition: affinegeometry.hh:565
ct ctype
Type used for coordinates.
Definition: affinegeometry.hh:486
FieldVector< ctype, coorddimension > GlobalCoordinate
Type for coordinate vector in world space.
Definition: affinegeometry.hh:498
const JacobianInverseTransposed & jacobianInverseTransposed([[maybe_unused]] const LocalCoordinate &local) const
Obtain the transposed of the Jacobian's inverse.
Definition: affinegeometry.hh:638
ctype integrationElement([[maybe_unused]] const LocalCoordinate &local) const
Obtain the integration element.
Definition: affinegeometry.hh:610
bool affine() const
Always true: this is an affine geometry.
Definition: affinegeometry.hh:550
Volume volume() const
Obtain the volume of the element.
Definition: affinegeometry.hh:616
const JacobianTransposed & jacobianTransposed([[maybe_unused]] const LocalCoordinate &local) const
Obtain the transposed of the Jacobian.
Definition: affinegeometry.hh:627
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:123