Compare commits

1 Commits

Author SHA1 Message Date
fc61442b68 Made my own equally wrong QR factorization
Some checks failed
Merge-Checker / build_and_test (pull_request) Failing after 20s
2025-06-02 21:19:17 -04:00
8 changed files with 308 additions and 334 deletions

View File

@@ -76,7 +76,7 @@
"clangd.enable": true,
"C_Cpp.dimInactiveRegions": false,
"editor.defaultFormatter": "xaver.clang-format",
"clangd.inactiveRegions.useBackgroundHighlight": false,
"clangd.inactiveRegions.useBackgroundHighlight": true,
"clangd.arguments": [
"--compile-commands-dir=${workspaceFolder}/build"
],

View File

@@ -2,11 +2,8 @@
This matrix math library is focused on embedded development and avoids any heap memory allocation unless you explicitly ask for it.
It uses templates to pre-allocate matrices on the stack.
# Building
1. Initialize the repositiory with the command:
```bash
cmake -S . -B build -G Ninja
```
2. Go into the build folder and run `ninja`
3. That's it. You can test out the build by running `./unit-tests/matrix-tests`
There are still several operations that are works in progress such as:
- Add a function to calculate eigenvalues/vectors
- Add a function to compute RREF
- Add a function for SVD decomposition
- Add a function for LQ decomposition

View File

