dune-geometry  3.0-git
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 <functional>
8 #include <iterator>
9 #include <limits>
10 #include <vector>
11 
12 #include <dune/common/fmatrix.hh>
13 #include <dune/common/fvector.hh>
14 #include <dune/common/typetraits.hh>
15 
17 #include <dune/geometry/type.hh>
20 
21 namespace Dune
22 {
23 
24  // External Forward Declarations
25  // -----------------------------
26 
27  template< class ctype, int dim >
28  class ReferenceElement;
29 
30  template< class ctype, int dim >
31  struct ReferenceElements;
32 
33 
34 
35  // MultiLinearGeometryTraits
36  // -------------------------
37 
47  template< class ct >
49  {
69 
71  static ct tolerance () { return ct( 16 ) * std::numeric_limits< ct >::epsilon(); }
72 
137  template< int mydim, int cdim >
139  {
140  typedef std::vector< FieldVector< ct, cdim > > Type;
141  };
142 
156  template< int dim >
158  {
159  static const bool v = false;
160  static const unsigned int topologyId = ~0u;
161  };
162  };
163 
164 
165 
166  // MultiLinearGeometry
167  // -------------------
168 
189  template< class ct, int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct > >
191  {
193 
194  public:
196  typedef ct ctype;
197 
199  static const int mydimension= mydim;
201  static const int coorddimension = cdim;
202 
204  typedef FieldVector< ctype, mydimension > LocalCoordinate;
206  typedef FieldVector< ctype, coorddimension > GlobalCoordinate;
207 
209  typedef FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed;
210 
213 
216 
217  private:
218  static const bool hasSingleGeometryType = Traits::template hasSingleGeometryType< mydimension >::v;
219 
220  protected:
221  typedef typename Traits::MatrixHelper MatrixHelper;
222  typedef typename std::conditional< hasSingleGeometryType, std::integral_constant< unsigned int, Traits::template hasSingleGeometryType< mydimension >::topologyId >, unsigned int >::type TopologyId;
223 
225 
226  public:
236  template< class Corners >
238  const Corners &corners )
239  : refElement_( &refElement ),
240  corners_( corners )
241  {}
242 
252  template< class Corners >
254  const Corners &corners )
255  : refElement_( &ReferenceElements::general( gt ) ),
256  corners_( corners )
257  {}
258 
260  bool affine () const
261  {
263  return affine( jt );
264  }
265 
267  Dune::GeometryType type () const { return GeometryType( toUnsignedInt(topologyId()), mydimension ); }
268 
270  int corners () const { return refElement().size( mydimension ); }
271 
273  GlobalCoordinate corner ( int i ) const
274  {
275  assert( (i >= 0) && (i < corners()) );
276  return std::cref(corners_).get()[ i ];
277  }
278 
280  GlobalCoordinate center () const { return global( refElement().position( 0, 0 ) ); }
281 
288  GlobalCoordinate global ( const LocalCoordinate &local ) const
289  {
290  using std::begin;
291 
292  auto cit = begin(std::cref(corners_).get());
294  global< false >( topologyId(), std::integral_constant< int, mydimension >(), cit, ctype( 1 ), local, ctype( 1 ), y );
295  return y;
296  }
297 
310  LocalCoordinate local ( const GlobalCoordinate &globalCoord ) const
311  {
312  const ctype tolerance = Traits::tolerance();
313  LocalCoordinate x = refElement().position( 0, 0 );
314  LocalCoordinate dx;
315  do
316  {
317  // Newton's method: DF^n dx^n = F^n, x^{n+1} -= dx^n
318  const GlobalCoordinate dglobal = (*this).global( x ) - globalCoord;
319  MatrixHelper::template xTRightInvA< mydimension, coorddimension >( jacobianTransposed( x ), dglobal, dx );
320  x -= dx;
321  } while( dx.two_norm2() > tolerance );
322  return x;
323  }
324 
339  ctype integrationElement ( const LocalCoordinate &local ) const
340  {
341  return MatrixHelper::template sqrtDetAAT< mydimension, coorddimension >( jacobianTransposed( local ) );
342  }
343 
352  ctype volume () const
353  {
354  return integrationElement( refElement().position( 0, 0 ) ) * refElement().volume();
355  }
356 
367  {
368  using std::begin;
369 
371  auto cit = begin(std::cref(corners_).get());
372  jacobianTransposed< false >( topologyId(), std::integral_constant< int, mydimension >(), cit, ctype( 1 ), local, ctype( 1 ), jt );
373  return jt;
374  }
375 
382  JacobianInverseTransposed jacobianInverseTransposed ( const LocalCoordinate &local ) const;
383 
384  friend const ReferenceElement &referenceElement ( const MultiLinearGeometry &geometry ) { return geometry.refElement(); }
385 
386  protected:
387  const ReferenceElement &refElement () const { return *refElement_; }
388 
390  {
391  return topologyId( std::integral_constant< bool, hasSingleGeometryType >() );
392  }
393 
394  template< bool add, int dim, class CornerIterator >
395  static void global ( TopologyId topologyId, std::integral_constant< int, dim >,
396  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
397  const ctype &rf, GlobalCoordinate &y );
398  template< bool add, class CornerIterator >
399  static void global ( TopologyId topologyId, std::integral_constant< int, 0 >,
400  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
401  const ctype &rf, GlobalCoordinate &y );
402 
403  template< bool add, int rows, int dim, class CornerIterator >
404  static void jacobianTransposed ( TopologyId topologyId, std::integral_constant< int, dim >,
405  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
406  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt );
407  template< bool add, int rows, class CornerIterator >
408  static void jacobianTransposed ( TopologyId topologyId, std::integral_constant< int, 0 >,
409  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
410  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt );
411 
412  template< int dim, class CornerIterator >
413  static bool affine ( TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt );
414  template< class CornerIterator >
415  static bool affine ( TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt );
416 
417  bool affine ( JacobianTransposed &jacobianT ) const
418  {
419  using std::begin;
420 
421  auto cit = begin(std::cref(corners_).get());
422  return affine( topologyId(), std::integral_constant< int, mydimension >(), cit, jacobianT );
423  }
424 
425  private:
426  // The following methods are needed to convert the return type of topologyId to
427  // unsigned int with g++-4.4. It has problems casting integral_constant to the
428  // integral type.
429  static unsigned int toUnsignedInt(unsigned int i) { return i; }
430  template<unsigned int v>
431  static unsigned int toUnsignedInt(std::integral_constant<unsigned int,v> i) { return v; }
432  TopologyId topologyId ( std::integral_constant< bool, true > ) const { return TopologyId(); }
433  unsigned int topologyId ( std::integral_constant< bool, false > ) const { return refElement().type().id(); }
434 
435  const ReferenceElement *refElement_;
436  typename Traits::template CornerStorage< mydimension, coorddimension >::Type corners_;
437  };
438 
439 
440 
441  // MultiLinearGeometry::JacobianInverseTransposed
442  // ----------------------------------------------
443 
444  template< class ct, int mydim, int cdim, class Traits >
445  class MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed
446  : public FieldMatrix< ctype, coorddimension, mydimension >
447  {
448  typedef FieldMatrix< ctype, coorddimension, mydimension > Base;
449 
450  public:
451  void setup ( const JacobianTransposed &jt )
452  {
453  detInv_ = MatrixHelper::template rightInvA< mydimension, coorddimension >( jt, static_cast< Base & >( *this ) );
454  }
455 
457  {
458  detInv_ = MatrixHelper::template sqrtDetAAT< mydimension, coorddimension >( jt );
459  }
460 
461  ctype det () const { return ctype( 1 ) / detInv_; }
462  ctype detInv () const { return detInv_; }
463 
464  private:
465  ctype detInv_;
466  };
467 
468 
469 
482  template< class ct, int mydim, int cdim, class Traits = MultiLinearGeometryTraits< ct > >
484  : public MultiLinearGeometry< ct, mydim, cdim, Traits >
485  {
488 
489  protected:
491 
492  public:
494 
495  typedef typename Base::ctype ctype;
496 
497  using Base::mydimension;
498  using Base::coorddimension;
499 
502 
505 
506  template< class CornerStorage >
507  CachedMultiLinearGeometry ( const ReferenceElement &referenceElement, const CornerStorage &cornerStorage )
508  : Base( referenceElement, cornerStorage ),
509  affine_( Base::affine( jacobianTransposed_ ) ),
510  jacobianInverseTransposedComputed_( false ),
511  integrationElementComputed_( false )
512  {}
513 
514  template< class CornerStorage >
515  CachedMultiLinearGeometry ( Dune::GeometryType gt, const CornerStorage &cornerStorage )
516  : Base( gt, cornerStorage ),
517  affine_( Base::affine( jacobianTransposed_ ) ),
518  jacobianInverseTransposedComputed_( false ),
519  integrationElementComputed_( false )
520  {}
521 
523  bool affine () const { return affine_; }
524 
525  using Base::corner;
526 
528  GlobalCoordinate center () const { return global( refElement().position( 0, 0 ) ); }
529 
536  GlobalCoordinate global ( const LocalCoordinate &local ) const
537  {
538  if( affine() )
539  {
541  jacobianTransposed_.umtv( local, global );
542  return global;
543  }
544  else
545  return Base::global( local );
546  }
547 
560  LocalCoordinate local ( const GlobalCoordinate &global ) const
561  {
562  if( affine() )
563  {
564  LocalCoordinate local;
565  if( jacobianInverseTransposedComputed_ )
566  jacobianInverseTransposed_.mtv( global - corner( 0 ), local );
567  else
568  MatrixHelper::template xTRightInvA< mydimension, coorddimension >( jacobianTransposed_, global - corner( 0 ), local );
569  return local;
570  }
571  else
572  return Base::local( global );
573  }
574 
589  ctype integrationElement ( const LocalCoordinate &local ) const
590  {
591  if( affine() )
592  {
593  if( !integrationElementComputed_ )
594  {
595  jacobianInverseTransposed_.setupDeterminant( jacobianTransposed_ );
596  integrationElementComputed_ = true;
597  }
598  return jacobianInverseTransposed_.detInv();
599  }
600  else
601  return Base::integrationElement( local );
602  }
603 
605  ctype volume () const
606  {
607  if( affine() )
608  return integrationElement( refElement().position( 0, 0 ) ) * refElement().volume();
609  else
610  return Base::volume();
611  }
612 
623  {
624  if( affine() )
625  return jacobianTransposed_;
626  else
627  return Base::jacobianTransposed( local );
628  }
629 
637  {
638  if( affine() )
639  {
640  if( !jacobianInverseTransposedComputed_ )
641  {
642  jacobianInverseTransposed_.setup( jacobianTransposed_ );
643  jacobianInverseTransposedComputed_ = true;
644  integrationElementComputed_ = true;
645  }
646  return jacobianInverseTransposed_;
647  }
648  else
649  return Base::jacobianInverseTransposed( local );
650  }
651 
652  protected:
653  using Base::refElement;
654 
655  private:
656  mutable JacobianTransposed jacobianTransposed_;
657  mutable JacobianInverseTransposed jacobianInverseTransposed_;
658 
659  mutable bool affine_ : 1;
660 
661  mutable bool jacobianInverseTransposedComputed_ : 1;
662  mutable bool integrationElementComputed_ : 1;
663  };
664 
665 
666 
667  // Implementation of MultiLinearGeometry
668  // -------------------------------------
669 
670  template< class ct, int mydim, int cdim, class Traits >
671  inline typename MultiLinearGeometry< ct, mydim, cdim, Traits >::JacobianInverseTransposed
673  {
675  jit.setup( jacobianTransposed( local ) );
676  return jit;
677  }
678 
679 
680  template< class ct, int mydim, int cdim, class Traits >
681  template< bool add, int dim, class CornerIterator >
683  ::global ( TopologyId topologyId, std::integral_constant< int, dim >,
684  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
685  const ctype &rf, GlobalCoordinate &y )
686  {
687  const ctype xn = df*x[ dim-1 ];
688  const ctype cxn = ctype( 1 ) - xn;
689 
690  if( GenericGeometry::isPrism( toUnsignedInt(topologyId), mydimension, mydimension-dim ) )
691  {
692  // apply (1-xn) times mapping for bottom
693  global< add >( topologyId, std::integral_constant< int, dim-1 >(), cit, df, x, rf*cxn, y );
694  // apply xn times mapping for top
695  global< true >( topologyId, std::integral_constant< int, dim-1 >(), cit, df, x, rf*xn, y );
696  }
697  else
698  {
699  assert( GenericGeometry::isPyramid( toUnsignedInt(topologyId), mydimension, mydimension-dim ) );
700  // apply (1-xn) times mapping for bottom (with argument x/(1-xn))
701  if( cxn > Traits::tolerance() || cxn < -Traits::tolerance() )
702  global< add >( topologyId, std::integral_constant< int, dim-1 >(), cit, df/cxn, x, rf*cxn, y );
703  else
704  global< add >( topologyId, std::integral_constant< int, dim-1 >(), cit, df, x, ctype( 0 ), y );
705  // apply xn times the tip
706  y.axpy( rf*xn, *cit );
707  ++cit;
708  }
709  }
710 
711  template< class ct, int mydim, int cdim, class Traits >
712  template< bool add, class CornerIterator >
714  ::global ( TopologyId topologyId, std::integral_constant< int, 0 >,
715  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
716  const ctype &rf, GlobalCoordinate &y )
717  {
718  const GlobalCoordinate &origin = *cit;
719  ++cit;
720  for( int i = 0; i < coorddimension; ++i )
721  y[ i ] = (add ? y[ i ] + rf*origin[ i ] : rf*origin[ i ]);
722  }
723 
724 
725  template< class ct, int mydim, int cdim, class Traits >
726  template< bool add, int rows, int dim, class CornerIterator >
728  ::jacobianTransposed ( TopologyId topologyId, std::integral_constant< int, dim >,
729  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
730  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt )
731  {
732  assert( rows >= dim );
733 
734  const ctype xn = df*x[ dim-1 ];
735  const ctype cxn = ctype( 1 ) - xn;
736 
737  auto cit2( cit );
738  if( GenericGeometry::isPrism( toUnsignedInt(topologyId), mydimension, mydimension-dim ) )
739  {
740  // apply (1-xn) times Jacobian for bottom
741  jacobianTransposed< add >( topologyId, std::integral_constant< int, dim-1 >(), cit2, df, x, rf*cxn, jt );
742  // apply xn times Jacobian for top
743  jacobianTransposed< true >( topologyId, std::integral_constant< int, dim-1 >(), cit2, df, x, rf*xn, jt );
744  // compute last row as difference between top value and bottom value
745  global< add >( topologyId, std::integral_constant< int, dim-1 >(), cit, df, x, -rf, jt[ dim-1 ] );
746  global< true >( topologyId, std::integral_constant< int, dim-1 >(), cit, df, x, rf, jt[ dim-1 ] );
747  }
748  else
749  {
750  assert( GenericGeometry::isPyramid( toUnsignedInt(topologyId), mydimension, mydimension-dim ) );
751  /*
752  * In the pyramid case, we need a transformation Tb: B -> R^n for the
753  * base B \subset R^{n-1}. The pyramid transformation is then defined as
754  * T: P \subset R^n -> R^n
755  * (x, xn) |-> (1-xn) Tb(x*) + xn t (x \in R^{n-1}, xn \in R)
756  * with the tip of the pyramid mapped to t and x* = x/(1-xn)
757  * the projection of (x,xn) onto the base.
758  *
759  * For the Jacobi matrix DT we get
760  * DT = ( A | b )
761  * with A = DTb(x*) (n x n-1 matrix)
762  * and b = dT/dxn (n-dim column vector).
763  * Furthermore
764  * b = -Tb(x*) + t + \sum_i dTb/dx_i(x^*) x_i/(1-xn)
765  *
766  * Note that both A and b are not defined in the pyramid tip (x=0, xn=1)!
767  * Indeed for B the unit square, Tb mapping B to the quadrilateral given
768  * by the vertices (0,0,0), (2,0,0), (0,1,0), (1,1,0) and t=(0,0,1), we get
769  *
770  * T(x,y,xn) = ( x(2-y/(1-xn)), y, xn )
771  * / 2-y/(1-xn) -x 0 \
772  * DT(x,y,xn) = | 0 1 0 |
773  * \ 0 0 1 /
774  * which is not continuous for xn -> 1, choose for example
775  * x=0, y=1-xn, xn -> 1 --> DT -> diag(1,1,1)
776  * x=1-xn, y=0, xn -> 1 --> DT -> diag(2,1,1)
777  *
778  * However, for Tb affine-linear, Tb(y) = My + y0, DTb = M:
779  * A = M
780  * b = -M x* - y0 + t + \sum_i M_i x_i/(1-xn)
781  * = -M x* - y0 + t + M x*
782  * = -y0 + t
783  * which is continuous for xn -> 1. Note that this b is also given by
784  * b = -Tb(0) + t + \sum_i dTb/dx_i(0) x_i/1
785  * that is replacing x* by 1 and 1-xn by 1 in the formular above.
786  *
787  * For xn -> 1, we can thus set x*=0, "1-xn"=1 (or anything != 0) and get
788  * the right result in case Tb is affine-linear.
789  */
790 
791  /* The second case effectively results in x* = 0 */
792  ctype dfcxn = (cxn > Traits::tolerance() || cxn < -Traits::tolerance()) ? ctype(df / cxn) : ctype(0);
793 
794  // initialize last row
795  // b = -Tb(x*)
796  // (b = -Tb(0) = -y0 in case xn -> 1 and Tb affine-linear)
797  global< add >( topologyId, std::integral_constant< int, dim-1 >(), cit, dfcxn, x, -rf, jt[ dim-1 ] );
798  // b += t
799  jt[ dim-1 ].axpy( rf, *cit );
800  ++cit;
801  // apply Jacobian for bottom (with argument x/(1-xn)) and correct last row
802  if( add )
803  {
804  FieldMatrix< ctype, dim-1, coorddimension > jt2;
805  // jt2 = dTb/dx_i(x*)
806  jacobianTransposed< false >( topologyId, std::integral_constant< int, dim-1 >(), cit2, dfcxn, x, rf, jt2 );
807  // A = dTb/dx_i(x*) (jt[j], j=0..dim-1)
808  // b += \sum_i dTb/dx_i(x*) x_i/(1-xn) (jt[dim-1])
809  // (b += 0 in case xn -> 1)
810  for( int j = 0; j < dim-1; ++j )
811  {
812  jt[ j ] += jt2[ j ];
813  jt[ dim-1 ].axpy( dfcxn*x[ j ], jt2[ j ] );
814  }
815  }
816  else
817  {
818  // jt = dTb/dx_i(x*)
819  jacobianTransposed< false >( topologyId, std::integral_constant< int, dim-1 >(), cit2, dfcxn, x, rf, jt );
820  // b += \sum_i dTb/dx_i(x*) x_i/(1-xn)
821  for( int j = 0; j < dim-1; ++j )
822  jt[ dim-1 ].axpy( dfcxn*x[ j ], jt[ j ] );
823  }
824  }
825  }
826 
827  template< class ct, int mydim, int cdim, class Traits >
828  template< bool add, int rows, class CornerIterator >
830  ::jacobianTransposed ( TopologyId topologyId, std::integral_constant< int, 0 >,
831  CornerIterator &cit, const ctype &df, const LocalCoordinate &x,
832  const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt )
833  {
834  ++cit;
835  }
836 
837 
838 
839  template< class ct, int mydim, int cdim, class Traits >
840  template< int dim, class CornerIterator >
842  ::affine ( TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt )
843  {
844  const GlobalCoordinate &orgBottom = *cit;
845  if( !affine( topologyId, std::integral_constant< int, dim-1 >(), cit, jt ) )
846  return false;
847  const GlobalCoordinate &orgTop = *cit;
848 
849  if( GenericGeometry::isPrism( toUnsignedInt(topologyId), mydimension, mydimension-dim ) )
850  {
851  JacobianTransposed jtTop;
852  if( !affine( topologyId, std::integral_constant< int, dim-1 >(), cit, jtTop ) )
853  return false;
854 
855  // check whether both jacobians are identical
856  ctype norm( 0 );
857  for( int i = 0; i < dim-1; ++i )
858  norm += (jtTop[ i ] - jt[ i ]).two_norm2();
859  if( norm >= Traits::tolerance() )
860  return false;
861  }
862  else
863  ++cit;
864  jt[ dim-1 ] = orgTop - orgBottom;
865  return true;
866  }
867 
868  template< class ct, int mydim, int cdim, class Traits >
869  template< class CornerIterator >
871  ::affine ( TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt )
872  {
873  ++cit;
874  return true;
875  }
876 
877 } // namespace Dune
878 
879 #endif // #ifndef DUNE_GEOMETRY_MULTILINEARGEOMETRY_HH
A unique label for each type of element that can occur in a grid.
Definition: affinegeometry.hh:19
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:150
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:168
This class provides access to geometric and topological properties of a reference element.
Definition: referenceelements.hh:29
int size(int c) const
number of subentities of codimension c
Definition: referenceelements.hh:80
ctype volume() const
obtain the volume of the reference element
Definition: referenceelements.hh:186
const FieldVector< ctype, dim > & position(int i, int c) const
position of the barycenter of entity (i,c)
Definition: referenceelements.hh:148
const GeometryType & type(int i, int c) const
obtain the type of subentity (i,c)
Definition: referenceelements.hh:130
Class providing access to the singletons of the reference elements.
Definition: referenceelements.hh:455
Definition: matrixhelper.hh:36
default traits class for MultiLinearGeometry
Definition: multilineargeometry.hh:49
static ct tolerance()
tolerance to numerical algorithms
Definition: multilineargeometry.hh:71
GenericGeometry::MatrixHelper< GenericGeometry::DuneCoordTraits< ct > > MatrixHelper
helper structure containing some matrix routines
Definition: multilineargeometry.hh:68
template specifying the storage for the corners
Definition: multilineargeometry.hh:139
std::vector< FieldVector< ct, cdim > > Type
Definition: multilineargeometry.hh:140
will there be only one geometry type for a dimension?
Definition: multilineargeometry.hh:158
static const unsigned int topologyId
Definition: multilineargeometry.hh:160
static const bool v
Definition: multilineargeometry.hh:159
generic geometry implementation based on corner coordinates
Definition: multilineargeometry.hh:191
static void global(TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, GlobalCoordinate &y)
Definition: multilineargeometry.hh:683
static const int mydimension
geometry dimension
Definition: multilineargeometry.hh:199
Dune::ReferenceElement< ctype, mydimension > ReferenceElement
type of reference element
Definition: multilineargeometry.hh:212
Dune::GeometryType type() const
obtain the name of the reference element
Definition: multilineargeometry.hh:267
Traits::MatrixHelper MatrixHelper
Definition: multilineargeometry.hh:221
FieldVector< ctype, coorddimension > GlobalCoordinate
type of global coordinates
Definition: multilineargeometry.hh:206
static void jacobianTransposed(TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt)
Definition: multilineargeometry.hh:830
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian
Definition: multilineargeometry.hh:366
const ReferenceElement & refElement() const
Definition: multilineargeometry.hh:387
GlobalCoordinate global(const LocalCoordinate &local) const
evaluate the mapping
Definition: multilineargeometry.hh:288
GlobalCoordinate center() const
obtain the centroid of the mapping's image
Definition: multilineargeometry.hh:280
GlobalCoordinate corner(int i) const
obtain coordinates of the i-th corner
Definition: multilineargeometry.hh:273
Dune::ReferenceElements< ctype, mydimension > ReferenceElements
Definition: multilineargeometry.hh:224
ct ctype
coordinate type
Definition: multilineargeometry.hh:196
static void jacobianTransposed(TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, FieldMatrix< ctype, rows, cdim > &jt)
Definition: multilineargeometry.hh:728
static const int coorddimension
coordinate dimension
Definition: multilineargeometry.hh:201
int corners() const
obtain number of corners of the corresponding reference element
Definition: multilineargeometry.hh:270
TopologyId topologyId() const
Definition: multilineargeometry.hh:389
std::conditional< hasSingleGeometryType, std::integral_constant< unsigned int, Traits::template hasSingleGeometryType< mydimension >::topologyId >, unsigned int >::type TopologyId
Definition: multilineargeometry.hh:222
static void global(TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, const ctype &df, const LocalCoordinate &x, const ctype &rf, GlobalCoordinate &y)
Definition: multilineargeometry.hh:714
FieldVector< ctype, mydimension > LocalCoordinate
type of local coordinates
Definition: multilineargeometry.hh:204
ctype volume() const
obtain the volume of the mapping's image
Definition: multilineargeometry.hh:352
MultiLinearGeometry(const ReferenceElement &refElement, const Corners &corners)
constructor
Definition: multilineargeometry.hh:237
static bool affine(TopologyId topologyId, std::integral_constant< int, 0 >, CornerIterator &cit, JacobianTransposed &jt)
Definition: multilineargeometry.hh:871
static bool affine(TopologyId topologyId, std::integral_constant< int, dim >, CornerIterator &cit, JacobianTransposed &jt)
Definition: multilineargeometry.hh:842
MultiLinearGeometry(Dune::GeometryType gt, const Corners &corners)
constructor
Definition: multilineargeometry.hh:253
ctype integrationElement(const LocalCoordinate &local) const
obtain the integration element
Definition: multilineargeometry.hh:339
FieldMatrix< ctype, mydimension, coorddimension > JacobianTransposed
type of jacobian transposed
Definition: multilineargeometry.hh:209
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian's inverse
Definition: multilineargeometry.hh:672
friend const ReferenceElement & referenceElement(const MultiLinearGeometry &geometry)
Definition: multilineargeometry.hh:384
bool affine() const
is this mapping affine?
Definition: multilineargeometry.hh:260
bool affine(JacobianTransposed &jacobianT) const
Definition: multilineargeometry.hh:417
Definition: multilineargeometry.hh:447
void setup(const JacobianTransposed &jt)
Definition: multilineargeometry.hh:451
ctype det() const
Definition: multilineargeometry.hh:461
ctype detInv() const
Definition: multilineargeometry.hh:462
void setupDeterminant(const JacobianTransposed &jt)
Definition: multilineargeometry.hh:456
Implement a MultiLinearGeometry with additional caching.
Definition: multilineargeometry.hh:485
GlobalCoordinate global(const LocalCoordinate &local) const
evaluate the mapping
Definition: multilineargeometry.hh:536
Base::ReferenceElement ReferenceElement
Definition: multilineargeometry.hh:493
bool affine() const
is this mapping affine?
Definition: multilineargeometry.hh:523
CachedMultiLinearGeometry(const ReferenceElement &referenceElement, const CornerStorage &cornerStorage)
Definition: multilineargeometry.hh:507
const ReferenceElement & refElement() const
Definition: multilineargeometry.hh:387
Base::MatrixHelper MatrixHelper
Definition: multilineargeometry.hh:490
Base::LocalCoordinate LocalCoordinate
Definition: multilineargeometry.hh:500
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian
Definition: multilineargeometry.hh:622
GlobalCoordinate corner(int i) const
obtain coordinates of the i-th corner
Definition: multilineargeometry.hh:273
CachedMultiLinearGeometry(Dune::GeometryType gt, const CornerStorage &cornerStorage)
Definition: multilineargeometry.hh:515
ctype volume() const
obtain the volume of the mapping's image
Definition: multilineargeometry.hh:605
ctype integrationElement(const LocalCoordinate &local) const
obtain the integration element
Definition: multilineargeometry.hh:589
Base::ctype ctype
Definition: multilineargeometry.hh:495
Base::JacobianInverseTransposed JacobianInverseTransposed
Definition: multilineargeometry.hh:504
Base::JacobianTransposed JacobianTransposed
Definition: multilineargeometry.hh:503
Base::GlobalCoordinate GlobalCoordinate
Definition: multilineargeometry.hh:501
GlobalCoordinate center() const
obtain the centroid of the mapping's image
Definition: multilineargeometry.hh:528
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
obtain the transposed of the Jacobian's inverse
Definition: multilineargeometry.hh:636
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:25
unsigned int id() const
Return the topology id the type.
Definition: type.hh:326