type traits: need to finalize set of type aliases that all types must have

[fnrizzi]

what is the problem ?

We currently lack a well-defined specification/definition for what a "trait" should be for any container-like/algebra type (vectors, MV, matrices, etc). Need to fix this. In the various components of pressio, we currently have this usage of traits:

• type_traits: this is where things are defined so obviously matters
• expressions: only access scalar_type and size_type of traits
• ops: this is where traits are used the most and most completely. this makes sense because ops are specialized for types. However, this approach is not entirely solid because here is a place where we actually know what the concrete types are here, so have common names to query types is not particularly useful. We could in fact extract the info we need from named types directly: for example, if we (partial) specialize an operation for Eigen, we know that we using Eigen, so we could extract from that Eigen type the info needed. One reason to keep the current approach is that if a library changes names for things, we do NOT have to change everything.
• qr: similar usage to expressions
• solvers_linear, nonlinear: none
• ode_steppers: only scalar_type
• ode_advancers: none
• rom: scalar_type and size_type

These are my observations:

• for the pressio components that are supposed to be "really generic" (solvers_nonlinear, ode, rom) we are okeish, I think, because we have written things in a generic way in terms of ops. So I don't think these would need to change much.
• in ops, however, we are missing "fall back" cases for arbitrary types that we don't know anything about (more on this below) but that we could do something about if we had decent traits.
• once we revise traits, the components mostly affected will be expressions and ops

goal/discussion

we need to have a minimal (easily expandable) well-defined set of traits that we use to reason about types.

For example, we could imrpove things by making the distinction between value_type and element_type like C++ span does: https://en.cppreference.com/w/cpp/container/span And we need to make the traits richer to support a wider specialization flexibility.

Also, long term potentially we want pressio to handle operations on arbitrary types as fall-back scenario when for example users don't provide ops specializations for these arbitrary types. At the very least we need both things that describe things like value type as well as memory features (e.g. shared mem obj, distributed). If we do this well enough, I envision a scenario like this:

// pressio code
namespace pressio{ namespace ops{
 // specialize  gemv for Eigen 
 // specialize  gemv for Kokkos 
 // ... all specializations for libraries we support

 // fallback case: y = A^T*x
 // this would be called only if no other specialization is pickedup
 template<class T> 
 std::enable_if_t< 
   pressio::is_shared_mem_something && 
   is_host_accessible<T> && 
   is_rank_2<T> && 
   subscriptable<T> >
 product(pressio::transpose, const T & v, Eigen::Vector x, Eigen::Vector & y)
 {
   // we know it is sharedmem so we can potentially handle this 
   // even if not super efficiently. if user wants effieciency they have to specialize product for their types.
 }
}}

// user code
namespace pressio{
struct Traits<MyCustomMatrixType>{
  using scalar_type = float;
  // this naming is bad but just to give an idea 
  static constexpr MemEnum objMemSomething = pressio::MemEnum::SharedMemory; // or DistributedMemory
}
}
int main(){ 
  MyCustomMatrixType A(...);
  Eigen::Vector romState (...).
  pressio::rom::something_that_uses_gemv( ...., A, romState);
}

Let's say we have that trait. Then, if the user provides specializations of ops to operate on that great, otherwise we could envision specializing ops to handle that where we can do basic operations because a sufficient number of features is given by the user in the trait to do so. The more we allow the user to refine the trait, the better we can handle types in a generic fashion.

trait levels

The more information user provides in the traits, the better we can handle things (specialize ops):

• level0: the most basic set of traits that are absolutely needed (this is most likely things like value_type, etc) - not enough for us to do any operation ourselves on these types, so user will need to provide implementations
• level1: level0 + info about the type syntax (e.g. is_subscriptable, etc)
• level2: level1 + the trait provides info about the memory layout (left, right, something else) - need to figure out a pressio convention to specify this
• level3: level2 + something else even more specialized

solution plan

I like to solve this with these steps:

• propose a draft of what we think a trait should look like: we need both a minimal set and a more refined one
• draft what trait classes would look like for e.g. Eigen matrix, a Tpetra MV, Kokkos view rank2 and something else maybe
• evaluate and refine

Immediate items to work on:

• [ ] figure out what traits are used and what are NOT used at all so we can simpilfy
• [ ] put in the issue what the trait class looks like (after removing the unused traits) for a vector: what i mean here is to literally write down what the complete trait class looks like for a vector.

  struct Traits<some vector>{
  using scalar_type = 
  // what else is here?
  };

  NOTE that the above is different than just saying this: basic traits - defined as static constexpr in ContainersSharedTraits. That is just an impl detail, we need to see what the trait looks like.

• [ ] checkout what traits Trilinos and PETCs libraries are using, and how their naming convention differs from mdspan and std containes in general
• [ ] based on that we decide how to change names

[mzuzek]

Current set of container type traits

Generic traits

• basic traits - defined as static constexpr in ContainersSharedTraits:
  • PackageIdentifier package_identifier
  • bool is_shared_mem
  • bool is_distributed
  • int rank
• ordinal - defined as typedefs in OrdinalTrait:
  • ordinal_type
  • size_type
• scalar traits - defined as typedefs in ScalarTrait:
  • scalar_type
  • reference_type
  • const_reference_type
• static/dynamic rank traits - defined as static constexpr in AllocTrait:
  • bool is_static
  • derived bool is_dynamic = !is_static;
