dune-geometry  2.3.1-rc1
multilineargeometry.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_MULTILINEARGEOMETRY_HH
4 #define DUNE_GEOMETRY_MULTILINEARGEOMETRY_HH
5 
6 #include <cassert>
7 #include <limits>
8 #include <vector>
9 
10 #include <dune/common/fmatrix.hh>
11 #include <dune/common/fvector.hh>
12 #include <dune/common/typetraits.hh>
13 
15 #include <dune/geometry/type.hh>
18 
19 namespace Dune
20 {
21 
22  // External Forward Declarations
23  // -----------------------------
24 
25  template< class ctype, int dim >
26  class ReferenceElement;
27 
28  template< class ctype, int dim >
29  struct ReferenceElements;
30 
31 
32 
33  // MultiLinearGeometryTraits
34  // -------------------------
35 
45  template< class ct >
47  {
67 
69  static ct tolerance () { return ct( 16 ) * std::numeric_limits< ct >::epsilon(); }
70 
94  template< int mydim, int cdim >
96  {
97  typedef std::vector< FieldVector< ct, cdim > > Type;
98  };
99 
113  template< int dim >
115  {
116  static const bool v = false;
117  static const unsigned int topologyId = ~0u;
118  };
119  };
120 
121 
122 
123  // MultiLinearGeometry
124  // -------------------
125 
149  template< class ct, int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct > >
151  {
153 
154  public:
156  typedef ct ctype;
157 
159  static const int mydimension= mydim;
161  static const int coorddimension = cdim;
162 
167  typedef FieldVector< ctype, mydimension > LocalCoordinate;
170  typedef FieldVector< ctype, coorddimension > GlobalCoordinate;
171 
173  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
174 
177 
182 
185 
186  private:
187  static const bool hasSingleGeometryType = Traits::template hasSingleGeometryType< mydimension >::v;
188 
189  protected:
190  typedef typename Traits::MatrixHelper MatrixHelper;
191  typedef typename conditional< hasSingleGeometryType, integral_constant< unsigned int, Traits::template hasSingleGeometryType< mydimension >::topologyId >, unsigned int >::type TopologyId;
192 
194 
195  private:
196  typedef typename Traits::template CornerStorage< mydimension, coorddimension >::Type::const_iterator CornerIterator;
197 
198  public:
208  template< class Corners >
210  const Corners &corners )
211  : refElement_( &refElement ),
212  corners_( corners )
213  {}
214 
224  template< class Corners >
226  const Corners &corners )
227  : refElement_( &ReferenceElements::general( gt ) ),
228  corners_( corners )
229  {}
230 
232  bool affine () const
233  {
234  CornerIterator cit = corners_.begin();
235  return affine( topologyId(), integral_constant< int, mydimension >(), cit, jacobianTransposed_ );
236  }
237 
240 
242  int corners () const { return refElement().size( mydimension ); }
243 
245  GlobalCoordinate corner ( int i ) const
246  {
247  assert( (i >= 0) && (i < corners()) );
248  return corners_[ i ];
249  }
250 
252  GlobalCoordinate center () const { return global( refElement().position( 0, 0 ) ); }
253 
261  {
262  CornerIterator cit = corners_.begin();
264  global< false >( topologyId(), integral_constant< int, mydimension >(), cit, ctype( 1 ), local, ctype( 1 ), y );
265  return y;
266  }
267 
281  {
282  const ctype tolerance = Traits::tolerance();
283  LocalCoordinate x = refElement().position( 0, 0 );
284  LocalCoordinate dx;
285  do
286  {
287  // Newton's method: DF^n dx^n = F^n, x^{n+1} -= dx^n
288  const GlobalCoordinate dglobal = (*this).global( x ) - global;
289  MatrixHelper::template xTRightInvA< mydimension, coorddimension >( jacobianTransposed( x ), dglobal, dx );
290  x -= dx;
291  } while( dx.two_norm2() > tolerance );
292  return x;
293  }
294 
310  {
311  return MatrixHelper::template sqrtDetAAT< mydimension, coorddimension >( jacobianTransposed( local ) );
312  }
313 
322  ctype volume () const
323  {
324  return integrationElement( refElement().position( 0, 0 ) ) * refElement().volume();
325  }
326 
337  {
338  CornerIterator cit = corners_.begin();
339  jacobianTransposed< false >( topologyId(), integral_constant< int, mydimension >(), cit, ctype( 1 ), local, ctype( 1 ), jacobianTransposed_ );
340  return jacobianTransposed_;
341  }
342 
349  const JacobianInverseTransposed &jacobianInverseTransposed ( const LocalCoordinate &local ) const;
350 
351  protected:
352  const ReferenceElement &refElement () const { return *refElement_; }
353 
355  {
356  return topologyId( integral_constant< bool, hasSingleGeometryType >() );
357  }
358 
359  template< bool add, int dim >
360  static void global ( TopologyId topologyId, integral_constant< int, dim >,
361  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
362  const ctype &rf, GlobalCoordinate &y );
363  template< bool add >
364  static void global ( TopologyId topologyId, integral_constant< int, 0 >,
365  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
366  const ctype &rf, GlobalCoordinate &y );
367 
368  template< bool add, int rows, int dim >
369  static void jacobianTransposed ( TopologyId topologyId, integral_constant< int, dim >,
370  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
371  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt );
372  template< bool add, int rows >
373  static void jacobianTransposed ( TopologyId topologyId, integral_constant< int, 0 >,
374  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
375  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt );
376 
377  template< int dim >
378  static bool affine ( TopologyId topologyId, integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt );
379  static bool affine ( TopologyId topologyId, integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt );
380 
381  protected:
382  TopologyId topologyId ( integral_constant< bool, true > ) const { return TopologyId(); }
383  unsigned int topologyId ( integral_constant< bool, false > ) const { return refElement().type().id(); }
384 
387 
388  private:
389  const ReferenceElement *refElement_;
390  typename Traits::template CornerStorage< mydimension, coorddimension >::Type corners_;
391  };
392 
393 
394 
395  // MultiLinearGeometry::JacobianInverseTransposed
396  // ----------------------------------------------
397 
398  template< class ct, int mydim, int cdim, class Traits >
399  class MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed
400  : public FieldMatrix< ctype, coorddimension, mydimension >
401  {
402  typedef FieldMatrix< ctype, coorddimension, mydimension > Base;
403 
404  public:
405  void setup ( const JacobianTransposed &jt )
406  {
407  detInv_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jt, static_cast< Base & >( *this ) );
408  }
409 
411  {
412  detInv_ = MatrixHelper::template sqrtDetAAT< mydimension, coorddimension >( jt );
413  }
414 
415  ctype det () const { return ctype( 1 ) / detInv_; }
416  ctype detInv () const { return detInv_; }
417 
418  private:
419  ctype detInv_;
420  };
421 
422 
423 
436  template< class ct, int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct > >
438  : public MultiLinearGeometry< ct, mydim, cdim, Traits >
439  {
442 
443  protected:
445 
446  public:
448 
449  typedef typename Base::ctype ctype;
450 
451  using Base::mydimension;
452  using Base::coorddimension;
453 
456 
459 
460  template< class CornerStorage >
461  CachedMultiLinearGeometry ( const ReferenceElement &refElement, const CornerStorage &cornerStorage )
462  : Base( refElement, cornerStorage ),
463  affine_( Base::affine() ),
464  jacobianInverseTransposedComputed_( false ),
465  integrationElementComputed_( false )
466  {}
467 
468  template< class CornerStorage >
469  CachedMultiLinearGeometry ( Dune::GeometryType gt, const CornerStorage &cornerStorage )
470  : Base( gt, cornerStorage ),
471  affine_( Base::affine() ),
472  jacobianInverseTransposedComputed_( false ),
473  integrationElementComputed_( false )
474  {}
475 
477  bool affine () const { return affine_; }
478 
479  using Base::corner;
480 
482  GlobalCoordinate center () const { return global( refElement().position( 0, 0 ) ); }
483 
491  {
492  if( affine() )
493  {
495  jacobianTransposed_.umtv( local, global );
496  return global;
497  }
498  else
499  return Base::global( local );
500  }
501 
515  {
516  if( affine() )
517  {
519  if( jacobianInverseTransposedComputed_ )
520  jacobianInverseTransposed_.mtv( global - corner( 0 ), local );
521  else
522  MatrixHelper::template xTRightInvA< mydimension, coorddimension >( jacobianTransposed_, global - corner( 0 ), local );
523  return local;
524  }
525  else
526  return Base::local( global );
527  }
528 
544  {
545  if( affine() )
546  {
547  if( !integrationElementComputed_ )
548  {
550  integrationElementComputed_ = true;
551  }
553  }
554  else
555  return Base::integrationElement( local );
556  }
557 
559  ctype volume () const
560  {
561  if( affine() )
562  return integrationElement( refElement().position( 0, 0 ) ) * refElement().volume();
563  else
564  return Base::volume();
565  }
566 
577  {
578  if( affine() )
579  return jacobianTransposed_;
580  else
581  return Base::jacobianTransposed( local );
582  }
583 
591  {
592  if( affine() )
593  {
594  if( !jacobianInverseTransposedComputed_ )
595  {
597  jacobianInverseTransposedComputed_ = true;
598  integrationElementComputed_ = true;
599  }
601  }
602  else
603  return Base::jacobianInverseTransposed( local );
604  }
605 
606  protected:
607  using Base::refElement;
608 
611 
612  private:
613  mutable bool affine_ : 1;
614 
615  mutable bool jacobianInverseTransposedComputed_ : 1;
616  mutable bool integrationElementComputed_ : 1;
617  };
618 
619 
620 
621  // Implementation of MultiLinearGeometry
622  // -------------------------------------
623 
624  template< class ct, int mydim, int cdim, class Traits >
625  inline const typename MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed &
628  {
629  jacobianInverseTransposed_.setup( jacobianTransposed( local ) );
630  return jacobianInverseTransposed_;
631  }
632 
633 
634  template< class ct, int mydim, int cdim, class Traits >
635  template< bool add, int dim >
637  ::global ( TopologyId topologyId, integral_constant< int, dim >,
638  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
639  const ctype &rf, GlobalCoordinate &y )
640  {
641  const ctype xn = df*x[ dim-1 ];
642  const ctype cxn = ctype( 1 ) - xn;
643 
644  if( GenericGeometry::isPrism( topologyId, mydimension, mydimension-dim ) )
645  {
646  // apply (1-xn) times mapping for bottom
647  global< add >( topologyId, integral_constant< int, dim-1 >(), cit, df, x, rf*cxn, y );
648  // apply xn times mapping for top
649  global< true >( topologyId, integral_constant< int, dim-1 >(), cit, df, x, rf*xn, y );
650  }
651  else
652  {
653  assert( GenericGeometry::isPyramid( topologyId, mydimension, mydimension-dim ) );
654  // apply (1-xn) times mapping for bottom (with argument x/(1-xn))
655  if( cxn > Traits::tolerance() || cxn < -Traits::tolerance() )
656  global< add >( topologyId, integral_constant< int, dim-1 >(), cit, df/cxn, x, rf*cxn, y );
657  else
658  global< add >( topologyId, integral_constant< int, dim-1 >(), cit, df, x, ctype( 0 ), y );
659  // apply xn times the tip
660  y.axpy( rf*xn, *cit );
661  ++cit;
662  }
663  }
664 
665  template< class ct, int mydim, int cdim, class Traits >
666  template< bool add >
668  ::global ( TopologyId topologyId, integral_constant< int, 0 >,
669  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
670  const ctype &rf, GlobalCoordinate &y )
671  {
672  const GlobalCoordinate &origin = *cit;
673  ++cit;
674  for( int i = 0; i < coorddimension; ++i )
675  y[ i ] = (add ? y[ i ] + rf*origin[ i ] : rf*origin[ i ]);
676  }
677 
678 
679  template< class ct, int mydim, int cdim, class Traits >
680  template< bool add, int rows, int dim >
682  ::jacobianTransposed ( TopologyId topologyId, integral_constant< int, dim >,
683  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
684  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt )
685  {
686  assert( rows >= dim );
687 
688  const ctype xn = df*x[ dim-1 ];
689  const ctype cxn = ctype( 1 ) - xn;
690 
691  CornerIterator cit2( cit );
692  if( GenericGeometry::isPrism( topologyId, mydimension, mydimension-dim ) )
693  {
694  // apply (1-xn) times Jacobian for bottom
695  jacobianTransposed< add >( topologyId, integral_constant< int, dim-1 >(), cit2, df, x, rf*cxn, jt );
696  // apply xn times Jacobian for top
697  jacobianTransposed< true >( topologyId, integral_constant< int, dim-1 >(), cit2, df, x, rf*xn, jt );
698  // compute last row as difference between top value and bottom value
699  global< add >( topologyId, integral_constant< int, dim-1 >(), cit, df, x, -rf, jt[ dim-1 ] );
700  global< true >( topologyId, integral_constant< int, dim-1 >(), cit, df, x, rf, jt[ dim-1 ] );
701  }
702  else
703  {
704  assert( GenericGeometry::isPyramid( topologyId, mydimension, mydimension-dim ) );
705  /*
706  * In the pyramid case, we need a transformation Tb: B -> R^n for the
707  * base B \subset R^{n-1}. The pyramid transformation is then defined as
708  * T: P \subset R^n -> R^n
709  * (x, xn) |-> (1-xn) Tb(x*) + xn t (x \in R^{n-1}, xn \in R)
710  * with the tip of the pyramid mapped to t and x* = x/(1-xn)
711  * the projection of (x,xn) onto the base.
712  *
713  * For the Jacobi matrix DT we get
714  * DT = ( A | b )
715  * with A = DTb(x*) (n x n-1 matrix)
716  * and b = dT/dxn (n-dim column vector).
717  * Furthermore
718  * b = -Tb(x*) + t + \sum_i dTb/dx_i(x^*) x_i/(1-xn)
719  *
720  * Note that both A and b are not defined in the pyramid tip (x=0, xn=1)!
721  * Indeed for B the unit square, Tb mapping B to the quadrilateral given
722  * by the vertices (0,0,0), (2,0,0), (0,1,0), (1,1,0) and t=(0,0,1), we get
723  *
724  * T(x,y,xn) = ( x(2-y/(1-xn)), y, xn )
725  * / 2-y/(1-xn) -x 0 \
726  * DT(x,y,xn) = | 0 1 0 |
727  * \ 0 0 1 /
728  * which is not continuous for xn -> 1, choose for example
729  * x=0, y=1-xn, xn -> 1 --> DT -> diag(1,1,1)
730  * x=1-xn, y=0, xn -> 1 --> DT -> diag(2,1,1)
731  *
732  * However, for Tb affine-linear, Tb(y) = My + y0, DTb = M:
733  * A = M
734  * b = -M x* - y0 + t + \sum_i M_i x_i/(1-xn)
735  * = -M x* - y0 + t + M x*
736  * = -y0 + t
737  * which is continuous for xn -> 1. Note that this b is also given by
738  * b = -Tb(0) + t + \sum_i dTb/dx_i(0) x_i/1
739  * that is replacing x* by 1 and 1-xn by 1 in the formular above.
740  *
741  * For xn -> 1, we can thus set x*=0, "1-xn"=1 (or anything != 0) and get
742  * the right result in case Tb is affine-linear.
743  */
744 
745  /* The second case effectively results in x* = 0 */
746  ctype dfcxn = (cxn > Traits::tolerance() || cxn < -Traits::tolerance()) ? df / cxn : 0;
747 
748  // initialize last row
749  // b = -Tb(x*)
750  // (b = -Tb(0) = -y0 in case xn -> 1 and Tb affine-linear)
751  global< add >( topologyId, integral_constant< int, dim-1 >(), cit, dfcxn, x, -rf, jt[ dim-1 ] );
752  // b += t
753  jt[ dim-1 ].axpy( rf, *cit );
754  ++cit;
755  // apply Jacobian for bottom (with argument x/(1-xn)) and correct last row
756  if( add )
757  {
758  FieldMatrix< ctype, dim-1, coorddimension > jt2;
759  // jt2 = dTb/dx_i(x*)
760  jacobianTransposed< false >( topologyId, integral_constant< int, dim-1 >(), cit2, dfcxn, x, rf, jt2 );
761  // A = dTb/dx_i(x*) (jt[j], j=0..dim-1)
762  // b += \sum_i dTb/dx_i(x*) x_i/(1-xn) (jt[dim-1])
763  // (b += 0 in case xn -> 1)
764  for( int j = 0; j < dim-1; ++j )
765  {
766  jt[ j ] += jt2[ j ];
767  jt[ dim-1 ].axpy( dfcxn*x[ j ], jt2[ j ] );
768  }
769  }
770  else
771  {
772  // jt = dTb/dx_i(x*)
773  jacobianTransposed< false >( topologyId, integral_constant< int, dim-1 >(), cit2, dfcxn, x, rf, jt );
774  // b += \sum_i dTb/dx_i(x*) x_i/(1-xn)
775  for( int j = 0; j < dim-1; ++j )
776  jt[ dim-1 ].axpy( dfcxn*x[ j ], jt[ j ] );
777  }
778  }
779  }
780 
781  template< class ct, int mydim, int cdim, class Traits >
782  template< bool add, int rows >
784  ::jacobianTransposed ( TopologyId topologyId, integral_constant< int, 0 >,
785  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
786  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt )
787  {
788  ++cit;
789  }
790 
791 
792 
793  template< class ct, int mydim, int cdim, class Traits >
794  template< int dim >
796  ::affine ( TopologyId topologyId, integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt )
797  {
798  const GlobalCoordinate &orgBottom = *cit;
799  if( !affine( topologyId, integral_constant< int, dim-1 >(), cit, jt ) )
800  return false;
801  const GlobalCoordinate &orgTop = *cit;
802 
803  if( GenericGeometry::isPrism( topologyId, mydimension, mydimension-dim ) )
804  {
805  JacobianTransposed jtTop;
806  if( !affine( topologyId, integral_constant< int, dim-1 >(), cit, jtTop ) )
807  return false;
808 
809  // check whether both jacobians are identical
810  ctype norm( 0 );
811  for( int i = 0; i < dim-1; ++i )
812  norm += (jtTop[ i ] - jt[ i ]).two_norm2();
813  if( norm >= Traits::tolerance() )
814  return false;
815  }
816  else
817  ++cit;
818  jt[ dim-1 ] = orgTop - orgBottom;
819  return true;
820  }
821 
822  template< class ct, int mydim, int cdim, class Traits >
824  ::affine ( TopologyId topologyId, integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt )
825  {
826  ++cit;
827  return true;
828  }
829 
830 } // namespace Dune
831 
832 #endif // #ifndef DUNE_GEOMETRY_MULTILINEARGEOMETRY_HH
Traits::MatrixHelper MatrixHelper
Definition: multilineargeometry.hh:190
CachedMultiLinearGeometry(const ReferenceElement &refElement, const CornerStorage &cornerStorage)
Definition: multilineargeometry.hh:461
Dune::GeometryType type() const
obtain the name of the reference element
Definition: multilineargeometry.hh:239
JacobianInverseTransposed jacobianInverseTransposed_
Definition: multilineargeometry.hh:386
static const bool v
Definition: multilineargeometry.hh:116
will there be only one geometry type for a dimension?
Definition: multilineargeometry.hh:114
TopologyId topologyId() const
Definition: multilineargeometry.hh:354
MultiLinearGeometry(const ReferenceElement &refElement, const Corners &corners)
constructor
Definition: multilineargeometry.hh:209
Definition: matrixhelper.hh:33
LocalCoordinate local(const GlobalCoordinate &global) const
evaluate the inverse mapping
Definition: multilineargeometry.hh:514
Implement a MultiLinearGeometry with additional caching.
Definition: multilineargeometry.hh:437
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
type of jacobian transposed
Definition: multilineargeometry.hh:173
Base::ctype ctype
Definition: multilineargeometry.hh:449
GlobalCoordinate center() const
obtain the centroid of the mapping's image
Definition: multilineargeometry.hh:252
bool affine() const
is this mapping affine?
Definition: multilineargeometry.hh:477
unsigned int topologyId(integral_constant< bool, false >) const
Definition: multilineargeometry.hh:383
bool isPrism(unsigned int topologyId, int dim, int codim=0)
check whether a prism construction was used to create a given codimension
Definition: topologytypes.hh:169
default traits class for MultiLinearGeometry
Definition: multilineargeometry.hh:46
int size(int c) const
number of subentities of codimension c
Definition: referenceelements.hh:86
Dune::ReferenceElement< ctype, mydimension > ReferenceElement
type of reference element
Definition: multilineargeometry.hh:184
A unique label for each type of element that can occur in a grid.
ctype volume() const
obtain the volume of the mapping's image
Definition: multilineargeometry.hh:559
ctype volume() const
obtain the volume of the reference element
Definition: referenceelements.hh:280
void setupDeterminant(const JacobianTransposed &jt)
Definition: multilineargeometry.hh:410
const ReferenceElement & refElement() const
Definition: multilineargeometry.hh:352
FieldVector< ctype, mydimension > LocalCoordinate
type of user data
Definition: multilineargeometry.hh:168
static const unsigned int topologyId
Definition: multilineargeometry.hh:117
This class provides access to geometric and topological properties of a reference element...
Definition: affinegeometry.hh:25
const FieldVector< ctype, dim > & position(int i, int c) const
position of the barycenter of entity (i,c)
Definition: referenceelements.hh:154
static const int coorddimension
coordinate dimension
Definition: multilineargeometry.hh:161
Base::JacobianTransposed JacobianTransposed
Definition: multilineargeometry.hh:457
static ct tolerance()
tolerance to numerical algorithms
Definition: multilineargeometry.hh:69
bool isPyramid(unsigned int topologyId, int dim, int codim=0)
check whether a pyramid construction was used to create a given codimension
Definition: topologytypes.hh:151
void setup(const JacobianTransposed &jt)
Definition: multilineargeometry.hh:405
ct ctype
coordinate type
Definition: multilineargeometry.hh:156
CachedMultiLinearGeometry(Dune::GeometryType gt, const CornerStorage &cornerStorage)
Definition: multilineargeometry.hh:469
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian
Definition: multilineargeometry.hh:336
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian's inverse
Definition: multilineargeometry.hh:627
Base::MatrixHelper MatrixHelper
Definition: multilineargeometry.hh:444
GlobalCoordinate corner(int i) const
obtain coordinates of the i-th corner
Definition: multilineargeometry.hh:245
GenericGeometry::MatrixHelper< GenericGeometry::DuneCoordTraits< ct > > MatrixHelper
helper structure containing some matrix routines
Definition: multilineargeometry.hh:66
unsigned int id() const
Return the topology id the type.
Definition: type.hh:327
Class providing access to the singletons of the reference elements. Special methods are available for...
Definition: affinegeometry.hh:28
GlobalCoordinate center() const
obtain the centroid of the mapping's image
Definition: multilineargeometry.hh:482
MultiLinearGeometry(Dune::GeometryType gt, const Corners &corners)
constructor
Definition: multilineargeometry.hh:225
Dune::ReferenceElements< ctype, mydimension > ReferenceElements
Definition: multilineargeometry.hh:193
GlobalCoordinate global(const LocalCoordinate &local) const
evaluate the mapping
Definition: multilineargeometry.hh:490
int corners() const
obtain number of corners of the corresponding reference element
Definition: multilineargeometry.hh:242
FieldVector< ctype, coorddimension > GlobalCoordinate
type of global coordinates
Definition: multilineargeometry.hh:170
ctype det() const
Definition: multilineargeometry.hh:415
TopologyId topologyId(integral_constant< bool, true >) const
Definition: multilineargeometry.hh:382
Base::ReferenceElement ReferenceElement
Definition: multilineargeometry.hh:447
Base::LocalCoordinate LocalCoordinate
Definition: multilineargeometry.hh:454
LocalCoordinate local(const GlobalCoordinate &global) const
evaluate the inverse mapping
Definition: multilineargeometry.hh:280
Definition: multilineargeometry.hh:399
ctype integrationElement(const LocalCoordinate &local) const
obtain the integration element
Definition: multilineargeometry.hh:543
Base::GlobalCoordinate GlobalCoordinate
Definition: multilineargeometry.hh:455
Base::JacobianInverseTransposed JacobianInverseTransposed
Definition: multilineargeometry.hh:458
GlobalCoordinate global(const LocalCoordinate &local) const
evaluate the mapping
Definition: multilineargeometry.hh:260
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:24
JacobianInverseTransposed Jacobian
For backward-compatibility, export the type JacobianInverseTransposed as Jacobian.
Definition: multilineargeometry.hh:176
ctype integrationElement(const LocalCoordinate &local) const
obtain the integration element
Definition: multilineargeometry.hh:309
ctype detInv() const
Definition: multilineargeometry.hh:416
const JacobianInverseTransposed & jacobianInverseTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian's inverse
Definition: multilineargeometry.hh:590
const GeometryType & type(int i, int c) const
obtain the type of subentity (i,c)
Definition: referenceelements.hh:136
template specifying the storage for the corners
Definition: multilineargeometry.hh:95
conditional< hasSingleGeometryType, integral_constant< unsigned int, Traits::template hasSingleGeometryType< mydimension >::topologyId >, unsigned int >::type TopologyId
Definition: multilineargeometry.hh:191
const JacobianTransposed & jacobianTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian
Definition: multilineargeometry.hh:576
generic geometry implementation based on corner coordinates
Definition: multilineargeometry.hh:150
bool affine() const
is this mapping affine?
Definition: multilineargeometry.hh:232
static const int mydimension
geometry dimension
Definition: multilineargeometry.hh:159
ctype volume() const
obtain the volume of the mapping's image
Definition: multilineargeometry.hh:322
std::vector< FieldVector< ct, cdim > > Type
Definition: multilineargeometry.hh:97
JacobianTransposed jacobianTransposed_
Definition: multilineargeometry.hh:385