type.hh

Go to the documentation of this file.

// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
// SPDX-FileCopyrightInfo: Copyright © DUNE Project contributors, see file LICENSE.md in module root
// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
#ifndef DUNE_GEOMETRY_TYPE_HH
#define DUNE_GEOMETRY_TYPE_HH
 
#include <cassert>
#include <cstdint>
 
#include <string>
#include <type_traits>
 
#include <dune/common/exceptions.hh>
#include <dune/common/typetraits.hh>
 
namespace Dune
{
 
  namespace Impl
  {
 
    enum TopologyConstruction { pyramidConstruction = 0, prismConstruction = 1 };
 
    // Dynamic Topology Properties
    // ---------------------------
 
    inline static unsigned int numTopologies ( int dim ) noexcept
    {
      return (1u << dim);
    }
 
    inline bool static isPyramid ( unsigned int topologyId, int dim, int codim = 0 ) noexcept
    {
      assert( (dim > 0) && (topologyId < numTopologies( dim )) );
      assert( (0 <= codim) && (codim < dim) );
      return (((topologyId & ~1) & (1u << (dim-codim-1))) == 0);
    }
 
    inline static bool isPrism ( unsigned int topologyId, int dim, int codim = 0 ) noexcept
    {
      assert( (dim > 0) && (topologyId < numTopologies( dim )) );
      assert( (0 <= codim) && (codim < dim) );
      return (( (topologyId | 1) & (1u << (dim-codim-1))) != 0);
    }
 
    inline static unsigned int baseTopologyId ( unsigned int topologyId, int dim, int codim = 1 ) noexcept
    {
      assert( (dim >= 0) && (topologyId < numTopologies( dim )) );
      assert( (0 <= codim) && (codim <= dim) );
      return topologyId & ((1u << (dim-codim)) - 1);
    }
 
  } // namespace Impl
 
  // GeometryType
  // -------------
 
  class GeometryType
  {
  public:
 
    enum
    BasicType {
      simplex,       
      cube,          
      pyramid,       
      prism,         
      extended,      
      none           
    };
 
  private:
 
    unsigned char dim_;
 
    bool none_;
 
    unsigned int topologyId_;
 
    // Internal type used for the Id. The exact nature of this type is kept
    // as an implementation detail on purpose. We use a scoped enum here because scoped enums
    // can be used as template parameters, but are not implicitly converted to other integral
    // types by the compiler. That way, we avoid unfortunate implicit conversion chains, e.g.
    // people trying to work with GlobalGeometryTypeIndex, but forgetting to actually call
    // GlobalGeometryTypeIndex::index(gt) and just using gt directly.
    enum class IdType : std::uint64_t
    {};
 
  public:
 
    using Id = IdType;
 
    constexpr operator Id() const
    {
      // recreate the exact storage layout that this class is using, making conversion
      // extremely cheap
      std::uint64_t id = dim_ | (std::uint64_t(none_) << 8) | (std::uint64_t(topologyId_) << 32);
      return static_cast<Id>(id);
    }
 
    constexpr Id toId() const
    {
      return static_cast<Id>(*this);
    }
 
    constexpr GeometryType(Id id)
      : dim_(static_cast<std::uint64_t>(id) & 0xFF)
      , none_(static_cast<std::uint64_t>(id) & 0x100)
      , topologyId_(static_cast<std::uint64_t>(id) >> 32)
    {}
 
 
    constexpr GeometryType ()
      : dim_(0), none_(true), topologyId_(0)
    {}
 
    constexpr GeometryType(unsigned int topologyId, unsigned int dim, bool isNone)
      : dim_(dim), none_(isNone), topologyId_(topologyId)
    {}
 
    constexpr GeometryType(unsigned int topologyId, unsigned int dim)
      : dim_(dim), none_(false), topologyId_(topologyId)
    {}
 
    template<class TopologyType,
      class = std::void_t<decltype(TopologyType::dimension), decltype(TopologyType::id)>>
    explicit GeometryType([[maybe_unused]] TopologyType t)
      : dim_(TopologyType::dimension), none_(false), topologyId_(TopologyType::id)
    {}
 
    constexpr bool isVertex() const {
      return dim_==0;
    }
 
    constexpr bool isLine() const {
      return dim_==1;
    }
 
    constexpr bool isTriangle() const {
      return ! none_ && dim_==2 && (topologyId_ | 1) == 0b0001;
    }
 
    constexpr bool isQuadrilateral() const {
      return ! none_ && dim_==2 && (topologyId_ | 1) == 0b0011;
    }
 
    constexpr bool isTetrahedron() const {
      return ! none_ && dim_==3 && (topologyId_ | 1) == 0b0001;
    }
 
    constexpr bool isPyramid() const {
      return ! none_ && dim_==3 && (topologyId_ | 1) == 0b0011;
    }
 
    constexpr bool isPrism() const {
      return ! none_ && dim_==3 && (topologyId_ | 1) == 0b0101;
    }
 
    constexpr bool isHexahedron() const {
      return ! none_ && dim_==3 && (topologyId_ | 1) == 0b0111;
    }
 
    constexpr bool isSimplex() const {
      return ! none_ && (topologyId_ | 1) == 1;
    }
 
    constexpr bool isCube() const {
      return ! none_ && ((topologyId_ ^ ((1 << dim_)-1)) >> 1 == 0);
    }
 
