12#ifndef DUMUX_DISCRETIZATION_FACECENTERED_DIAMOND_GEOMETRY_HELPER_HH 
   13#define DUMUX_DISCRETIZATION_FACECENTERED_DIAMOND_GEOMETRY_HELPER_HH 
   17#include <dune/common/reservedvector.hh> 
   18#include <dune/common/fvector.hh> 
   19#include <dune/geometry/multilineargeometry.hh> 
   20#include <dune/geometry/type.hh> 
   32    template< 
int mydim, 
int cdim >
 
   35        using Type = Dune::ReservedVector< Dune::FieldVector< ct, cdim >, (1<<mydim)>;
 
 
 
   41template<Dune::GeometryType::Id gt>
 
   47    static constexpr Dune::GeometryType 
type()
 
   48    { 
return Dune::GeometryTypes::triangle; }
 
 
   50    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
   51    static constexpr std::array<std::array<Key, 3>, 3> 
keys = {{
 
 
 
   61    static constexpr Dune::GeometryType 
type()
 
   62    { 
return Dune::GeometryTypes::triangle; }
 
 
   64    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
   65    static constexpr std::array<std::array<Key, 3>, 4> 
keys = {{
 
 
 
   76    static constexpr Dune::GeometryType 
type()
 
   77    { 
return Dune::GeometryTypes::tetrahedron; }
 
 
   79    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
   80    static constexpr std::array<std::array<Key, 4>, 4> 
keys = {{
 
 
 
   91    static constexpr Dune::GeometryType 
type()
 
   92    { 
return Dune::GeometryTypes::pyramid; }
 
 
   94    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
   95    static constexpr std::array<std::array<Key, 5>, 6> 
keys = {{
 
   96        { 
Key{4, 3}, 
Key{0, 3}, 
Key{6, 3}, 
Key{2, 3}, 
Key{0, 0} },
 
   97        { 
Key{1, 3}, 
Key{5, 3}, 
Key{3, 3}, 
Key{7, 3}, 
Key{0, 0} },
 
   98        { 
Key{4, 3}, 
Key{5, 3}, 
Key{0, 3}, 
Key{1, 3}, 
Key{0, 0} },
 
   99        { 
Key{2, 3}, 
Key{3, 3}, 
Key{6, 3}, 
Key{7, 3}, 
Key{0, 0} },
 
  100        { 
Key{0, 3}, 
Key{1, 3}, 
Key{2, 3}, 
Key{3, 3}, 
Key{0, 0} },
 
 
 
  105template<Dune::GeometryType::Id gt>
 
  111    static constexpr Dune::GeometryType 
type()
 
  112    { 
return Dune::GeometryTypes::line; }
 
 
  114    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
  115    static constexpr std::array<std::array<Key, 2>, 3> 
keys = {{
 
 
 
  125    static constexpr Dune::GeometryType 
type()
 
  126    { 
return Dune::GeometryTypes::line; }
 
 
  128    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
  129    static constexpr std::array<std::array<Key, 2>, 4> 
keys = {{
 
 
 
  140    static constexpr Dune::GeometryType 
type()
 
  141    { 
return Dune::GeometryTypes::triangle; }
 
 
  143    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
  144    static constexpr std::array<std::array<Key, 3>, 6> 
keys = {{
 
 
 
  157    static constexpr Dune::GeometryType 
type()
 
  158    { 
return Dune::GeometryTypes::triangle; }
 
 
  160    using Key = std::pair<std::uint8_t, std::uint8_t>; 
 
  161    static constexpr std::array<std::array<Key, 3>, 12> 
keys = {{
 
 
 
  178template<
class S, 
class Geo, 
class KeyArray, std::size_t... I>
 
  181    using Dune::referenceElement;
 
  182    const auto ref = referenceElement(geo);
 
  184    return { geo.global(ref.position(key[I].first, key[I].second))... };
 
 
  188template<
class S, 
class Geo, 
class T, std::
size_t N, 
class Indices = std::make_index_sequence<N>>
 
  195template<
class S, 
class Geo, std::size_t... ii>
 
  198    using Dune::referenceElement;
 
  199    const auto ref = referenceElement(geo);
 
  201    return { geo.global(ref.position(ref.subEntity(i, 1, ii, Geo::mydimension), Geo::mydimension))... };
 
 
  205template<
class S, std::
size_t numCorners, 
class Geo>
 
  211template<
class IndexType, Dune::GeometryType::Id gt>
 
  214template<
class IndexType>
 
  217    static constexpr std::array<std::array<IndexType, 2>, 3> 
pairs = {{
 
  218        {0, 1}, {0, 2}, {1, 2}
 
 
 
  222template<
class IndexType>
 
  225    static constexpr std::array<std::array<IndexType, 2>, 4> 
pairs = {{
 
  226        {0, 2}, {1, 2}, {0, 3}, {1, 3}
 
 
 
  230template<
class IndexType>
 
  233    static constexpr std::array<std::array<IndexType, 2>, 6> 
pairs = {{
 
  234        {0, 1}, {0, 2}, {0, 3}, {1, 2}, {1, 3}, {2, 3},
 
 
 
  238template<
class IndexType>
 
  241    static constexpr std::array<std::array<IndexType, 2>, 12> 
pairs = {{
 
  242        {0, 2}, {1, 2}, {0, 3}, {1, 3}, {0, 4}, {1, 4},
 
  243        {2, 4}, {3, 4}, {0, 5}, {1, 5}, {2, 5}, {3, 5}
 
 
 
 
  253template <
class Gr
idView, 
class ScvType, 
class ScvfType>
 
  256    using Element = 
typename GridView::template Codim<0>::Entity;
 
  258    static constexpr auto dim = GridView::dimension;
 
  259    static constexpr auto dimWorld = GridView::dimensionworld;
 
  261    using ScvCornerStorage = 
typename ScvType::Traits::CornerStorage;
 
  262    using LocalIndexType = 
typename ScvType::Traits::LocalIndexType;
 
  263    using ScvfCornerStorage = 
typename ScvfType::Traits::CornerStorage;
 
  266    using ctype = 
typename GridView::ctype;
 
  267    using GlobalPosition = 
typename Dune::FieldVector<ctype, GridView::dimensionworld>;
 
  277        const auto type = geo_.type();
 
  278        if (type == Dune::GeometryTypes::triangle)
 
  283        else if (type == Dune::GeometryTypes::quadrilateral)
 
  288        else if (type == Dune::GeometryTypes::tetrahedron)
 
  293        else if (type == Dune::GeometryTypes::hexahedron)
 
  299            DUNE_THROW(Dune::NotImplemented, 
"Scv geometries for type " << type);
 
 
  305        const auto type = geo_.type();
 
  306        if (type == Dune::GeometryTypes::triangle)
 
  311        else if (type == Dune::GeometryTypes::quadrilateral)
 
  316        else if (type == Dune::GeometryTypes::tetrahedron)
 
  321        else if (type == Dune::GeometryTypes::hexahedron)
 
  327            DUNE_THROW(Dune::NotImplemented, 
"Scvf geometries for type " << type);
 
 
  333        const auto type = geo_.type();
 
  334        if (type == Dune::GeometryTypes::triangle || type == Dune::GeometryTypes::quadrilateral)
 
  336        else if (type == Dune::GeometryTypes::tetrahedron)
 
  338        else if (type == Dune::GeometryTypes::hexahedron)
 
  341            DUNE_THROW(Dune::NotImplemented, 
"Boundary scvf geometries for type " << type);
 
 
  346        return geo_.global(referenceElement(geo_).position(localFacetIndex, 1));
 
 
  351        const auto type = geo_.type();
 
  352        if (type == Dune::GeometryTypes::triangle)
 
  354        else if (type == Dune::GeometryTypes::quadrilateral)
 
  356        else if (type == Dune::GeometryTypes::tetrahedron)
 
  358        else if (type == Dune::GeometryTypes::hexahedron)
 
  361            DUNE_THROW(Dune::NotImplemented, 
"Inside outside scv pairs for type " << type);
 
 
  367        return referenceElement(geo_).size(2);
 
 
  373        return referenceElement(geo_).size(1);
 
 
  376    template<
int d = dimWorld, std::enable_if_t<(d==3), 
int> = 0>
 
  377    GlobalPosition 
normal(
const ScvfCornerStorage& p, 
const std::array<LocalIndexType, 2>& scvPair)
 
  382        const auto ref = referenceElement(geo_);
 
  383        const auto v = facetCenter_(scvPair[1], ref) - facetCenter_(scvPair[0], ref);
 
 
  392    template<
int d = dimWorld, std::enable_if_t<(d==2), 
int> = 0>
 
  393    GlobalPosition 
normal(
const ScvfCornerStorage& p, 
const std::array<LocalIndexType, 2>& scvPair)
 
  396        const auto t = p[1] - p[0];
 
  397        GlobalPosition 
normal({-t[1], t[0]});
 
  400        const auto ref = referenceElement(geo_);
 
  401        const auto v = facetCenter_(scvPair[1], ref) - facetCenter_(scvPair[0], ref);
 
 
  414    template<
class RefElement>
 
  415    GlobalPosition facetCenter_(
unsigned int localFacetIndex, 
const RefElement& ref)
 const 
  417        return geo_.global(ref.position(localFacetIndex, 1));
 
  420    const typename Element::Geometry& geo_; 
 
 
GlobalPosition facetCenter(unsigned int localFacetIndex) const
Definition discretization/facecentered/diamond/geometryhelper.hh:344
ScvfCornerStorage getScvfCorners(unsigned int localEdgeIndex) const
Create a corner storage with the scvf corners for a given edge (codim-2) index.
Definition discretization/facecentered/diamond/geometryhelper.hh:303
std::size_t numInteriorScvf()
number of interior sub control volume faces (number of codim-2 entities)
Definition discretization/facecentered/diamond/geometryhelper.hh:365
std::array< LocalIndexType, 2 > getInsideOutsideScvForScvf(unsigned int localEdgeIndex)
Definition discretization/facecentered/diamond/geometryhelper.hh:349
ScvCornerStorage getScvCorners(unsigned int localFacetIndex) const
Create a corner storage with the scv corners for a given face (codim-1) index.
Definition discretization/facecentered/diamond/geometryhelper.hh:275
DiamondGeometryHelper(const typename Element::Geometry &geo)
Definition discretization/facecentered/diamond/geometryhelper.hh:270
GlobalPosition normal(const ScvfCornerStorage &p, const std::array< LocalIndexType, 2 > &scvPair)
Definition discretization/facecentered/diamond/geometryhelper.hh:377
const Element::Geometry & elementGeometry() const
Definition discretization/facecentered/diamond/geometryhelper.hh:410
ScvfCornerStorage getBoundaryScvfCorners(unsigned int localFacetIndex) const
Create the sub control volume face geometries on the boundary for a given face index.
Definition discretization/facecentered/diamond/geometryhelper.hh:331
std::size_t numScv()
number of sub control volumes (number of codim-1 entities)
Definition discretization/facecentered/diamond/geometryhelper.hh:371
Dune::FieldVector< Scalar, 3 > crossProduct(const Dune::FieldVector< Scalar, 3 > &vec1, const Dune::FieldVector< Scalar, 3 > &vec2)
Cross product of two vectors in three-dimensional Euclidean space.
Definition math.hh:671
Define some often used mathematical functions.
Definition discretization/facecentered/diamond/geometryhelper.hh:39
S keyToCornerStorage(const Geo &geo, const std::array< T, N > &key)
Definition discretization/facecentered/diamond/geometryhelper.hh:189
S keyToCornerStorageImpl(const Geo &geo, const KeyArray &key, std::index_sequence< I... >)
Definition discretization/facecentered/diamond/geometryhelper.hh:179
S boundaryCornerStorageImpl(const Geo &geo, unsigned int i, std::index_sequence< ii... >)
Definition discretization/facecentered/diamond/geometryhelper.hh:196
S boundaryCornerStorage(const Geo &geo, unsigned int i)
Definition discretization/facecentered/diamond/geometryhelper.hh:206
Definition common/pdesolver.hh:24
static constexpr std::array< std::array< IndexType, 2 >, 4 > pairs
Definition discretization/facecentered/diamond/geometryhelper.hh:225
static constexpr std::array< std::array< IndexType, 2 >, 6 > pairs
Definition discretization/facecentered/diamond/geometryhelper.hh:233
static constexpr std::array< std::array< IndexType, 2 >, 12 > pairs
Definition discretization/facecentered/diamond/geometryhelper.hh:241
static constexpr std::array< std::array< IndexType, 2 >, 3 > pairs
Definition discretization/facecentered/diamond/geometryhelper.hh:217
Definition discretization/facecentered/diamond/geometryhelper.hh:212
static constexpr std::array< std::array< Key, 5 >, 6 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:95
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:91
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:94
static constexpr std::array< std::array< Key, 3 >, 4 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:65
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:61
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:64
static constexpr std::array< std::array< Key, 4 >, 4 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:80
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:79
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:76
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:47
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:50
static constexpr std::array< std::array< Key, 3 >, 3 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:51
Definition discretization/facecentered/diamond/geometryhelper.hh:42
static constexpr std::array< std::array< Key, 3 >, 12 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:161
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:157
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:160
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:125
static constexpr std::array< std::array< Key, 2 >, 4 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:129
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:128
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:143
static constexpr std::array< std::array< Key, 3 >, 6 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:144
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:140
std::pair< std::uint8_t, std::uint8_t > Key
Definition discretization/facecentered/diamond/geometryhelper.hh:114
static constexpr std::array< std::array< Key, 2 >, 3 > keys
Definition discretization/facecentered/diamond/geometryhelper.hh:115
static constexpr Dune::GeometryType type()
Definition discretization/facecentered/diamond/geometryhelper.hh:111
Definition discretization/facecentered/diamond/geometryhelper.hh:106
Definition discretization/facecentered/diamond/geometryhelper.hh:34
Dune::ReservedVector< Dune::FieldVector< ct, cdim >,(1<< mydim)> Type
Definition discretization/facecentered/diamond/geometryhelper.hh:35
Traits for an efficient corner storage for fc diamond method.
Definition discretization/facecentered/diamond/geometryhelper.hh:29