matrix.hpp

#include <array>
#include <cassert>
#pragma once
template <size_t M, size_t N, typename T = float>
struct Matrix
{
    std::array<std::array<T, N>, M> lattice;

    constexpr Matrix sum(const Matrix<M, N, T> &rhs) const
    {
        Matrix m = *this;
        for (size_t i = 0; i < M; i++)
        {
            for (size_t j = 0; j < N; j++)
            {
                m.lattice[i][j] += rhs.lattice[i][j];
            }
        }
        return m;
    }

    template <size_t R>
    constexpr Matrix product(const Matrix<N, R, T> &rhs) const
    {
        Matrix m = *this;
        for (size_t i = 0; i < M; i++)
        {
            for (size_t j = 0; j < N; j++)
            {
                T counter = 0;
                for (size_t k = 0; k < N; k++)
                {
                    counter += m.lattice[i][k] * rhs.lattice[k][j];
                }
                m.lattice[i][j] = counter;
            }
        }
        return m;
    }

    constexpr Matrix scalar_product(const T &mul) const
    {
        Matrix m = *this;
        for (size_t i = 0; i < M; i++)
        {
            for (size_t j = 0; j < N; j++)
            {
                m.lattice[i][j] *= mul;
            }
        }
        return m;
    }

    constexpr Matrix swap_rows(size_t i, size_t j) const{
        assert(i <= M && j <= M);
        Matrix m = *this;
        if(i == j) return m;
        for(size_t k = 0; k < N; k++){
            T tmp = m.lattice[i][k];
            m.lattice[i][k] = m.lattice[j][k];  // std::swap is not constexpr, thanks committee
            m.lattice[j][k] = tmp;
        }
        return m;
    }

    constexpr static Matrix Zero()
    {
        return Matrix{};
    }
    
    constexpr static Matrix Identity()
    {
        static_assert(M == N, "Identity matrix must be square.");
        Matrix m{};
        for (size_t i = 0; i < N; i++)
        {
            m.lattice[i][i] = 1;
        }
        return m;
    }


};