4

I'm trying to write a program to find the determinant of an NxN matrix, and in doing so, learning about std::array and templates for the first time.

I came up with this implementation, based on the method of cofactors and minor.

#include <array>
#include <iostream>

template <int N>
using Matrix = std::array<std::array<double, N>, N>;

template <int N>
Matrix<N - 1> subMatrix(Matrix<N> matrix, size_t focusRow, size_t focusColumn) {
    Matrix<N - 1> returnMatrix;

    int subMatrixRow = 0, subMatrixColumn = 0;
    static const int matrixSize = matrix.size();

    for (size_t matrixRow = 0; matrixRow < matrixSize; matrixRow++) {
        for (size_t matrixColumn = 0; matrixColumn < matrixSize; matrixColumn++) {
            if (matrixRow != focusRow && matrixColumn != focusColumn) {
                returnMatrix[subMatrixRow][subMatrixColumn++] = matrix[matrixRow][matrixColumn];

                if (subMatrixColumn == matrixSize - 1) {
                    subMatrixColumn = 0;
                    subMatrixRow++;
                }
            }
        }
    }

    return returnMatrix;
}

template <int N>
double getDeterminant(Matrix<N> matrix) {
    static const int matrixSize = matrix.size();
    double determinant = 0;

    if (matrixSize == 1) {
        determinant = matrix[0][0];
    }

    else if (matrixSize == 2) {
        determinant = (matrix[0][0] * matrix[1][1]) - (matrix[0][1] * matrix[1][0]);
    }

    else if (matrixSize > 0) {
        int sign = 1;

        for (size_t dimension = 0; dimension < N; dimension++) {
            determinant += sign * matrix[0][dimension] * getDeterminant(subMatrix(matrix, 0, dimension));
            sign = -sign;
        }
    }

    else {
        throw std::invalid_argument("expected square matrix");
    }

    return determinant;
}

int main(int argc, char const* argv[]) {
    static const int length = 3;

    Matrix<length> matrix = {{{1, 2, 3},
                              {4, 5, 6},
                              {7, 8, 9}}};

    printf("determinant = %.2f\n", getDeterminant(matrix));
}

However, when I run make, I get an error I'm not quite able to resolve. How to proceed?

vscode ➜ /workspaces/c-cpp-mirror/11032023-arrays $ make determinant
g++     determinant.cpp   -o determinant
determinant.cpp: In function ‘int main(int, const char**)’:
determinant.cpp:68:57: error: no matching function for call to ‘getDeterminant(Matrix<3>&)’
   68 |     printf("determinant = %.2f\n", getDeterminant(matrix));
      |                                                         ^
