MatrixCSR.hpp

Go to the documentation of this file.

// This code is based on Jet framework.
// Copyright (c) 2018 Doyub Kim
// CubbyFlow is voxel-based fluid simulation engine for computer games.
// Copyright (c) 2020 CubbyFlow Team
// Core Part: Chris Ohk, Junwoo Hwang, Jihong Sin, Seungwoo Yoo
// AI Part: Dongheon Cho, Minseo Kim
// We are making my contributions/submissions to this project solely in our
// personal capacity and are not conveying any rights to any intellectual
// property of any third parties.

#ifndef CUBBYFLOW_MATRIX_CSR_HPP
#define CUBBYFLOW_MATRIX_CSR_HPP

#include <Core/Matrix/Matrix.hpp>

namespace CubbyFlow
{
template <typename T>
class MatrixCSR;

template <typename T, typename ME>
class MatrixCSRMatrixMul
    : public MatrixExpression<T, MATRIX_SIZE_DYNAMIC, MATRIX_SIZE_DYNAMIC,
                              MatrixCSRMatrixMul<T, ME>>
{
 public:
    MatrixCSRMatrixMul(const MatrixCSR<T>& m1, const ME& m2);

    [[nodiscard]] size_t GetRows() const;

    [[nodiscard]] size_t GetCols() const;

    T operator()(size_t i, size_t j) const;

 private:
    const MatrixCSR<T>& m_m1;
    const ME& m_m2;
    const T* const m_nnz;
    const size_t* const m_rp;
    const size_t* const m_ci;
};

template <typename T>
class MatrixCSR final
    : public MatrixExpression<T, MATRIX_SIZE_DYNAMIC, MATRIX_SIZE_DYNAMIC,
                              MatrixCSR<T>>
{
 public:
    static_assert(
        std::is_floating_point<T>::value,
        "MatrixCSR only can be instantiated with floating point types");

    struct Element
    {
        size_t i;
        size_t j;
        T value;

        Element();

        Element(size_t i, size_t j, const T& value);
    };

    using NonZeroContainerType = std::vector<T>;
    using NonZeroIterator = typename NonZeroContainerType::iterator;
    using ConstNonZeroIterator = typename NonZeroContainerType::const_iterator;

    using IndexContainerType = std::vector<size_t>;
    using IndexIterator = IndexContainerType::iterator;
    using ConstIndexIterator = IndexContainerType::const_iterator;

    MatrixCSR();

    MatrixCSR(const std::initializer_list<std::initializer_list<T>>& lst,
              T epsilon = std::numeric_limits<T>::epsilon());

    template <size_t R, size_t C, typename E>
    MatrixCSR(const MatrixExpression<T, R, C, E>& other,
              T epsilon = std::numeric_limits<T>::epsilon());

    ~MatrixCSR() = default;

    MatrixCSR(const MatrixCSR& other);

    MatrixCSR(MatrixCSR&& other) noexcept;

    MatrixCSR& operator=(const MatrixCSR& other);

    MatrixCSR& operator=(MatrixCSR&& other) noexcept;

    void Clear();

    void Set(const T& s);

    void Set(const MatrixCSR& other);

    void Reserve(size_t rows, size_t cols, size_t numNonZeros);

    void Compress(const std::initializer_list<std::initializer_list<T>>& lst,
                  T epsilon = std::numeric_limits<T>::epsilon());

    template <size_t R, size_t C, typename E>
    void Compress(const MatrixExpression<T, R, C, E>& other,
                  T epsilon = std::numeric_limits<T>::epsilon());

    void AddElement(size_t i, size_t j, const T& value);

    void AddElement(const Element& element);

    void AddRow(const NonZeroContainerType& nonZeros,
                const IndexContainerType& columnIndices);

    void SetElement(size_t i, size_t j, const T& value);

    void SetElement(const Element& element);

    [[nodiscard]] bool IsEqual(const MatrixCSR& other) const;

    [[nodiscard]] bool IsSimilar(
        const MatrixCSR& other,
        double tol = std::numeric_limits<double>::epsilon()) const;

    [[nodiscard]] bool IsSquare() const;

    [[nodiscard]] Vector2UZ Size() const;

    [[nodiscard]] size_t GetRows() const;

    [[nodiscard]] size_t GetCols() const;

    [[nodiscard]] size_t NumberOfNonZeros() const;

    [[nodiscard]] const T& NonZero(size_t i) const;

    [[nodiscard]] T& NonZero(size_t i);

    [[nodiscard]] const size_t& RowPointer(size_t i) const;

    [[nodiscard]] const size_t& ColumnIndex(size_t i) const;

    [[nodiscard]] T* NonZeroData();

    [[nodiscard]] const T* NonZeroData() const;

    [[nodiscard]] const size_t* RowPointersData() const;

    [[nodiscard]] const size_t* ColumnIndicesData() const;

    [[nodiscard]] NonZeroIterator NonZeroBegin();

    [[nodiscard]] ConstNonZeroIterator NonZeroBegin() const;

    [[nodiscard]] NonZeroIterator NonZeroEnd();

    [[nodiscard]] ConstNonZeroIterator NonZeroEnd() const;

    [[nodiscard]] IndexIterator RowPointersBegin();

    [[nodiscard]] ConstIndexIterator RowPointersBegin() const;

    [[nodiscard]] IndexIterator RowPointersEnd();

    [[nodiscard]] ConstIndexIterator RowPointersEnd() const;

    [[nodiscard]] IndexIterator ColumnIndicesBegin();

    [[nodiscard]] ConstIndexIterator ColumnIndicesBegin() const;

    [[nodiscard]] IndexIterator ColumnIndicesEnd();

    [[nodiscard]] ConstIndexIterator ColumnIndicesEnd() const;

    [[nodiscard]] MatrixCSR Add(const T& s) const;

    [[nodiscard]] MatrixCSR Add(const MatrixCSR& m) const;

    [[nodiscard]] MatrixCSR Sub(const T& s) const;

    [[nodiscard]] MatrixCSR Sub(const MatrixCSR& m) const;

    [[nodiscard]] MatrixCSR Mul(const T& s) const;

    template <size_t R, size_t C, typename ME>
    [[nodiscard]] MatrixCSRMatrixMul<T, ME> Mul(
        const MatrixExpression<T, R, C, ME>& m) const;

    [[nodiscard]] MatrixCSR Div(const T& s) const;

    [[nodiscard]] MatrixCSR RAdd(const T& s) const;

    [[nodiscard]] MatrixCSR RAdd(const MatrixCSR& m) const;

    [[nodiscard]] MatrixCSR RSub(const T& s) const;

    [[nodiscard]] MatrixCSR RSub(const MatrixCSR& m) const;

    [[nodiscard]] MatrixCSR RMul(const T& s) const;

    [[nodiscard]] MatrixCSR RDiv(const T& s) const;

    void IAdd(const T& s);

    void IAdd(const MatrixCSR& m);

    void ISub(const T& s);

    void ISub(const MatrixCSR& m);

    void IMul(const T& s);

    template <size_t R, size_t C, typename ME>
    void IMul(const MatrixExpression<T, R, C, ME>& m);

    void IDiv(const T& s);

    [[nodiscard]] T Sum() const;

    [[nodiscard]] T Avg() const;

    [[nodiscard]] T Min() const;

    [[nodiscard]] T Max() const;

    [[nodiscard]] T AbsMin() const;

    [[nodiscard]] T AbsMax() const;

    [[nodiscard]] T Trace() const;

    template <typename U>
    MatrixCSR<U> CastTo() const;

    template <size_t R, size_t C, typename ME>
    MatrixCSR& operator=(const MatrixExpression<T, R, C, ME>& m);

    MatrixCSR& operator+=(const T& s);

    MatrixCSR& operator+=(const MatrixCSR& m);

    MatrixCSR& operator-=(const T& s);

    MatrixCSR& operator-=(const MatrixCSR& m);

    MatrixCSR& operator*=(const T& s);

    template <size_t R, size_t C, typename ME>
    MatrixCSR& operator*=(const MatrixExpression<T, R, C, ME>& m);

    MatrixCSR& operator/=(const T& s);

    T operator()(size_t i, size_t j) const;

    bool operator==(const MatrixCSR& m) const;

    bool operator!=(const MatrixCSR& m) const;

    static MatrixCSR<T> MakeIdentity(size_t m);

 private:
    [[nodiscard]] size_t HasElement(size_t i, size_t j) const;

    template <typename Op>
    MatrixCSR BinaryOp(const MatrixCSR& m, Op op) const;

    Vector2UZ m_size;
    NonZeroContainerType m_nonZeros;
    IndexContainerType m_rowPointers;
    IndexContainerType m_columnIndices;
};

typedef MatrixCSR<float> MatrixCSRF;

typedef MatrixCSR<double> MatrixCSRD;
}  // namespace CubbyFlow

#include <Core/Matrix/MatrixCSR-Impl.hpp>

#endif