TensorIndices.h

Go to the documentation of this file.

#ifndef LARDATA_UTILITIES_TENSORINDICES_H
#define LARDATA_UTILITIES_TENSORINDICES_H


// C/C++ standard libraries
#include <cstddef> // std::size_t
#include <stdexcept> // std::out_of_range
#include <string> // std::to_string()
#include <type_traits> // std::enable_if


namespace util {
  
  
  struct TensorIndicesBasicTypes {
    
    using Index_t = std::ptrdiff_t;
    
    using DimSize_t = std::size_t;
    
    using LinIndex_t = std::size_t;
    
  }; // TensorIndicesBasicTypes
  
  
  // This is the class declaration; it has a specialisation for rank 1 and
  // a generic implementation (valid for rank larger than 1).
  template <unsigned int RANK>
  class TensorIndices;
  
  
  namespace details {
    template <unsigned int RANK, unsigned int DIM>
    struct ExtractTensorDimension;
  } // namespace details

  
  template <>
  class TensorIndices<1U> {
      public:
    
    using Index_t    = TensorIndicesBasicTypes::Index_t   ;
    
    using DimSize_t  = TensorIndicesBasicTypes::DimSize_t ;
    
    using LinIndex_t = TensorIndicesBasicTypes::LinIndex_t;
    
    static constexpr unsigned int rank() { return 1U; }
    
    
    TensorIndices(DimSize_t dim): dimSize(dim) {}
    
    template <
      typename ITER,
      typename = std::enable_if_t
        <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, void>
      >
    TensorIndices(ITER dimIter): TensorIndices(*dimIter) {}
    
    
    LinIndex_t operator() (Index_t index) const { return index; }
    
    template <typename ITER>
    std::enable_if_t
      <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, LinIndex_t>
    operator() (ITER indexIter) const { return *indexIter; }
    
    LinIndex_t at(Index_t index) const { return checkOuterIndex(index); }
    
    template <typename ITER>
    std::enable_if_t
      <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, LinIndex_t>
    at(ITER indexIter) const
      { return checkOuterIndex(*indexIter); }
    
    bool has(Index_t index) const
      { return (index >= 0) && ((DimSize_t) index < size()); }
    
    template <typename ITER>
    std::enable_if_t
      <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, bool>
    has(ITER indexIter) const
      { return (*indexIter >= 0) && ((DimSize_t) *indexIter < size()); }
    
    
    template <unsigned int DIM = 0>
    DimSize_t size() const { return (DIM == 0)? totalSize(): dim<DIM>(); }
    
    template <unsigned int DIM>
    DimSize_t dim() const
      {
        static_assert(DIM == 0, "Invalid dimension requested");
        return dimSize;
      } // dim()
    
    template <unsigned int DIM>
    bool hasIndex(Index_t index) const
      { 
        static_assert(DIM == 0, "Invalid dimension requested");
        return has(index);
      }
    
    bool hasLinIndex(LinIndex_t linIndex) const
      { return has((Index_t) linIndex); }
    
    bool operator== (TensorIndices<1U> const& t) const
      { return t.size() == size(); }
    bool operator!= (TensorIndices<1U> const& t) const
      { return t.size() != size(); }
    
      protected:
    Index_t checkOuterIndex(Index_t index) const
      {
        if (has(index)) return index; // good
        throw std::out_of_range("Requested index " + std::to_string(index)
          + " for a dimension of size " + std::to_string(size()));
      }
    
      protected:
    DimSize_t dim0() const { return dimSize; }
    
      private:
    template <unsigned int R, unsigned int D>
    friend struct details::ExtractTensorDimension;
    
    DimSize_t dimSize; 
    
    DimSize_t totalSize() const { return dim0(); }
    
  }; // class TensorIndices<1>
  
  
  template <unsigned int RANK>
  class TensorIndices: private TensorIndices<1U> {
    static_assert(RANK > 1, "TensorIndices must have rank 1 or higher");
    
    using Base_t = TensorIndices<1U>; 
    
      public:
    
    using Index_t    = Base_t::Index_t   ;
    
    using DimSize_t  = Base_t::DimSize_t ;
    
    using LinIndex_t = Base_t::LinIndex_t;
    
    
    static constexpr unsigned int rank() { return RANK; }
    
    using MinorTensor_t = TensorIndices<rank() - 1>;
    
    
    template <typename... OTHERDIMS>
    TensorIndices(DimSize_t first, OTHERDIMS... others)
      : Base_t(first), m(others...), totSize(dim0() * m.size()) {}
    
    
    // only if ITER can be dereferenced into something convertible to DimSize_t
    template <
      typename ITER,
      typename = std::enable_if_t
        <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, void>
      >
    TensorIndices(ITER dimIter)
      : Base_t(dimIter), m(++dimIter), totSize(dim0() * m.size()) {}
    
    
    template <typename... OTHERINDICES>
    LinIndex_t operator() (Index_t first, OTHERINDICES... others) const
      { return first * minor().size() + minor()(others...); }
    
    
    template <typename ITER>
    std::enable_if_t
      <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, LinIndex_t>
    operator() (ITER indexIter) const
      {
        auto const baseSize = (*indexIter) * minor().size();
        return baseSize + minor()(++indexIter);
      }
    
    
    template <typename... OTHERINDICES>
    LinIndex_t at(Index_t first, OTHERINDICES... others) const
      {
        return Base_t::checkOuterIndex(first) * minor().size()
          + minor().at(others...);
      }
    