determinant.cpp:31:8: note: candidate: ‘template<int N> double getDeterminant(Matrix<N>)’
   31 | double getDeterminant(Matrix<N> matrix) {
      |        ^~~~~~~~~~~~~~
determinant.cpp:31:8: note:   template argument deduction/substitution failed:
determinant.cpp:68:57: note:   mismatched types ‘int’ and ‘long unsigned int’
   68 |     printf("determinant = %.2f\n", getDeterminant(matrix));
      |                                                         ^
make: *** [<builtin>: determinant] Error 1

Here's a smaller program that still exhibits the same behaviour:

#include <array>
#include <iostream>

template <int N>
using Matrix = std::array<std::array<double, N>, N>;

template <int N>
void printMatrix(Matrix<N> matrix) {}

int main(int argc, char const* argv[]) {
    static const int length = 4;

    Matrix<length> matrix = {{{7, -2, 2, 1},
                              {3, 1, -5, 2},
                              {2, 2, -5, 3},
                              {3, 2, 5, 1}}};

    printMatrix(matrix);
}
eccentricOrange
  • 876
  • 5
  • 18
  • @FrançoisAndrieux Yeah, I just realised that after posting too. Have edited, thanks! – eccentricOrange Mar 12 '23 at 11:45
  • 1
    `std::array`'s size argument type is `std::size_t` and not `int`. Mixing the two might be messing with template argument deduction, and since no arguments are provided explicitly the compiler can't find a function it can call. – François Andrieux Mar 12 '23 at 11:48
  • This is not minimal, see [mre]. – Passer By Mar 12 '23 at 11:48
  • Clang has a better diagnostic: _"candidate template ignored: substitution failure: deduced non-type template argument does not have the same type as the corresponding template parameter ('unsigned long' vs 'int')"_ – Passer By Mar 12 '23 at 11:49
  • @PasserBy this is fair critique. Please allow a few minutes, I'll try to minimize the example. – eccentricOrange Mar 12 '23 at 11:52
  • @PasserBy I've added a separate smaller example. However, I didn't remove the original code because it prompted an important point in the accepted answer (regarding `constexpr`) which was vital to solving the problem. Please advise if this seems okay – eccentricOrange Mar 12 '23 at 12:55
  • The new example is much better. Leaving the original there is correct, otherwise the answer would not make sense. – Passer By Mar 12 '23 at 14:18

1 Answers1

4

There are multiple fundamental issues with the shown code.

The first one is thatstd::array's 2nd template parameter is not an int. It's a std::size_t.

Put this fact in a blender, together with a using alias, swirl it around, and you end up in, basically, asking a suffering C++ compiler to deduce an int from a size_t. It's an apple vs oranges problem.

To fix this first problem, all template template<int N> declarations must be changed to template<size_t n>.

But when this happens your compiler will be grumpy, because you will be instructing your suffering C++ compiler to create a matrix with more values than there are atoms in our shared universe:

    else if (matrixSize > 0) {
        Matrix<N - 1> temp;

When Matrix is 0 this will end up being a pretty large Matrix, won't you say? Even though this branch will not be taken, this still has to compile, and as a template your suffering C++ compiler needs to, basically, evaluate it. The chances of that are not very likely.

So, you need to fix that, and change all of these to if constexpr.

Then, the code will compile.

#include <array>
#include <iostream>

template <std::size_t N>
using Matrix = std::array<std::array<double, N>, N>;

template <std::size_t N>
Matrix<N - 1> subMatrix(Matrix<N> matrix, size_t focusRow, size_t focusColumn) {
    Matrix<N - 1> returnMatrix;

    int subMatrixRow = 0, subMatrixColumn = 0;
    static const int matrixSize = matrix.size();

    for (size_t matrixRow = 0; matrixRow < matrixSize; matrixRow++) {
        for (size_t matrixColumn = 0; matrixColumn < matrixSize; matrixColumn++) {
            if (matrixRow != focusRow && matrixColumn != focusColumn) {
                returnMatrix[subMatrixRow][subMatrixColumn++] = matrix[matrixRow][matrixColumn];

                if (subMatrixColumn == matrixSize - 1) {
                    subMatrixColumn = 0;
                    subMatrixRow++;
                }
            }
        }
    }

    return returnMatrix;
}

template <std::size_t N>
double getDeterminant(Matrix<N> matrix) {
    static const int matrixSize = matrix.size();
    double determinant;

    if constexpr (matrixSize == 1) {
        determinant = matrix[0][0];
    }

    else if constexpr (matrixSize == 2) {
        determinant = (matrix[0][0] * matrix[1][1]) - (matrix[0][1] * matrix[1][0]);
    }

    else if constexpr (matrixSize > 0) {
        Matrix<N - 1> temp;
        int sign = 1;

        for (size_t dimension = 0; dimension < N; dimension++) {
            temp = subMatrix(matrix, 0, dimension);
            determinant += sign * matrix[0][dimension] * getDeterminant(temp);
            sign = -sign;
        }
    }

    else {
        throw std::invalid_argument("expected square matrix");
    }

    return determinant;
}

int main(int argc, char const* argv[]) {
    static const int length = 3;

    Matrix<length> matrix = {{{1, 2, 3},
                         {4, 5, 6},
                         {7, 8, 9}}};

    printf("determinant = %.2f\n", getDeterminant(matrix));
}
Sam Varshavchik
  • 114,536
  • 5
  • 94
  • 148
  • Thanks for the detailed explanation. Indeed, I had tried `std::size_t` first, but it gave [these exceptions](https://pastebin.com/vaV4FYR4) and https://stackoverflow.com/a/63372812/10907280 directed me to use `int`. Your explanation with pointing to `N - 1` helps a lot here. – eccentricOrange Mar 12 '23 at 12:11
  • 1
    Yes. I didn't immediately see the 2nd issue. After fixing the first compilation error this exact error came up. In C++ there can be different reasons for the same error. It takes time to learn how to understand compiler errors. Don't take the results of a search blindly without reading it and fully understanding how it applies to one's own code. – Sam Varshavchik Mar 12 '23 at 14:35