@@ -1,10 +1,3 @@
// This #ifndef section makes clangd happy so that it can properly do type hints
// in this file
#ifndef MATRIX_H_
#define MATRIX_H_
#include "Matrix.hpp"
#endif
#ifdef MATRIX_H_ // since the .cpp file has to be included by the .hpp file this
// will evaluate to true
#include "Matrix.hpp"
@@ -13,6 +6,12 @@
#include <cmath>
#include <cstdlib>
#include <cstring>
#include <type_traits>
template <uint8_t rows, uint8_t columns>
Matrix<rows, columns>::Matrix(float value) {
this->Fill(value);
}
template <uint8_t rows, uint8_t columns>
Matrix<rows, columns>::Matrix(const std::array<float, rows * columns> &array) {
@@ -26,14 +25,6 @@ Matrix<rows, columns>::Matrix(Args... args) {
static_cast<uint16_t>(columns)};
std::initializer_list<float> initList{static_cast<float>(args)...};
// if there is only one value, we actually want to do a fill
if (sizeof...(args) == 1) {
this->Fill(*initList.begin());
}
static_assert(sizeof...(args) == arraySize || sizeof...(args) == 1,
"You did not provide the right amount of initializers for this "
"matrix size");
// choose whichever buffer size is smaller for the copy length
uint32_t minSize =
std::min(arraySize, static_cast<uint16_t>(initList.size()));
@@ -41,13 +32,11 @@ Matrix<rows, columns>::Matrix(Args... args) {
}
template <uint8_t rows, uint8_t columns>
Matrix<rows, columns> Matrix<rows, columns>::Identity() {
Matrix<rows, columns> identityMatrix{0};
uint32_t minDimension = std::min(rows, columns);
for (uint8_t idx{0}; idx < minDimension; idx++) {
identityMatrix[idx][idx] = 1;
void Matrix<rows, columns>::Identity() {
this->Fill(0);
for (uint8_t idx{0}; idx < rows; idx++) {
this->matrix[idx * columns + idx] = 1;
}
return identityMatrix;
}
template <uint8_t rows, uint8_t columns>
@@ -568,19 +557,16 @@ void Matrix<rows, columns>::EigenQR(Matrix<rows, rows> &eigenVectors,
uint32_t maxIterations,
float tolerance) const {
static_assert(rows > 1, "Matrix size must be > 1 for QR iteration");
static_assert(rows == columns, "Matrix size must be square for QR iteration");
Matrix<rows, rows> Ak = *this; // Copy original matrix
Matrix<rows, rows> QQ{Matrix<rows, rows>::Identity()};
Matrix<rows, rows> shift{0};
Matrix<rows, rows> QQ{};
QQ.Identity();
for (uint32_t iter = 0; iter < maxIterations; ++iter) {
Matrix<rows, rows> Q, R;
Ak.QRDecomposition(Q, R);
// // QR shift lets us "attack" the first diagonal to speed up the algorithm
// shift = Matrix<rows, rows>::Identity() * Ak[rows - 1][rows - 1];
(Ak - shift).QRDecomposition(Q, R);
Ak = R * Q + shift;
Ak = R * Q;
QQ = QQ * Q;
// Check convergence: off-diagonal norm

View File

@@ -1,4 +1,5 @@
#pragma once
#ifndef MATRIX_H_
#define MATRIX_H_
#include <array>
#include <cstdint>
@@ -18,6 +19,11 @@ public:
*/
Matrix() = default;
/**
* @brief Create a matrix but fill all of its entries with one value
*/
Matrix(float value);
/**
* @brief Initialize a matrix with an array
*/
@@ -34,9 +40,9 @@ public:
template <typename... Args> Matrix(Args... args);
/**
* @brief Create an identity matrix
* @brief set the matrix diagonals to 1 and all other values to 0
*/
static Matrix<rows, columns> Identity();
void Identity();
/**
* @brief Set all elements in this to value
@@ -252,6 +258,6 @@ private:
void setMatrixToArray(const std::array<float, rows * columns> &array);
};
#ifndef MATRIX_H_
#include "Matrix.cpp"
#endif // MATRIX_H_

View File

@@ -2,11 +2,12 @@
#define QUATERNION_H_
#include "Matrix.hpp"
class Quaternion : public Matrix<1, 4> {
class Quaternion : public Matrix<1, 4>
{
public:
Quaternion() : Matrix<1, 4>() {}
Quaternion(float w, float v1, float v2, float v3)
: Matrix<1, 4>(w, v1, v2, v3) {}
Quaternion(float fillValue) : Matrix<1, 4>(fillValue) {}
Quaternion(float w, float v1, float v2, float v3) : Matrix<1, 4>(w, v1, v2, v3) {}
Quaternion(const Quaternion &q) : Matrix<1, 4>(q.w, q.v1, q.v2, q.v3) {}
Quaternion(const Matrix<1, 4> &matrix) : Matrix<1, 4>(matrix) {}
Quaternion(const std::array<float, 4> &array) : Matrix<1, 4>(array) {}

View File

@@ -10,61 +10,41 @@
#include <cmath>
#include <iostream>
// Helper functions
template <uint8_t rows, uint8_t columns>
float matrixSum(const Matrix<rows, columns> &matrix) {
float sum = 0;
for (uint32_t i = 0; i < rows * columns; i++) {
float number = matrix.ToArray()[i];
sum += number * number;
}
return std::sqrt(sum);
}
template <uint8_t rows, uint8_t columns>
void printLabeledMatrix(const std::string &label,
const Matrix<rows, columns> &matrix) {
std::string strBuf = "";
matrix.ToString(strBuf);
std::cout << label << ":\n" << strBuf << std::endl;
}
TEST_CASE("Initialization", "Matrix") {
SECTION("Array Initialization") {
std::array<float, 4> arr2{5, 6, 7, 8};
Matrix<2, 2> mat2{arr2};
// array initialization
REQUIRE(mat2.Get(0, 0) == 5);
REQUIRE(mat2.Get(0, 1) == 6);
REQUIRE(mat2.Get(1, 0) == 7);
REQUIRE(mat2.Get(1, 1) == 8);
}
SECTION("Argument Pack Initialization") {
Matrix<2, 2> mat1{1, 2, 3, 4};
// template pack initialization
REQUIRE(mat1.Get(0, 0) == 1);
REQUIRE(mat1.Get(0, 1) == 2);
REQUIRE(mat1.Get(1, 0) == 3);
REQUIRE(mat1.Get(1, 1) == 4);
}
SECTION("Single Argument Pack Initialization") {
Matrix<2, 2> mat1{2};
// template pack initialization
REQUIRE(mat1.Get(0, 0) == 2);
REQUIRE(mat1.Get(0, 1) == 2);
REQUIRE(mat1.Get(1, 0) == 2);
REQUIRE(mat1.Get(1, 1) == 2);
}
}
TEST_CASE("Elementary Matrix Operations", "Matrix") {
std::array<float, 4> arr2{5, 6, 7, 8};
Matrix<2, 2> mat1{1, 2, 3, 4};
Matrix<2, 2> mat2{arr2};
Matrix<2, 2> mat3{};
SECTION("Initialization") {
// array initialization
REQUIRE(mat1.Get(0, 0) == 1);
REQUIRE(mat1.Get(0, 1) == 2);
REQUIRE(mat1.Get(1, 0) == 3);
REQUIRE(mat1.Get(1, 1) == 4);
// empty initialization
REQUIRE(mat3.Get(0, 0) == 0);
REQUIRE(mat3.Get(0, 1) == 0);
REQUIRE(mat3.Get(1, 0) == 0);
REQUIRE(mat3.Get(1, 1) == 0);
// template pack initialization
REQUIRE(mat2.Get(0, 0) == 5);
REQUIRE(mat2.Get(0, 1) == 6);
REQUIRE(mat2.Get(1, 0) == 7);
REQUIRE(mat2.Get(1, 1) == 8);
// large matrix
Matrix<255, 255> mat6{};
mat6.Fill(4);
for (uint8_t row{0}; row < 255; row++) {
for (uint8_t column{0}; column < 255; column++) {
REQUIRE(mat6.Get(row, column) == 4);
}
}
}
SECTION("Fill") {
mat1.Fill(0);
REQUIRE(mat1.Get(0, 0) == 0);
@@ -86,6 +66,10 @@ TEST_CASE("Elementary Matrix Operations", "Matrix") {
}
SECTION("Addition") {
std::string strBuf1 = "";
mat1.ToString(strBuf1);
std::cout << "Matrix 1:\n" << strBuf1 << std::endl;
mat1.Add(mat2, mat3);
REQUIRE(mat3.Get(0, 0) == 6);
@@ -379,58 +363,18 @@ TEST_CASE("Elementary Matrix Operations", "Matrix") {
}
}
TEST_CASE("Identity Matrix", "Matrix") {
SECTION("Square Matrix") {
Matrix<5, 5> matrix = Matrix<5, 5>::Identity();
uint32_t oneColumnIndex{0};
for (uint32_t row = 0; row < 5; row++) {
for (uint32_t column = 0; column < 5; column++) {
float value = matrix[row][column];
if (oneColumnIndex == column) {
REQUIRE_THAT(value, Catch::Matchers::WithinRel(1.0f, 1e-6f));
} else {
REQUIRE_THAT(value, Catch::Matchers::WithinRel(0.0f, 1e-6f));
}
}
oneColumnIndex++;
}
}
SECTION("Wide Matrix") {
Matrix<2, 5> matrix = Matrix<2, 5>::Identity();
uint32_t oneColumnIndex{0};
for (uint32_t row = 0; row < 2; row++) {
for (uint32_t column = 0; column < 5; column++) {
float value = matrix[row][column];
if (oneColumnIndex == column && row < 3) {
REQUIRE_THAT(value, Catch::Matchers::WithinRel(1.0f, 1e-6f));
} else {
REQUIRE_THAT(value, Catch::Matchers::WithinRel(0.0f, 1e-6f));
}
}
oneColumnIndex++;
}
}
SECTION("Tall Matrix") {
Matrix<5, 2> matrix = Matrix<5, 2>::Identity();
uint32_t oneColumnIndex{0};
for (uint32_t row = 0; row < 5; row++) {
for (uint32_t column = 0; column < 2; column++) {
float value = matrix[row][column];
if (oneColumnIndex == column) {
REQUIRE_THAT(value, Catch::Matchers::WithinRel(1.0f, 1e-6f));
} else {
REQUIRE_THAT(value, Catch::Matchers::WithinRel(0.0f, 1e-6f));
}
}
oneColumnIndex++;
}
template <uint8_t rows, uint8_t columns>
float matrixSum(const Matrix<rows, columns> &matrix) {
float sum = 0;
for (uint32_t i = 0; i < rows * columns; i++) {
float number = matrix.ToArray()[i];
sum += number * number;
}
return std::sqrt(sum);
}
// TODO: Add test for scalar division
TEST_CASE("Euclidean Norm", "Matrix") {
SECTION("2x2 Normalize") {
@@ -479,48 +423,48 @@ TEST_CASE("Euclidean Norm", "Matrix") {
}
TEST_CASE("QR Decompositions", "Matrix") {
SECTION("2x2 QRDecomposition") {
Matrix<2, 2> A{1.0f, 2.0f, 3.0f, 4.0f};
Matrix<2, 2> Q{}, R{};
A.QRDecomposition(Q, R);
// SECTION("2x2 QRDecomposition") {
// Matrix<2, 2> A{1.0f, 2.0f, 3.0f, 4.0f};
// Matrix<2, 2> Q{}, R{};
// A.QRDecomposition(Q, R);
// Check that Q * R ≈ A
Matrix<2, 2> QR{};
Q.Mult(R, QR);
for (int i = 0; i < 2; ++i) {
for (int j = 0; j < 2; ++j) {
REQUIRE_THAT(QR[i][j], Catch::Matchers::WithinRel(A[i][j], 1e-4f));
}
}
// // Check that Q * R ≈ A
// Matrix<2, 2> QR{};
// Q.Mult(R, QR);
// for (int i = 0; i < 2; ++i) {
// for (int j = 0; j < 2; ++j) {
// REQUIRE_THAT(QR[i][j], Catch::Matchers::WithinRel(A[i][j], 1e-4f));
// }
// }
// Check that Q is orthonormal: Qᵀ * Q ≈ I
Matrix<2, 2> Qt = Q.Transpose();
Matrix<2, 2> QtQ{};
Qt.Mult(Q, QtQ);
for (int i = 0; i < 2; ++i) {
for (int j = 0; j < 2; ++j) {
if (i == j)
REQUIRE_THAT(QtQ[i][j], Catch::Matchers::WithinRel(1.0f, 1e-4f));
else
REQUIRE_THAT(QtQ[i][j], Catch::Matchers::WithinAbs(0.0f, 1e-4f));
}
}
// // Check that Q is orthonormal: Qᵀ * Q ≈ I
// Matrix<2, 2> Qt = Q.Transpose();
// Matrix<2, 2> QtQ{};
// Qt.Mult(Q, QtQ);
// for (int i = 0; i < 2; ++i) {
// for (int j = 0; j < 2; ++j) {
// if (i == j)
// REQUIRE_THAT(QtQ[i][j], Catch::Matchers::WithinRel(1.0f, 1e-4f));
// else
// REQUIRE_THAT(QtQ[i][j], Catch::Matchers::WithinAbs(0.0f, 1e-4f));
// }
// }
// Optional: R should be upper triangular
REQUIRE(std::fabs(R[1][0]) < 1e-4f);
// // Optional: R should be upper triangular
// REQUIRE(std::fabs(R[1][0]) < 1e-4f);
// check that all Q values are correct
REQUIRE_THAT(Q[0][0], Catch::Matchers::WithinRel(0.3162f, 1e-4f));
REQUIRE_THAT(Q[0][1], Catch::Matchers::WithinRel(0.94868f, 1e-4f));
REQUIRE_THAT(Q[1][0], Catch::Matchers::WithinRel(0.94868f, 1e-4f));
REQUIRE_THAT(Q[1][1], Catch::Matchers::WithinRel(-0.3162f, 1e-4f));
// // check that all Q values are correct
// REQUIRE_THAT(Q[0][0], Catch::Matchers::WithinRel(0.3162f, 1e-4f));
// REQUIRE_THAT(Q[0][1], Catch::Matchers::WithinRel(0.94868f, 1e-4f));
// REQUIRE_THAT(Q[1][0], Catch::Matchers::WithinRel(0.94868f, 1e-4f));
// REQUIRE_THAT(Q[1][1], Catch::Matchers::WithinRel(-0.3162f, 1e-4f));
// check that all R values are correct
REQUIRE_THAT(R[0][0], Catch::Matchers::WithinRel(3.16228f, 1e-4f));
REQUIRE_THAT(R[0][1], Catch::Matchers::WithinRel(4.42719f, 1e-4f));
REQUIRE_THAT(R[1][0], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(R[1][1], Catch::Matchers::WithinRel(0.63246f, 1e-4f));
}
// // check that all R values are correct
// REQUIRE_THAT(R[0][0], Catch::Matchers::WithinRel(3.16228f, 1e-4f));
// REQUIRE_THAT(R[0][1], Catch::Matchers::WithinRel(4.42719f, 1e-4f));
// REQUIRE_THAT(R[1][0], Catch::Matchers::WithinRel(0.0f, 1e-4f));
// REQUIRE_THAT(R[1][1], Catch::Matchers::WithinRel(0.63246f, 1e-4f));
// }
SECTION("3x3 QRDecomposition") {
// this symmetrix tridiagonal matrix is well behaved for testing
@@ -529,6 +473,13 @@ TEST_CASE("QR Decompositions", "Matrix") {
Matrix<3, 3> Q{}, R{};
A.QRDecomposition(Q, R);
std::string strBuf1 = "";
Q.ToString(strBuf1);
std::cout << "Q:\n" << strBuf1 << std::endl;
strBuf1 = "";
R.ToString(strBuf1);
std::cout << "R:\n" << strBuf1 << std::endl;
// Check that Q * R ≈ A
Matrix<3, 3> QR{};
QR = Q * R;
@@ -539,13 +490,11 @@ TEST_CASE("QR Decompositions", "Matrix") {
}
// Check that Qᵀ * Q ≈ I
// Since the rank of this matrix is 2, only the top left 2x2 sub-matrix will
// equal I.
Matrix<3, 3> Qt = Q.Transpose();
Matrix<3, 3> QtQ{};
QtQ = Qt * Q;
for (int i = 0; i < 2; ++i) {
for (int j = 0; j < 2; ++j) {
for (int i = 0; i < 3; ++i) {
for (int j = 0; j < 3; ++j) {
if (i == j)
REQUIRE_THAT(QtQ[i][j], Catch::Matchers::WithinRel(1.0f, 1e-4f));
else
@@ -559,6 +508,28 @@ TEST_CASE("QR Decompositions", "Matrix") {
REQUIRE(std::fabs(R[i][j]) < 1e-4f);
}
}
// check that all Q values are correct
REQUIRE_THAT(Q[0][0], Catch::Matchers::WithinRel(0.1231f, 1e-4f));
REQUIRE_THAT(Q[0][1], Catch::Matchers::WithinRel(0.904534f, 1e-4f));
REQUIRE_THAT(Q[0][2], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(Q[1][0], Catch::Matchers::WithinRel(0.49237f, 1e-4f));
REQUIRE_THAT(Q[1][1], Catch::Matchers::WithinRel(0.301511f, 1e-4f));
REQUIRE_THAT(Q[1][2], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(Q[2][0], Catch::Matchers::WithinRel(0.86164f, 1e-4f));
REQUIRE_THAT(Q[2][1], Catch::Matchers::WithinRel(-0.30151f, 1e-4f));
REQUIRE_THAT(Q[2][2], Catch::Matchers::WithinRel(0.0f, 1e-4f));
// check that all R values are correct
REQUIRE_THAT(R[0][0], Catch::Matchers::WithinRel(8.124038f, 1e-4f));
REQUIRE_THAT(R[0][1], Catch::Matchers::WithinRel(9.60114f, 1e-4f));
REQUIRE_THAT(R[0][2], Catch::Matchers::WithinRel(11.07823f, 1e-4f));
REQUIRE_THAT(R[1][0], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(R[1][1], Catch::Matchers::WithinRel(0.90453f, 1e-4f));
REQUIRE_THAT(R[1][2], Catch::Matchers::WithinRel(1.80907f, 1e-4f));
REQUIRE_THAT(R[2][0], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(R[2][1], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(R[2][2], Catch::Matchers::WithinRel(1.0f, 1e-4f));
}
SECTION("4x2 QRDecomposition") {
@@ -600,41 +571,42 @@ TEST_CASE("QR Decompositions", "Matrix") {
}
}
TEST_CASE("Eigenvalues and Vectors", "Matrix") {
SECTION("2x2 Eigen") {
Matrix<2, 2> A{1.0f, 2.0f, 3.0f, 4.0f};
Matrix<2, 2> vectors{};
Matrix<2, 1> values{};
// TEST_CASE("Eigenvalues and Vectors", "Matrix") {
// SECTION("2x2 Eigen") {
// Matrix<2, 2> A{1.0f, 2.0f, 3.0f, 4.0f};
// Matrix<2, 2> vectors{};
// Matrix<2, 1> values{};
A.EigenQR(vectors, values, 1000000, 1e-20f);
// A.EigenQR(vectors, values, 1000000, 1e-20f);
REQUIRE_THAT(vectors[0][0], Catch::Matchers::WithinRel(0.41597f, 1e-4f));
REQUIRE_THAT(vectors[1][0], Catch::Matchers::WithinRel(0.90938f, 1e-4f));
REQUIRE_THAT(values[0][0], Catch::Matchers::WithinRel(5.372282f, 1e-4f));
REQUIRE_THAT(values[1][0], Catch::Matchers::WithinRel(-0.372281f, 1e-4f));
}
// REQUIRE_THAT(vectors[0][0], Catch::Matchers::WithinRel(0.41597f, 1e-4f));
// REQUIRE_THAT(vectors[1][0], Catch::Matchers::WithinRel(0.90938f, 1e-4f));
// REQUIRE_THAT(values[0][0], Catch::Matchers::WithinRel(5.372282f, 1e-4f));
// REQUIRE_THAT(values[1][0], Catch::Matchers::WithinRel(-0.372281f,
// 1e-4f));
// }
SECTION("3x3 Rank Defficient Eigen") {
SKIP("Skipping this because QR decomposition isn't ready for it");
// this symmetrix tridiagonal matrix is well behaved for testing
Matrix<3, 3> A{1, 2, 3, 4, 5, 6, 7, 8, 9};
// SECTION("3x3 Eigen") {
// // this symmetrix tridiagonal matrix is well behaved for testing
// Matrix<3, 3> A{1, 2, 3, 4, 5, 6, 7, 8, 9};
Matrix<3, 3> vectors{};
Matrix<3, 1> values{};
A.EigenQR(vectors, values, 1000000, 1e-8f);
// Matrix<3, 3> vectors{};
// Matrix<3, 1> values{};
// A.EigenQR(vectors, values, 1000000, 1e-8f);
std::string strBuf1 = "";
vectors.ToString(strBuf1);
std::cout << "Vectors:\n" << strBuf1 << std::endl;
strBuf1 = "";
values.ToString(strBuf1);
std::cout << "Values:\n" << strBuf1 << std::endl;
// std::string strBuf1 = "";
// vectors.ToString(strBuf1);
// std::cout << "Vectors:\n" << strBuf1 << std::endl;
// strBuf1 = "";
// values.ToString(strBuf1);
// std::cout << "Values:\n" << strBuf1 << std::endl;
REQUIRE_THAT(vectors[0][0], Catch::Matchers::WithinRel(0.23197f, 1e-4f));
REQUIRE_THAT(vectors[1][0], Catch::Matchers::WithinRel(0.525322f, 1e-4f));
REQUIRE_THAT(vectors[2][0], Catch::Matchers::WithinRel(0.81867f, 1e-4f));
REQUIRE_THAT(values[0][0], Catch::Matchers::WithinRel(-1.11684f, 1e-4f));
REQUIRE_THAT(values[1][0], Catch::Matchers::WithinRel(0.0f, 1e-4f));
REQUIRE_THAT(values[2][0], Catch::Matchers::WithinRel(16.1168f, 1e-4f));
}
}
// REQUIRE_THAT(vectors[0][0], Catch::Matchers::WithinRel(0.23197f, 1e-4f));
// REQUIRE_THAT(vectors[1][0], Catch::Matchers::WithinRel(0.525322f,
// 1e-4f)); REQUIRE_THAT(vectors[2][0], Catch::Matchers::WithinRel(0.81867f,
// 1e-4f)); REQUIRE_THAT(values[0][0], Catch::Matchers::WithinRel(-1.11684f,
// 1e-4f)); REQUIRE_THAT(values[1][0], Catch::Matchers::WithinRel(0.0f,
// 1e-4f)); REQUIRE_THAT(values[2][0], Catch::Matchers::WithinRel(16.1168f,
// 1e-4f));
// }
// }

View File

@@ -8,7 +8,6 @@
// any other libraries
#include <array>
#include <cmath>
#include <cstdint>
// basically re-run all of the matrix tests with huge matrices and time the
// results.
@@ -30,13 +29,13 @@ TEST_CASE("Timing Tests", "Matrix") {
Matrix<4, 4> mat5{};
SECTION("Addition") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat3 = mat1 + mat2;
}
}
SECTION("Subtraction") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat3 = mat1 - mat2;
}
}
@@ -48,19 +47,19 @@ TEST_CASE("Timing Tests", "Matrix") {
}
SECTION("Scalar Multiplication") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat3 = mat1 * 3;
}
}
SECTION("Element Multiply") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat1.ElementMultiply(mat2, mat3);
}
}
SECTION("Element Divide") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat1.ElementDivide(mat2, mat3);
}
}
@@ -69,59 +68,52 @@ TEST_CASE("Timing Tests", "Matrix") {
// what about matrices of 0,0 or 1,1?
// minor matrix for 2x2 matrix
Matrix<49, 49> minorMat1{};
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat1.MinorMatrix(minorMat1, 0, 0);
}
}
SECTION("Determinant") {
for (uint32_t i{0}; i < 1000000; i++) {
for (uint32_t i{0}; i < 100000; i++) {
float det1 = mat4.Det();
}
}
SECTION("Matrix of Minors") {
for (uint32_t i{0}; i < 1000000; i++) {
for (uint32_t i{0}; i < 100000; i++) {
mat4.MatrixOfMinors(mat5);
}
}
SECTION("Invert") {
for (uint32_t i{0}; i < 1000000; i++) {
for (uint32_t i{0}; i < 100000; i++) {
mat5 = mat4.Invert();
}
};
SECTION("Transpose") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat3 = mat1.Transpose();
}
}
SECTION("Normalize") {
for (uint32_t i{0}; i < 100000; i++) {
for (uint32_t i{0}; i < 10000; i++) {
mat3 = mat1 / mat1.EuclideanNorm();
}
}
SECTION("GET ROW") {
Matrix<1, 50> mat1Rows{};
for (uint32_t i{0}; i < 100000000; i++) {
for (uint32_t i{0}; i < 1000000; i++) {
mat1.GetRow(0, mat1Rows);
}
}
SECTION("GET COLUMN") {
Matrix<50, 1> mat1Columns{};
for (uint32_t i{0}; i < 100000000; i++) {
for (uint32_t i{0}; i < 1000000; i++) {
mat1.GetColumn(0, mat1Columns);
}
}
SECTION("QR Decomposition") {
Matrix<50, 50> Q, R{};
for (uint32_t i{0}; i < 500; i++) {
mat1.QRDecomposition(Q, R);
}
}
}

View File

@@ -1,36 +1,56 @@
Running matrix-timing-tests with timing
Randomness seeded to: 3567651885
1.857 s: Addition
1.857 s: Timing Tests
1.788 s: Subtraction
1.788 s: Timing Tests
1.929 s: Multiplication
1.929 s: Timing Tests
1.268 s: Scalar Multiplication
1.268 s: Timing Tests
1.798 s: Element Multiply
1.798 s: Timing Tests
1.802 s: Element Divide
1.803 s: Timing Tests
1.553 s: Minor Matrix
1.554 s: Timing Tests
1.009 s: Determinant
1.009 s: Timing Tests
4.076 s: Matrix of Minors
4.076 s: Timing Tests
1.066 s: Invert
1.066 s: Timing Tests
1.246 s: Transpose
1.246 s: Timing Tests
2.284 s: Normalize
2.284 s: Timing Tests
0.606 s: GET ROW
0.606 s: Timing Tests
24.629 s: GET COLUMN
24.630 s: Timing Tests
3.064 s: QR Decomposition
3.064 s: Timing Tests
Randomness seeded to: 2444679151
0.180 s: Addition
0.180 s: Timing Tests
0.177 s: Subtraction
0.177 s: Timing Tests
1.868 s: Multiplication
1.868 s: Timing Tests
0.127 s: Scalar Multiplication
0.127 s: Timing Tests
0.173 s: Element Multiply
0.173 s: Timing Tests
0.178 s: Element Divide
0.178 s: Timing Tests
0.172 s: Minor Matrix
0.172 s: Timing Tests
0.103 s: Determinant
0.103 s: Timing Tests
0.411 s: Matrix of Minors
0.411 s: Timing Tests
0.109 s: Invert
0.109 s: Timing Tests
0.122 s: Transpose
0.122 s: Timing Tests
0.190 s: Normalize
0.190 s: Timing Tests
0.006 s: GET ROW
0.006 s: Timing Tests
0.235 s: GET COLUMN
0.235 s: Timing Tests
===============================================================================
test cases: 1 | 1 passed
assertions: - none -
Command being timed: "build/unit-tests/matrix-timing-tests -d yes"
User time (seconds): 4.05
System time (seconds): 0.00
Percent of CPU this job got: 100%
Elapsed (wall clock) time (h:mm:ss or m:ss): 0:04.05
Average shared text size (kbytes): 0
Average unshared data size (kbytes): 0
Average stack size (kbytes): 0
Average total size (kbytes): 0
Maximum resident set size (kbytes): 3200
Average resident set size (kbytes): 0
Major (requiring I/O) page faults: 184
Minor (reclaiming a frame) page faults: 171
Voluntary context switches: 1
Involuntary context switches: 26
Swaps: 0
File system inputs: 12
File system outputs: 1
Socket messages sent: 0
Socket messages received: 0
Signals delivered: 0
Page size (bytes): 4096
Exit status: 0