vtkMatrixUtilities.h

Go to the documentation of this file.

/*=========================================================================
 
  Program:   Visualization Toolkit
  Module:    vtkMathPrivate.hxx
 
  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
  All rights reserved.
  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
 
     This software is distributed WITHOUT ANY WARRANTY; without even
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
     PURPOSE.  See the above copyright notice for more information.
 
=========================================================================*/
#ifndef vtkMatrixUtilities_h
#define vtkMatrixUtilities_h
 
#include <type_traits> // for type traits
 
namespace vtkMatrixUtilities
{
//=============================================================================
struct Layout
{
  struct Identity;
  /*
   * Input matrix is transposed, i.e. sorted in column-wise ordered.
   */
  struct Transpose;
 
  struct Diag;
};
 
namespace detail
{
// Extracting for STL-like containers
template <int ContainerTypeT, class ContainerT>
struct ScalarTypeExtractor
{
  typedef typename ContainerT::value_type value_type;
  static_assert(std::is_integral<value_type>::value || std::is_floating_point<value_type>::value,
    "value_type is not a numeric type");
};
 
// Extracting for C++ arrays
template <class ContainerT>
struct ScalarTypeExtractor<1, ContainerT>
{
  typedef typename std::remove_pointer<
    typename std::remove_all_extents<typename std::remove_pointer<ContainerT>::type>::type>::type
    value_type;
  static_assert(std::is_integral<value_type>::value || std::is_floating_point<value_type>::value,
    "value_type is not a numeric type");
};
} // namespace detail
 
//=============================================================================
template <class ContainerT>
struct ScalarTypeExtractor
{
private:
  typedef typename std::remove_reference<ContainerT>::type DerefContainer;
 
public:
  typedef typename detail::ScalarTypeExtractor<
    // This parameter equals 0 or 1
    std::is_array<DerefContainer>::value || std::is_pointer<DerefContainer>::value,
    ContainerT>::value_type value_type;
  static_assert(std::is_integral<value_type>::value || std::is_floating_point<value_type>::value,
    "value_type is not a numeric type");
};
 
//-----------------------------------------------------------------------------
template <class MatrixT>
static constexpr bool MatrixIs2DArray()
{
  typedef typename std::remove_extent<MatrixT>::type Row;
  typedef typename std::remove_extent<Row>::type Value;
  return std::is_array<MatrixT>::value && std::is_array<Row>::value && !std::is_array<Value>::value;
}
 
//-----------------------------------------------------------------------------
template <class MatrixT>
static constexpr bool MatrixIsPointerToPointer()
{
  typedef typename std::remove_pointer<MatrixT>::type Row;
  typedef typename std::remove_pointer<Row>::type Value;
  return std::is_pointer<MatrixT>::value && std::is_pointer<Row>::value &&
    !std::is_pointer<Value>::value;
}
 
//-----------------------------------------------------------------------------
template <class MatrixT>
static constexpr bool MatrixLayoutIs2D()
{
  typedef typename std::remove_pointer<MatrixT>::type RowPointer;
  typedef typename std::remove_extent<MatrixT>::type RowArray;
  typedef typename std::remove_pointer<MatrixT>::type ValuePointerPointer;
  typedef typename std::remove_extent<MatrixT>::type ValuePointerArray;
  typedef typename std::remove_pointer<MatrixT>::type ValueArrayPointer;
  typedef typename std::remove_extent<MatrixT>::type ValueArrayArray;
  return ((std::is_array<RowPointer>::value && !std::is_same<RowPointer, MatrixT>::value) ||
           std::is_pointer<RowPointer>::value || std::is_array<RowArray>::value ||
           (std::is_pointer<RowArray>::value && !std::is_same<RowArray, MatrixT>::value)) &&
    (!std::is_array<ValuePointerPointer>::value || !std::is_pointer<ValuePointerPointer>::value) &&
    (!std::is_array<ValueArrayPointer>::value || !std::is_pointer<ValueArrayPointer>::value) &&
    (!std::is_array<ValuePointerArray>::value || !std::is_pointer<ValuePointerArray>::value) &&
    (!std::is_array<ValueArrayArray>::value || !std::is_pointer<ValueArrayArray>::value);
}
 
namespace detail
{
// Class actually implementing matrix mapping.
template <int RowsT, int ColsT, class LayoutT>
struct Mapper;
 
// Specialization of the matrix mapper for when the layout is the identity
template <int RowsT, int ColsT>
struct Mapper<RowsT, ColsT, Layout::Identity>
{
  template <int RowT, int ColT>
  static constexpr int GetIndex()
  {
    static_assert(RowT >= 0 && RowT < RowsT, "RowT out of bounds");
    static_assert(ColT >= 0 && ColT < ColsT, "ColT out of bounds");
    return ColsT * RowT + ColT;
  }
};
 
template <int RowsT, int ColsT>
struct Mapper<RowsT, ColsT, Layout::Transpose>
{
  template <int RowT, int ColT>
  static constexpr int GetIndex()
  {
    static_assert(RowT >= 0 && RowT < RowsT, "RowT out of bounds");
    static_assert(ColT >= 0 && ColT < ColsT, "ColT out of bounds");
    return RowsT * ColT + RowT;
  }
};
} // namespace detail
 
//=============================================================================
template <int RowsT, int ColsT, class LayoutT = Layout::Identity>
struct Mapper
{
  template <int RowT, int ColT>
  static constexpr int GetIndex()
  {
    return detail::Mapper<RowsT, ColsT, LayoutT>::template GetIndex<RowT, ColT>();
  }
};
 
namespace detail
{
// Class implementing matrix wrapping.
template <int RowsT, int ColsT, class MatrixT, class LayoutT, bool MatrixLayoutIs2DT>
class Wrapper;
 
// Specializaion of matrix wrapping for matrices stored as 1D arrays
// in row-wise order
template <int RowsT, int ColsT, class MatrixT, class LayoutT>
class Wrapper<RowsT, ColsT, MatrixT, LayoutT, false>
{
private:
  using Scalar = typename vtkMatrixUtilities::ScalarTypeExtractor<MatrixT>::value_type;
 
public:
  template <int RowT, int ColT>
  static const Scalar& Get(const MatrixT& M)
  {
    return M[Mapper<RowsT, ColsT, LayoutT>::template GetIndex<RowT, ColT>()];
  }
 
