Dia/Kokkos/Morpheus_MatrixOperations_Impl.hpp

 
#ifndef MORPHEUS_DIA_KOKKOS_MATRIXOPERATIONS_IMPL_HPP
#define MORPHEUS_DIA_KOKKOS_MATRIXOPERATIONS_IMPL_HPP
 
#include <Morpheus_Exceptions.hpp>
#include <Morpheus_SpaceTraits.hpp>
#include <Morpheus_FormatTraits.hpp>
#include <Morpheus_FormatTags.hpp>
#include <Morpheus_Spaces.hpp>
 
namespace Morpheus {
namespace Impl {
 
template <typename ExecSpace, typename Matrix, typename Vector>
void update_diagonal(
    Matrix& A, const Vector& diagonal,
    typename std::enable_if_t<
        Morpheus::is_dia_matrix_format_container_v<Matrix> &&
        Morpheus::is_dense_vector_format_container_v<Vector> &&
        Morpheus::has_generic_backend_v<ExecSpace> &&
        Morpheus::has_access_v<ExecSpace, Matrix, Vector>>* = nullptr) {
  using execution_space  = typename ExecSpace::execution_space;
  using value_array_type = typename Matrix::value_array_type::value_array_type;
  using index_array_type = typename Matrix::index_array_type::value_array_type;
  using array_type       = typename Vector::value_array_type;
  using size_type        = typename Matrix::size_type;
  using index_type       = typename Matrix::index_type;
  using value_type       = typename Matrix::value_type;
  using range_policy =
      Kokkos::RangePolicy<Kokkos::IndexType<size_type>, execution_space>;
 
  value_array_type values           = A.values().view();
  index_array_type diagonal_offsets = A.diagonal_offsets().view();
  const array_type diagonal_view    = diagonal.const_view();
 
  size_type ndiag = A.values().ncols(), ncols = A.ncols();
 
  range_policy policy(0, A.nrows());
 
  Kokkos::parallel_for(
      policy, KOKKOS_LAMBDA(const size_type row) {
        for (size_type n = 0; n < ndiag; n++) {
          const index_type col = row + diagonal_offsets[n];
 
          if ((col >= 0 && col < (index_type)ncols) &&
              (col == (index_type)row)) {
            values(row, n) =
                (values(row, n) == value_type(0)) ? 0 : diagonal_view[col];
          }
        }
      });
}
 
template <typename ExecSpace, typename Matrix, typename Vector>
void get_diagonal(
    Matrix&, const Vector&,
    typename std::enable_if_t<
        Morpheus::is_dia_matrix_format_container_v<Matrix> &&
        Morpheus::is_dense_vector_format_container_v<Vector> &&
        Morpheus::has_generic_backend_v<ExecSpace> &&
        Morpheus::has_access_v<ExecSpace, Matrix, Vector>>* = nullptr) {
  throw Morpheus::NotImplementedException("get_diagonal not implemented yet");
}
 
template <typename ExecSpace, typename Matrix, typename SizeType,
          typename ValueType>
void set_value(Matrix&, SizeType, SizeType, ValueType,
               typename std::enable_if_t<
                   Morpheus::is_dia_matrix_format_container_v<Matrix> &&
                   Morpheus::has_generic_backend_v<ExecSpace> &&
                   Morpheus::has_access_v<ExecSpace, Matrix>>* = nullptr) {
  throw Morpheus::NotImplementedException("set_value not implemented yet");
}
 
template <typename ExecSpace, typename Matrix, typename IndexVector,
          typename ValueVector>
void set_values(
    Matrix&, typename IndexVector::value_type, const IndexVector,
    typename IndexVector::value_type, const IndexVector, const ValueVector,
    typename std::enable_if_t<
        Morpheus::is_dia_matrix_format_container_v<Matrix> &&
        Morpheus::is_dense_vector_format_container_v<IndexVector> &&
        Morpheus::is_dense_vector_format_container_v<ValueVector> &&
        Morpheus::has_generic_backend_v<ExecSpace> &&
        Morpheus::has_access_v<ExecSpace, Matrix, IndexVector, ValueVector>>* =
        nullptr) {
  throw Morpheus::NotImplementedException("set_values not implemented yet");
}
 
template <typename ExecSpace, typename Matrix, typename TransposeMatrix>
void transpose(
    const Matrix&, TransposeMatrix&,
    typename std::enable_if_t<
        Morpheus::is_dia_matrix_format_container_v<Matrix> &&
        Morpheus::is_dia_matrix_format_container_v<TransposeMatrix> &&
        Morpheus::has_generic_backend_v<ExecSpace> &&
        Morpheus::has_access_v<ExecSpace, Matrix, TransposeMatrix>>* =
        nullptr) {
  throw Morpheus::NotImplementedException("transpose not implemented yet");
}
 
}  // namespace Impl
}  // namespace Morpheus
 
#endif  // MORPHEUS_DIA_KOKKOS_MATRIXOPERATIONS_IMPL_HPP