dune-grid  2.9.0
common/geometry.hh
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
2 // SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4 // vi: set et ts=4 sw=2 sts=2:
5 #ifndef DUNE_GRID_COMMON_GEOMETRY_HH
6 #define DUNE_GRID_COMMON_GEOMETRY_HH
7 
12 #include <cassert>
13 
14 #include <dune/common/fmatrix.hh>
15 #include <dune/common/typetraits.hh>
16 #include <dune/common/transpose.hh>
17 #include <dune/common/std/type_traits.hh>
18 
19 #include <dune/geometry/referenceelements.hh>
20 
21 namespace Dune
22 {
23 
24  // External Forward Declarations
25  // -----------------------------
26 
27  template< int dim, int dimworld, class ct, class GridFamily >
28  class GridDefaultImplementation;
29 
30 
31 
32  //*****************************************************************************
33  //
34  // Geometry
35  // forwards the interface to the implementation
36  //
37  //*****************************************************************************
38 
69  template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp >
70  class Geometry
71  {
72  public:
78  typedef GeometryImp< mydim, cdim, GridImp > Implementation;
79 
91  const Implementation &impl () const { return realGeometry; }
92 
94  constexpr static int mydimension = mydim;
95 
97  constexpr static int coorddimension = cdim;
98 
100  typedef typename GridImp::ctype ctype;
101 
103  typedef FieldVector<ctype, mydim> LocalCoordinate;
104 
106  typedef FieldVector< ctype, cdim > GlobalCoordinate;
107 
109  typedef decltype(std::declval<Implementation>().volume()) Volume;
110 
121 
132 
133  private:
134 
135  template<class Implementation_T>
136  using JacobianInverseOfImplementation = decltype(typename Implementation_T::JacobianInverse{std::declval<Implementation_T>().jacobianInverse(std::declval<LocalCoordinate>())});
137 
138  using JacobianInverseDefault = decltype(transpose(std::declval<JacobianInverseTransposed>()));
139 
140  template<class Implementation_T>
141  using JacobianOfImplementation = decltype(typename Implementation_T::Jacobian{std::declval<Implementation_T>().jacobian(std::declval<LocalCoordinate>())});
142 
143  using JacobianDefault = decltype(transpose(std::declval<JacobianTransposed>()));
144 
145  auto jacobianImpl ( const LocalCoordinate &local, std::true_type /*implDetected*/ ) const
146  {
147  return impl().jacobian(local);
148  }
149 
150  [[deprecated("Geometry implementatons are required to provide a jacobian(local) method. The default implementation is deprecated and will be removed after release 2.9")]]
151  auto jacobianImpl ( const LocalCoordinate &local, std::false_type /*implNotDetected*/ ) const
152  {
153  return transpose(jacobianTransposed(local));
154  }
155 
156  auto jacobianInverseImpl ( const LocalCoordinate &local, std::true_type /*implDetected*/ ) const
157  {
158  return impl().jacobianInverse(local);
159  }
160 
161  [[deprecated("Geometry implementatons are required to provide a jacobianInverse(local) method. The default implementation is deprecated and will be removed after release 2.9")]]
162  auto jacobianInverseImpl ( const LocalCoordinate &local, std::false_type /*implNotDetected*/ ) const
163  {
164  return transpose(jacobianInverseTransposed(local));
165  }
166 
167  public:
168 
178  using JacobianInverse = Std::detected_or_t<JacobianInverseDefault, JacobianInverseOfImplementation, Implementation>;
179 
189  using Jacobian = Std::detected_or_t<JacobianDefault, JacobianOfImplementation, Implementation>;
190 
194  GeometryType type () const { return impl().type(); }
195 
197  bool affine() const { return impl().affine(); }
198 
205  int corners () const { return impl().corners(); }
206 
219  GlobalCoordinate corner ( int i ) const
220  {
221  return impl().corner( i );
222  }
223 
229  {
230  return impl().global( local );
231  }
232 
238  {
239  return impl().local( global );
240  }
241 
266  {
267  return impl().integrationElement( local );
268  }
269 
271  Volume volume () const
272  {
273  return impl().volume();
274  }
275 
287  {
288  return impl().center();
289  }
290 
303  {
304  return impl().jacobianTransposed( local );
305  }
306 
329  {
330  return impl().jacobianInverseTransposed(local);
331  }
332 
344  Jacobian jacobian ( const LocalCoordinate& local ) const
345  {
346  auto implDetected = Std::is_detected<JacobianOfImplementation, Implementation>{};
347  return jacobianImpl(local, implDetected);
348  }
349 
372  {
373  auto implDetected = Std::is_detected<JacobianInverseOfImplementation, Implementation>{};
374  return jacobianInverseImpl(local, implDetected);
375  }
376 
377  //===========================================================
381  //===========================================================
382 
384  explicit Geometry ( const Implementation &impl )
385  : realGeometry( impl )
386  {}
387 
389 
390  protected:
391 
393  };
394 
395 
396 
397  //************************************************************************
398  // GEOMETRY Default Implementations
399  //*************************************************************************
400  //
401  // --GeometryDefault
402  //
404  template<int mydim, int cdim, class GridImp, template<int,int,class> class GeometryImp>
406  {
407  public:
408  static const int mydimension = mydim;
409  static const int coorddimension = cdim;
410 
411  // save typing
412  typedef typename GridImp::ctype ctype;
413 
414  typedef FieldVector< ctype, mydim > LocalCoordinate;
415  typedef FieldVector< ctype, cdim > GlobalCoordinate;
416 
418  typedef ctype Volume;
419 
421  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
422 
424  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
425 
427  typedef FieldMatrix< ctype, mydim, cdim > JacobianInverse;
428 
430  typedef FieldMatrix< ctype, cdim, mydim > Jacobian;
431 
433  Volume volume () const
434  {
435  GeometryType type = asImp().type();
436 
437  // get corresponding reference element
438  auto refElement = referenceElement< ctype , mydim >(type);
439 
440  LocalCoordinate localBaryCenter ( 0 );
441  // calculate local bary center
442  const int corners = refElement.size(0,0,mydim);
443  for(int i=0; i<corners; ++i) localBaryCenter += refElement.position(i,mydim);
444  localBaryCenter *= (ctype) (1.0/corners);
445 
446  // volume is volume of reference element times integrationElement
447  return refElement.volume() * asImp().integrationElement(localBaryCenter);
448  }
449 
452  {
453  GeometryType type = asImp().type();
454 
455  // get corresponding reference element
456  auto refElement = referenceElement< ctype , mydim >(type);
457 
458  // center is (for now) the centroid of the reference element mapped to
459  // this geometry.
460  return asImp().global(refElement.position(0,0));
461  }
462 
464  Jacobian jacobian ( const LocalCoordinate& local ) const
465  {
466  return asImp().jacobianTransposed(local).transposed();
467  }
468 
471  {
472  return asImp().jacobianInverseTransposed(local).transposed();
473  }
474 
475  private:
477  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
478  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
479  }; // end GeometryDefault
480 
481  template<int cdim, class GridImp, template<int,int,class> class GeometryImp>
482  class GeometryDefaultImplementation<0,cdim,GridImp,GeometryImp>
483  {
484  // my dimension
485  constexpr static int mydim = 0;
486 
487  public:
488  static const int mydimension = mydim;
489  static const int coorddimension = cdim;
490 
491  // save typing
492  typedef typename GridImp::ctype ctype;
493 
494  typedef FieldVector< ctype, mydim > LocalCoordinate;
495  typedef FieldVector< ctype, cdim > GlobalCoordinate;
496  typedef ctype Volume;
497 
499  typedef FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed;
500 
502  typedef FieldMatrix< ctype, mydim, cdim > JacobianTransposed;
503 
505  typedef FieldMatrix< ctype, mydim, cdim > JacobianInverse;
506 
508  typedef FieldMatrix< ctype, cdim, mydim > Jacobian;
509 
511  FieldVector<ctype, cdim> global (const FieldVector<ctype, mydim>& local) const
512  {
513  return asImp().corner(0);
514  }
515 
517  FieldVector<ctype, mydim> local (const FieldVector<ctype, cdim>& ) const
518  {
519  return FieldVector<ctype, mydim>();
520  }
521 
523  Volume volume () const
524  {
525  return Volume(1.0);
526  }
527 
529  FieldVector<ctype, cdim> center () const
530  {
531  return asImp().corner(0);
532  }
533 
535  Jacobian jacobian ( const LocalCoordinate& local ) const
536  {
537  return asImp().jacobianTransposed(local).transposed();
538  }
539 
542  {
543  return asImp().jacobianInverseTransposed(local).transposed();
544  }
545 
546  private:
547  // Barton-Nackman trick
548  GeometryImp<mydim,cdim,GridImp>& asImp () {return static_cast<GeometryImp<mydim,cdim,GridImp>&>(*this);}
549  const GeometryImp<mydim,cdim,GridImp>& asImp () const {return static_cast<const GeometryImp<mydim,cdim,GridImp>&>(*this);}
550  }; // end GeometryDefault
551 
552 
553 
554  // referenceElement
555  // ----------------
556 
557  template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp>
559  -> decltype(referenceElement(geo,geo.impl()))
560  {
561  return referenceElement(geo,geo.impl());
562  }
563 
565 
580  template< int mydim, int cdim, class GridImp, template< int, int, class > class GeometryImp, typename Impl>
582  -> decltype(referenceElement<typename GridImp::ctype,mydim>(geo.type()))
583  {
585  return referenceElement<typename Geo::ctype,Geo::mydimension>(geo.type());
586  }
587 
588 } // namespace Dune
589 
590 #endif // DUNE_GRID_COMMON_GEOMETRY_HH
Include standard header files.
Definition: agrid.hh:60
auto referenceElement(const Geometry< mydim, cdim, GridImp, GeometryImp > &geo) -> decltype(referenceElement(geo, geo.impl()))
Definition: common/geometry.hh:558
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
Wrapper class for geometries.
Definition: common/geometry.hh:71
Volume integrationElement(const LocalCoordinate &local) const
Return the factor appearing in the integral transformation formula.
Definition: common/geometry.hh:265
Jacobian jacobian(const LocalCoordinate &local) const
Return the Jacobian.
Definition: common/geometry.hh:344
const Implementation & impl() const
access to the underlying implementation
Definition: common/geometry.hh:91
FieldVector< ctype, cdim > GlobalCoordinate
type of the global coordinates
Definition: common/geometry.hh:106
JacobianTransposed jacobianTransposed(const LocalCoordinate &local) const
Return the transposed of the Jacobian.
Definition: common/geometry.hh:302
GeometryType type() const
Return the type of the reference element. The type can be used to access the Dune::ReferenceElement.
Definition: common/geometry.hh:194
GlobalCoordinate global(const LocalCoordinate &local) const
Evaluate the map .
Definition: common/geometry.hh:228
JacobianInverse jacobianInverse(const LocalCoordinate &local) const
Return inverse of Jacobian.
Definition: common/geometry.hh:371
Geometry(const Implementation &impl)
copy constructor from implementation
Definition: common/geometry.hh:384
Implementation::JacobianTransposed JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:131
Implementation::JacobianInverseTransposed JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:120
GlobalCoordinate corner(int i) const
Obtain a corner of the geometry.
Definition: common/geometry.hh:219
decltype(std::declval< Implementation >().volume()) typedef Volume
Number type used for the geometry volume.
Definition: common/geometry.hh:109
int corners() const
Return the number of corners of the reference element.
Definition: common/geometry.hh:205
Volume volume() const
return volume of geometry
Definition: common/geometry.hh:271
Std::detected_or_t< JacobianInverseDefault, JacobianInverseOfImplementation, Implementation > JacobianInverse
type of jacobian inverse
Definition: common/geometry.hh:178
GridImp::ctype ctype
define type used for coordinates in grid module
Definition: common/geometry.hh:100
constexpr static int coorddimension
dimension of embedding coordinate system
Definition: common/geometry.hh:97
Implementation & impl()
access to the underlying implementation
Definition: common/geometry.hh:85
FieldVector< ctype, mydim > LocalCoordinate
type of local coordinates
Definition: common/geometry.hh:103
Implementation realGeometry
Definition: common/geometry.hh:392
GlobalCoordinate center() const
return center of geometry
Definition: common/geometry.hh:286
constexpr static int mydimension
geometry dimension
Definition: common/geometry.hh:94
bool affine() const
Return true if the geometry mapping is affine and false otherwise.
Definition: common/geometry.hh:197
auto referenceElement(const Geometry< mydim, cdim, GridImp, GeometryImp > &geo, const Impl &) -> decltype(referenceElement< typename GridImp::ctype, mydim >(geo.type()))
Second-level dispatch to select the correct reference element for a grid geometry.
Definition: common/geometry.hh:581
Std::detected_or_t< JacobianDefault, JacobianOfImplementation, Implementation > Jacobian
type of jacobian
Definition: common/geometry.hh:189
GeometryImp< mydim, cdim, GridImp > Implementation
type of underlying implementation
Definition: common/geometry.hh:78
JacobianInverseTransposed jacobianInverseTransposed(const LocalCoordinate &local) const
Return inverse of transposed of Jacobian.
Definition: common/geometry.hh:328
Default implementation for class Geometry.
Definition: common/geometry.hh:406
FieldVector< ctype, cdim > GlobalCoordinate
Definition: common/geometry.hh:415
static const int mydimension
Definition: common/geometry.hh:408
Volume volume() const
return volume of the geometry
Definition: common/geometry.hh:433
FieldVector< ctype, mydim > LocalCoordinate
Definition: common/geometry.hh:414
GlobalCoordinate center() const
return center of the geometry
Definition: common/geometry.hh:451
FieldMatrix< ctype, cdim, mydim > Jacobian
type of jacobian
Definition: common/geometry.hh:430
JacobianInverse jacobianInverse(const LocalCoordinate &local) const
Return inverse of Jacobian.
Definition: common/geometry.hh:470
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:424
Jacobian jacobian(const LocalCoordinate &local) const
Return the Jacobian.
Definition: common/geometry.hh:464
GridImp::ctype ctype
Definition: common/geometry.hh:412
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:421
static const int coorddimension
Definition: common/geometry.hh:409
FieldMatrix< ctype, mydim, cdim > JacobianInverse
type of jacobian inverse
Definition: common/geometry.hh:427
ctype Volume
Number type used for the geometry volume.
Definition: common/geometry.hh:418
GridImp::ctype ctype
Definition: common/geometry.hh:492
FieldVector< ctype, cdim > GlobalCoordinate
Definition: common/geometry.hh:495
FieldMatrix< ctype, cdim, mydim > JacobianInverseTransposed
type of jacobian inverse transposed
Definition: common/geometry.hh:499
FieldMatrix< ctype, mydim, cdim > JacobianInverse
type of jacobian inverse
Definition: common/geometry.hh:505
FieldMatrix< ctype, cdim, mydim > Jacobian
type of jacobian
Definition: common/geometry.hh:508
Jacobian jacobian(const LocalCoordinate &local) const
Return the Jacobian.
Definition: common/geometry.hh:535
JacobianInverse jacobianInverse(const LocalCoordinate &local) const
Return inverse of Jacobian.
Definition: common/geometry.hh:541
FieldVector< ctype, mydim > LocalCoordinate
Definition: common/geometry.hh:494
FieldMatrix< ctype, mydim, cdim > JacobianTransposed
type of jacobian transposed
Definition: common/geometry.hh:502
FieldVector< ctype, cdim > center() const
return center of the geometry
Definition: common/geometry.hh:529
Volume volume() const
return volume of the geometry
Definition: common/geometry.hh:523