  template <int RowT, int ColT>
  static Scalar& Get(MatrixT& M)
  {
    return M[Mapper<RowsT, ColsT, LayoutT>::template GetIndex<RowT, ColT>()];
  }
};
 
// Specialization for matrices stored as 2D arrays with an unchanged layout
template <int RowsT, int ColsT, class MatrixT>
class Wrapper<RowsT, ColsT, MatrixT, Layout::Identity, true>
{
private:
  using Scalar = typename vtkMatrixUtilities::ScalarTypeExtractor<MatrixT>::value_type;
 
public:
  template <int RowT, int ColT>
  static const Scalar& Get(const MatrixT& M)
  {
    return M[RowT][ColT];
  }
 
  template <int RowT, int ColT>
  static Scalar& Get(MatrixT& M)
  {
    return M[RowT][ColT];
  }
};
 
// Specialization for matrices stored as 2D arrays read as its transposed self.
template <int RowsT, int ColsT, class MatrixT>
class Wrapper<RowsT, ColsT, MatrixT, Layout::Transpose, true>
{
private:
  using Scalar = typename vtkMatrixUtilities::ScalarTypeExtractor<MatrixT>::value_type;
 
public:
  template <int RowT, int ColT>
  static const Scalar& Get(const MatrixT& M)
  {
    return M[ColT][RowT];
  }
 
  template <int RowT, int ColT>
  static Scalar& Get(MatrixT& M)
  {
    return M[ColT][RowT];
  }
};
 
// Specialization for diagonal matrices.
// Note: a diagonal matrix has to be stored in a 1D array.
template <int RowsT, int ColsT, class MatrixT>
class Wrapper<RowsT, ColsT, MatrixT, Layout::Diag, false>
{
private:
  using Scalar = typename vtkMatrixUtilities::ScalarTypeExtractor<MatrixT>::value_type;
 
  template <int RowT, int ColT>
  struct Helper
  {
    static constexpr Scalar ZERO = Scalar(0);
 
    static Scalar& Get(const MatrixT&) { return ZERO; }
  };
 
  template <int RowT>
  struct Helper<RowT, RowT>
  {
    static Scalar& Get(MatrixT& M) { return M[RowT]; }
 
    static const Scalar& Get(const MatrixT& M) { return M[RowT]; }
  };
 
public:
  template <int RowT, int ColT>
  const Scalar& Get(const MatrixT& M)
  {
    return Helper<RowT, ColT>::Get(M);
  }
 
  template <int RowT, int ColT>
  Scalar& Get(MatrixT& M)
  {
    return Helper<RowT, ColT>::Get(M);
  }
};
} // namespace detail
 
//=============================================================================
template <int RowsT, int ColsT, class MatrixT, class LayoutT = Layout::Identity>
class Wrapper
{
private:
  using Scalar = typename ScalarTypeExtractor<MatrixT>::value_type;
 
  static_assert(!MatrixLayoutIs2D<MatrixT>() || !std::is_same<LayoutT, Layout::Diag>::value,
    "A diagonal  matrix cannot be a 2D array");
 
public:
  template <int RowT, int ColT>
  static const Scalar& Get(const MatrixT& M)
  {
    return detail::Wrapper<RowsT, ColsT, MatrixT, LayoutT,
      MatrixLayoutIs2D<MatrixT>()>::template Get<RowT, ColT>(M);
  }
 
  template <int RowT, int ColT>
  static Scalar& Get(MatrixT& M)
  {
    return detail::Wrapper<RowsT, ColsT, MatrixT, LayoutT,
      MatrixLayoutIs2D<MatrixT>()>::template Get<RowT, ColT>(M);
  }
};
} // namespace vtkMatrixUtilities
#endif
 
// VTK-HeaderTest-Exclude: vtkMatrixUtilities.h