Single-file reduced Peanut example for blog

reduced_peanut_example.cpp

//
// This software is released under the MIT license.
//
// Copyright (c) 2022 Jino Park (http://github.com/pjessesco)
//
// Permission is hereby granted, free of charge, to any person obtaining a copy
// of this software and associated documentation files (the "Software"), to deal
// in the Software without restriction, including without limitation the rights
// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the Software is
// furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included in all
// copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
// SOFTWARE.
//

// Reduced Peanut library (https://github.com/pjessesco/peanut) for blog post.
// https://pjessesco.tistory.com/78 (Korean)
// https://pjessesco.tistory.com/79 (Korean)
// https://pjessesco.tistory.com/81 (Korean)

#include <iostream>
#include <array>
#include <chrono>

using Index = unsigned int;

template <typename E>
struct MatrixExpr{
    inline auto elem(Index r, Index c) const{
        return static_cast<const E&>(*this).elem(r, c);
    }

    inline auto& elem(Index r, Index c) {
        return static_cast<const E&>(*this).elem(r, c);
    }

    static constexpr Index row = E::row;
    static constexpr Index col = E::col;
};

template<typename T, Index Row, Index Col>
    requires std::is_arithmetic_v<T> && (Row > 0) && (Col > 0)
struct Matrix : public MatrixExpr<Matrix<T, Row, Col>>{
    using Type = T;
    static constexpr Index row = Row;
    static constexpr Index col = Col;

    template <typename ...TList> requires std::conjunction_v<std::is_same<T, TList>...>
    Matrix(TList ... tlist) : m_data{std::forward<T>(tlist)...} {}

    template<typename E>
    Matrix(const MatrixExpr<E> &expr) {
        for(Index r=0;r<Row;r++){
            for(Index c=0;c<Col;c++){
                elem(r, c) = expr.elem(r, c);
            }
        }
    }

    inline T elem(Index r, Index c) const{
        return m_data[Col*r+c];
    }

    inline T& elem(Index r, Index c) {
        return m_data[Col*r+c];
    }

    std::array<T, Row*Col> m_data;
};

template<typename E1, typename E2>
struct MatrixSum : public MatrixExpr<MatrixSum<E1, E2>> {
    MatrixSum(const E1 &x, const E2 &y) : x{x}, y{y} {}

    // Static polymorphism implementation of MatrixExpr
    inline auto elem(Index r, Index c) const {
        return x.elem(r, c) + y.elem(r, c);
    }

    using Type = typename E1::Type;
    static constexpr Index row = E1::row;
    static constexpr Index col = E1::col;

    const E1 &x;
    const E2 &y;
};

template<typename E1, typename E2>
MatrixSum<E1, E2> operator+(const MatrixExpr<E1> &x, const MatrixExpr<E2> &y) {
    return MatrixSum<E1, E2>(static_cast<const E1 &>(x), static_cast<const E2 &>(y));
}

template<typename E1, typename E2>
    requires(E1::col == E2::row)
struct MatrixMult : public MatrixExpr<MatrixMult<E1, E2>> {
    using Type = typename E1::Type;
    static constexpr Index row = E1::row;
    static constexpr Index col = E2::col;

    MatrixMult(const E1 &_x, const E2 &_y) : x{_x}, y{_y} {}

    // Static polymorphism implementation of MatrixExpr
    inline auto elem(Index r, Index c) const {
        Type ret = static_cast<Type>(0);
        for (Index i = 0; i < E1::col; i++) {
            ret += x.elem(r, i) * y.elem(i, c);
        }
        return ret;
    }

//#define USE_NAIVE_MATRIX_MULT

#if defined USE_NAIVE_MATRIX_MULT
    const E1 &x;
    const E2 &y;
#else
    const Matrix<Type, E1::row, E1::col> x;
    const Matrix<Type, E2::row, E2::col> y;
#endif
};

template<typename E1, typename E2>
    requires(E1::col == E2::row)
MatrixMult<E1, E2> operator*(const MatrixExpr<E1> &x, const MatrixExpr<E2> &y) {
    return MatrixMult<E1, E2>(static_cast<const E1 &>(x), static_cast<const E2 &>(y));
}

int main() {

    Matrix<int,5,5> mat1{1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5,1,2,3,4,5};
    Matrix<int,5,5> mat2{5,4,3,2,1,5,4,3,2,1,5,4,3,2,1,5,4,3,2,1,5,4,3,2,1};

    auto time1 = std::chrono::high_resolution_clock::now();
    for(int i=0;i<100000;i++){
        Matrix<int,2,2> mat = mat1*mat2*mat1*mat2*mat1*mat2*mat1*mat2;
    }
    auto time2 = std::chrono::high_resolution_clock::now();
    std::cout<<std::chrono::duration_cast<std::chrono::microseconds>(time2 - time1).count() << " [micro-s]"<<std::endl;

    return 0;
}