    template <typename ITER>
    std::enable_if_t
      <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, LinIndex_t>
    at(ITER indexIter) const
      {
        auto const baseSize
          = Base_t::checkOuterIndex(*indexIter) * minor().size();
        return baseSize + minor()(++indexIter);
      }
    
    
    template <typename... OTHERINDICES>
    bool has(Index_t first, OTHERINDICES... others) const
      { return Base_t::has(first) && minor().has(others...); }
    
    template <typename ITER>
    std::enable_if_t
      <std::is_convertible<decltype(*(ITER())), DimSize_t>::value, bool>
    has(ITER indexIter) const
      { return Base_t::has(*indexIter)? minor().has(++indexIter): false; }
    
    
    template <unsigned int DIM>
    DimSize_t dim() const
      { return details::ExtractTensorDimension<rank(), DIM>::dim(*this); }
    
    
    template <unsigned int DIM>
    bool hasIndex(Index_t index) const
      { return (index >= 0U) && ((DimSize_t) index < dim<DIM>()); }
    
    
    template <unsigned int DIM = 0>
    DimSize_t size() const
      { 
        return (DIM == 0)? totalSize(): 
          details::ExtractTensorDimension<rank(), DIM>::size(*this);
      }
    
    bool hasLinIndex(LinIndex_t linIndex) const
      { return (linIndex >= 0U) && ((DimSize_t) linIndex < size()); }
    
    
    MinorTensor_t const& minor() const { return m; }
    
    bool operator== (TensorIndices<RANK> const& t) const
      { return Base_t::operator==(t) && (minor() == t.minor()); }
    
    bool operator!= (TensorIndices<RANK> const& t) const
      { return Base_t::operator!=(t) || (minor() != t.minor()); }
    
    
    
      protected:
    // need to be friend of the dimension extractor
    template <unsigned int R, unsigned int D>
    friend struct details::ExtractTensorDimension;
    
    MinorTensor_t m;   
    DimSize_t totSize; 
    
    DimSize_t totalSize() const { return totSize; }
    
  }; // class TensorIndices<>
  
  
  
  template <
    unsigned int RANK1, unsigned int RANK2,
    typename = std::enable_if_t<(RANK1 != RANK2), bool>
    >
  bool operator== (TensorIndices<RANK1> const& a, TensorIndices<RANK2> const& b)
    { return false; }
  
  template <
    unsigned int RANK1, unsigned int RANK2,
    typename = std::enable_if_t<(RANK1 != RANK2), bool>
    >
  bool operator!= (TensorIndices<RANK1> const& a, TensorIndices<RANK2> const& b)
    { return true; }
  
  
  template <typename... DIMS>
  auto makeTensorIndices(DIMS... dims)
    {
      return TensorIndices<sizeof...(DIMS)>
        { TensorIndicesBasicTypes::DimSize_t(dims)... };
    }
  
  
  using MatrixIndices = TensorIndices<2U>;
  
  
} // namespace util


//------------------------------------------------------------------------------
//--- details implementation
//---
namespace util {
  namespace details {
    
    template <unsigned int RANK, unsigned int DIM>
    struct ExtractTensorDimension {
      static_assert(RANK > DIM, "Invalid dimension requested");
      static TensorIndicesBasicTypes::DimSize_t dim
        (TensorIndices<RANK> const& t)
        { return ExtractTensorDimension<RANK-1U, DIM-1U>::dim(t.minor()); }
      
      static TensorIndicesBasicTypes::DimSize_t size
        (TensorIndices<RANK> const& t)
        { return ExtractTensorDimension<RANK-1U, DIM-1U>::size(t.minor()); }
      
    }; // ExtractTensorDimension<>()
    
    template <unsigned int RANK>
    struct ExtractTensorDimension<RANK, 0U> {
      static_assert(RANK > 0, "Invalid rank 0 for TensorIndices");
      static TensorIndicesBasicTypes::DimSize_t dim
        (TensorIndices<RANK> const& t)
        { return t.dim0(); }
      
      static TensorIndicesBasicTypes::DimSize_t size
        (TensorIndices<RANK> const& t)
        { return t.size(); }
      
    }; // ExtractTensorDimension<RANK, 0>()

  } // namespace details
} // namespace util
  


#endif // LARDATA_UTILITIES_TENSORINDICES_H