Skip to content

perf(math): build matrix results without zero-filling them first #349

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
gabrielfrasantos merged 2 commits into from
Oct 5, 2026
Merged

perf(math): build matrix results without zero-filling them first #349

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
e5ea80d
perf(math): build matrix results without zero-filling them first
claude Oct 5, 2026
b54894a
fix(math): build the matrix product's result in a member function
claude Oct 5, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
63 changes: 43 additions & 20 deletions numerical/math/Matrix.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,15 +63,15 @@ namespace math

[[nodiscard]] OPTIMIZE_FOR_SPEED friend constexpr Matrix operator+(const Matrix& lhs, const Matrix& rhs)
{
Matrix result;
Matrix result{ Uninitialized{} };
for (size_type i = 0; i < size; ++i)
result.data[i] = lhs.data[i] + rhs.data[i];
return result;
}

[[nodiscard]] OPTIMIZE_FOR_SPEED friend constexpr Matrix operator-(const Matrix& lhs, const Matrix& rhs)
{
Matrix result;
Matrix result{ Uninitialized{} };
for (size_type i = 0; i < size; ++i)
result.data[i] = lhs.data[i] - rhs.data[i];
return result;
Expand All @@ -80,23 +80,12 @@ namespace math
template<size_t RhsCols>
[[nodiscard]] OPTIMIZE_FOR_SPEED friend constexpr Matrix<T, Rows, RhsCols> operator*(const Matrix& lhs, const Matrix<T, Cols, RhsCols>& rhs)
{
Matrix<T, Rows, RhsCols> result;
for (size_type i = 0; i < Rows; ++i)
{
for (size_type j = 0; j < RhsCols; ++j)
{
T sum{};
for (size_type k = 0; k < Cols; ++k)
sum += lhs.at(i, k) * rhs.at(k, j);
result.at(i, j) = sum;
}
}
return result;
return lhs.Multiply(rhs);
}

[[nodiscard]] OPTIMIZE_FOR_SPEED friend constexpr Matrix operator*(const Matrix& lhs, const T& scalar)
{
Matrix result;
Matrix result{ Uninitialized{} };
for (size_type i = 0; i < size; ++i)
result.data[i] = lhs.data[i] * scalar;
return result;
Expand All @@ -121,7 +110,19 @@ namespace math
[[nodiscard]] constexpr Matrix<T, Rows, 1> GetColumn(size_type col) const;

private:
std::array<T, Rows * Cols> data = {};
template<typename, size_t, size_t>
friend class Matrix;

struct Uninitialized
{};

constexpr explicit Matrix(Uninitialized) noexcept
{}

template<size_t RhsCols>
[[nodiscard]] constexpr Matrix<T, Rows, RhsCols> Multiply(const Matrix<T, Cols, RhsCols>& rhs) const;

std::array<T, Rows * Cols> data;
};

template<typename T, typename... U>
Expand All @@ -141,10 +142,13 @@ namespace math
}

template<typename T, size_t Rows, size_t Cols>
constexpr Matrix<T, Rows, Cols>::Matrix() noexcept = default;
constexpr Matrix<T, Rows, Cols>::Matrix() noexcept
: data{}
{}

template<typename T, size_t Rows, size_t Cols>
OPTIMIZE_FOR_SPEED constexpr Matrix<T, Rows, Cols>::Matrix(std::initializer_list<std::initializer_list<T>> init)
: data{}
{
size_t row = 0;
for (const auto& row_list : init)
Expand Down Expand Up @@ -263,11 +267,30 @@ namespace math
return *this;
}

template<typename T, size_t Rows, size_t Cols>
template<size_t RhsCols>
OPTIMIZE_FOR_SPEED constexpr Matrix<T, Rows, RhsCols>
Matrix<T, Rows, Cols>::Multiply(const Matrix<T, Cols, RhsCols>& rhs) const
{
Matrix<T, Rows, RhsCols> result{ typename Matrix<T, Rows, RhsCols>::Uninitialized{} };
for (size_type i = 0; i < Rows; ++i)
{
for (size_type j = 0; j < RhsCols; ++j)
{
T sum{};
for (size_type k = 0; k < Cols; ++k)
sum += at(i, k) * rhs.at(k, j);
result.at(i, j) = sum;
}
}
return result;
}

template<typename T, size_t Rows, size_t Cols>
OPTIMIZE_FOR_SPEED constexpr Matrix<T, Cols, Rows>
Matrix<T, Rows, Cols>::Transpose() const
{
Matrix<T, Cols, Rows> result;
Matrix<T, Cols, Rows> result{ typename Matrix<T, Cols, Rows>::Uninitialized{} };
for (size_type i = 0; i < Rows; ++i)
for (size_type j = 0; j < Cols; ++j)
result.at(j, i) = at(i, j);
Expand Down Expand Up @@ -328,7 +351,7 @@ namespace math
{
static_assert(BlockRows <= Rows && BlockCols <= Cols,
"Requested block exceeds source matrix dimensions");
Matrix<T, BlockRows, BlockCols> result;
Matrix<T, BlockRows, BlockCols> result{ typename Matrix<T, BlockRows, BlockCols>::Uninitialized{} };
for (size_type r = 0; r < BlockRows; ++r)
for (size_type c = 0; c < BlockCols; ++c)
result.at(r, c) = at(rowOffset + r, colOffset + c);
Expand All @@ -339,7 +362,7 @@ namespace math
[[nodiscard]] OPTIMIZE_FOR_SPEED constexpr Matrix<T, Rows, 1>
Matrix<T, Rows, Cols>::GetColumn(size_type col) const
{
Vector<T, Rows> result;
Vector<T, Rows> result{ typename Vector<T, Rows>::Uninitialized{} };
for (size_type r = 0; r < Rows; ++r)
result.at(r, 0) = at(r, col);
return result;
Expand Down
18 changes: 18 additions & 0 deletions numerical/math/test/TestMatrix.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -427,6 +427,24 @@ TEST_F(MatrixFloatShapeTest, three_by_three_compound_assignment_and_trace)
EXPECT_NEAR(square3.Transpose().at(0, 2), 7.0f, math::Tolerance<float>());
}

TEST_F(MatrixFloatShapeTest, operators_write_every_element_of_their_result)
{
constexpr math::Matrix<float, 3, 1> column{ { 1.0f }, { 2.0f }, { 3.0f } };
constexpr auto outer = column * column.Transpose();
constexpr auto sum = outer + outer;
constexpr auto difference = sum - outer;
constexpr auto scaled = difference * 0.5f;
constexpr auto block = scaled.GetBlock<2, 2>(1, 1);
constexpr auto lastColumn = scaled.GetColumn(2);

EXPECT_NEAR(outer.at(2, 1), 6.0f, math::Tolerance<float>());
EXPECT_NEAR(sum.at(0, 2), 6.0f, math::Tolerance<float>());
EXPECT_NEAR(difference.at(1, 1), 4.0f, math::Tolerance<float>());
EXPECT_NEAR(scaled.at(2, 2), 4.5f, math::Tolerance<float>());
EXPECT_NEAR(block.at(1, 0), 3.0f, math::Tolerance<float>());
EXPECT_NEAR(lastColumn.at(0, 0), 1.5f, math::Tolerance<float>());
}

TEST_F(MatrixFloatShapeTest, non_square_shapes_transpose_and_index)
{
math::Matrix<float, 1, 4> row{ { 1.0f, 2.0f, 3.0f, 4.0f } };
Expand Down
Loading