• matrix density trait - defined as static constexpr in MatrixDensityTrait:
  • bool bool is_sparse
  • derived bool is_dense = !is_sparse;
• layout trait - defined as static constexpr in MatrixLayoutTrait:
  • bool is_row_major
  • derived bool is_col_major = !is_row_major;

Package-specific traits

TPL-specific traits - like Kokkos host_mem_space_type or Trilinos communicator_type - are defined in type_traits/traits_tpl.hpp.

Note: Maybe some of these traits can be generalized (defined for all container types) to become useful in this context. I'm not sure what's the advantage over direct overloads (i.e. functions specialized for Kokkos containers) since these traits can't be used without prior type check anyway.

Type indetification

A container type will define one of listed static constexprs to identify itself:

• VectorIdentifier vector_identifier;
• static constexpr MultiVectorIdentifier multi_vector_identifier;
• static constexpr MatrixIdentifier matrix_identifier;
• static constexpr TensorIdentifier tensor_identifier;

Note: Same information could be carried with simpler and generic container_type enumeration (with values for vector, multivector, matrix, tensor, etc.) combined with already existing package_identifier.

[mzuzek]

Container Traits - idea

@fnrizzi Based on https://github.com/Pressio/pressio/pull/359. Essentially, pretty much only rank, scalar_type and ordinal_type (aliased as size_type) are widely used.

With the container type trait set described by following class:

template< int Rank, typename ScalarType, typename OrdinalType >
struct ContainerTraits
{
  using scalar_type  = ScalarType;
  using ordinal_type = OrdinalType;

  static constexpr int rank = Rank;
};

we could define vector and matrix traits like:

//---------- Eigen vector ----------//
template <typename T>
struct Traits<
    T,
    mpl::enable_if_t<is_vector_eigen<T>::value>
> : public ::pressio::impl::ContainerTraits<
      ::pressio::Vector,
      typename T::Scalar,
      typename T::StorageIndex
> {};

//---------- Kokkos dense matrix ----------//
template <typename T>
struct Traits<
  T,
  ::pressio::mpl::enable_if_t<is_dense_matrix_kokkos<T>::value>
> : public ::pressio::impl::ContainerTraits<
      ::pressio::Matrix,
      typename T::traits::value_type,
      typename T::traits::size_type
> {};

//---------- Tpetra multi-vector ----------//
template<typename T>
struct Traits<
  T,
  ::pressio::mpl::enable_if_t<
    is_multi_vector_tpetra<T>::value
    || is_multi_vector_tpetra_block<T>::value>
> : public ::pressio::impl::ContainerTraits<
      ::pressio::MultiVector,
      typename T::impl_scalar_type,
      typename T::local_ordinal_type
> {};

with the following remarks:

• Library specific traits that are rarely used and only after container type is already semi-resolved (e.g. to Kokkos vector) can be replaced with direct type properties (e.g. use T::memory_space instead of Traits<T>::memory_space);
• Identification like package_identifier (e.g. PackageIdentifier::Eigen) can be replaced with type detection helpers like is_eigen_vector<T>;
• Helpers can handle library specific traits that aren't defined for all containers anyway, e.g.: device_type<T>, execution_space<T>, have_matching_execution_space<T1, T2>, have_matching_device_type<T1, T2>.
• Helpers could be also used for reference_type and const_reference_type, as they seem to carry out expression related logic: https://github.com/Pressio/pressio/blob/a9eb57824bec75953581b6e6cc57d18c1aa9286a/include/pressio/expressions/impl/asdiagonalmatrix_traits.hpp#L86-L92

Edit 2022-06-30: plain int for rank, removed size_type and {const_}reference_type

[mzuzek]

Traits and their semantics in different libraries

@fnrizzi

Basic Set

trait	Pressio	Eigen	Kokkos	Tpetra	[PETSc]()
type of element value	scalar_type	Scalar	value_type const_value_type non_const_value_type	scalar_type	fixed ³
type of element index	ordinal_type	StorageIndex	arbitrary/templated ¹	local_ordinal_type global_ordinal_type ²	fixed ³
container rank	rank	-	rank rank_dynamic ⁴	-	-

¹ there is size_type which derives from memory_space::size_type and seems to be pointer-like byte index instead of logical element index ² for multi-node distributed computing ³ uses PetscScalar for values and PetscInt for indexes (configured via compile flags) ⁴ static vs dynamic allocation

Extended Set

Candidate traits:

trait	semantics	purpose
reference_type const_reference_type	type of reference to container value	TODO: explain what gains would it provide over auto-detection (in expressions, where these are used)

[fnrizzi]

Traits

• using pressio::Traits class

• simplify to only have scalar_type, rank, get rid of everything else pretty muhc

predicates based on trait class

• we decided on:

  • ::pressio::isvector{eigen, kokkos}: this identifies exactly a Kokkos View

  • ::pressio::is_expression_acting_on_kokkos<>

OPS

• we need to chagne all constraints to conform to the above

Stages of work

1. make list of pressio components and list all traits used inside each of them (most likely this is scalar_type and nothing else excpet for ops)

2. do changes as said above

[mzuzek]

@fnrizzi

Container identifying predicates

Here are some charts that may help navigating container predicates and their relationships. They can also help in completeness verification if (e.g. test coverage).

[fnrizzi]

thanks for making this!

[mzuzek]

@fnrizzi Would you see anything blocking v0.13.0 in this task at this point ?

[fnrizzi]

leaving this open for now