    constexpr bool isConical() const {
      return ! none_ && (((topologyId_ & ~1) & (1u << (dim_-1))) == 0);
    }
 
    constexpr bool isConical(const int& step) const {
      return ! none_ && (((topologyId_ & ~1) & (1u << step)) == 0);
    }
 
    constexpr bool isPrismatic() const {
      return ! none_ && (( (topologyId_ | 1) & (1u << (dim_-1))) != 0);
    }
 
    constexpr bool isPrismatic(const int& step) const {
      return ! none_ && (( (topologyId_ | 1) & (1u << step)) != 0);
    }
 
    constexpr bool isNone() const {
      return none_;
    }
 
    constexpr unsigned int dim() const {
      return dim_;
    }
 
    constexpr unsigned int id() const {
      return topologyId_;
    }
 
 
    constexpr bool operator==(const GeometryType& other) const {
      return ( ( none_ == other.none_ )
               && ( ( none_ == true )
                    || ( ( dim_ == other.dim_ )
                         && ( (topologyId_ >> 1) == (other.topologyId_ >> 1) )
                         )
                    )
               );
    }
 
    constexpr bool operator!=(const GeometryType& other) const {
      return ! ((*this)==other);
    }
 
    constexpr bool operator < (const GeometryType& other) const {
      return ( ( none_ < other.none_ )
               || ( !( other.none_ < none_ )
                    && ( ( dim_ < other.dim_ )
                         || ( (other.dim_ == dim_)
                              && ((topologyId_ >> 1) < (other.topologyId_ >> 1) )
                              )
                         )
                    )
               );
    }
 
  };
 
  inline std::ostream& operator<< (std::ostream& s, const GeometryType& a)
  {
    if (a.isSimplex())
    {
      s << "(simplex, " << a.dim() << ")";
      return s;
    }
    if (a.isCube())
    {
      s << "(cube, " << a.dim() << ")";
      return s;
    }
    if (a.isPyramid())
    {
      s << "(pyramid, 3)";
      return s;
    }
    if (a.isPrism())
    {
      s << "(prism, 3)";
      return s;
    }
    if (a.isNone())
    {
      s << "(none, " << a.dim() << ")";
      return s;
    }
    s << "(other [" << a.id() << "], " << a.dim() << ")";
    return s;
  }
 
 
 
  namespace GeometryTypes {
 
 
    inline constexpr GeometryType simplex(unsigned int dim)
    {
      return GeometryType(0,dim,false);
    }
 
 
    inline constexpr GeometryType cube(unsigned int dim)
    {
      return GeometryType(((dim>1) ? ((1 << dim) - 1) : 0),dim,false);
    }
 
 
    inline constexpr GeometryType none(unsigned int dim)
    {
      return GeometryType(0,dim,true);
    }
 
    inline constexpr GeometryType conicalExtension(const GeometryType& gt)
    {
      return GeometryType(gt.id(), gt.dim()+1, gt.isNone());
    }
 
    inline constexpr GeometryType prismaticExtension(const GeometryType& gt)
    {
      return GeometryType(gt.id() | ((1 << gt.dim())), gt.dim()+1, gt.isNone());
    }
 
 
    inline constexpr GeometryType vertex = GeometryType(0,0,false);
 
 
    inline constexpr GeometryType line = GeometryType(0,1,false);
 
 
    inline constexpr GeometryType triangle = simplex(2);
 
 
    inline constexpr GeometryType quadrilateral = cube(2);
 
 
    inline constexpr GeometryType tetrahedron = simplex(3);
 
 
    inline constexpr GeometryType pyramid = GeometryType(0b0011,3,false);
 
 
    inline constexpr GeometryType prism = GeometryType(0b0101,3,false);
 
 
    inline constexpr GeometryType hexahedron = cube(3);
 
  }
 
  namespace Impl
  {
 
    inline constexpr GeometryType getBase(const GeometryType& gt) {
      return GeometryType(gt.id() & ((1 << (gt.dim()-1))-1), gt.dim()-1, gt.isNone());
    }
 
 
    // IfGeometryType
    // ----------
 
    template< template< GeometryType::Id > class Operation, int dim, GeometryType::Id geometryId = GeometryTypes::vertex >
    struct IfGeometryType
    {
      static constexpr GeometryType geometry = geometryId;
      template< class... Args >
      static auto apply ( GeometryType gt, Args &&... args )
      {
        GeometryType lowerGeometry(gt.id() >>1 , gt.dim()-1, gt.isNone());
 
        if( gt.id() & 1 )
          return IfGeometryType< Operation, dim-1, GeometryTypes::prismaticExtension(geometry).toId() >::apply( lowerGeometry, std::forward< Args >( args )... );
        else
          return IfGeometryType< Operation, dim-1, GeometryTypes::conicalExtension(geometry).toId() >::apply( lowerGeometry, std::forward< Args >( args )... );
      }
    };
 
    template< template< GeometryType::Id > class Operation, GeometryType::Id geometryId >
    struct IfGeometryType< Operation, 0, geometryId>
    {
      template< class... Args >
      static auto apply ([[maybe_unused]] GeometryType gt, Args &&... args )
      {
        return Operation< geometryId >::apply( std::forward< Args >( args )... );
      }
    };
  } // namespace Impl
} // namespace Dune
 
#endif // DUNE_GEOMETRY_TYPE_HH