I wanted to do some exercise and came up with the idea of a good challenge (for my level of course). I tried to implement Laplace's algorithm for computing the determinant, recursively.

#include <math.h>
#include <stdlib.h>
#include <vector>
#include <iostream>
using namespace std;

int getMinimoCount = 0; //ignore these, just to keep track of the recursion.
int calcDetCount = 0;

void printMatrix ( vector< vector<double> > M) {
  //just does what it means
  int size = M.size();
  for( int i = 0; i < size; i++ ) {
    cout << "\t";
    for( int j = 0; j < size; j++ ) {
      cout << M[i][j] << "\t"; }
    cout << endl << endl << endl; }
  cout << endl;
}

vector< vector<double> > getMinimo( vector< vector<double> > src, int I, int J, int ordSrc ) {
  // Compute and return the minimum of the element I J
  // If the element is not in the Ith row or Jth column it will get copied to the minimum matrix
  getMinimoCount++;
  
  vector< vector<double> > minimo( ordSrc-1, vector<double> (ordSrc-1,0));

  int rowCont = 0;
  for( int i=0; i < ordSrc; i++)
  {
    int colCont = 0;
    if ( i != I ) { 
      
      for ( int j=0; j < ordSrc; j++)
      { 
        if ( j != J ) { 
          minimo[rowCont][colCont] = src[i][j];
          colCont++; }
      };
      
      rowCont++; }
  };
  return minimo;
}


double calcDet( vector< vector<double> > src, int ord) {
  // Here be recursion. 
  calcDetCount++;

  if ( ord == 2 ) {

    double mainDiag = src[0][0] * src[1][1];
    double negDiag = src[1][0] * src[0][1];

    return mainDiag - negDiag; }
  else {
    double det = 0;

    for( int J = 0; J < ord; J++) 
    {
      vector< vector<double> > min = getMinimo( src, 0, J, ord);
      
      if ( (J % 2) == 0 ) { det += src[0][J] * calcDet( min, ord-1); }
      else { det -= src[0][J] * calcDet( min, ord-1); }
      
    };
    
    return det;
  }
}

int main() {

  // Just some UI to gather the matrix. not really convinced of this.
  int ord;
  cout << "############## MATRIX DET ##############" << endl << endl;
  cout << " Matrix order: "; cin >> ord; cout << endl;
  vector <vector<double> > mainMatrix( ord, vector<double> (ord, 0));

  cout << """ insert values one row at time. Top to bottom:\n\n""";
  for ( int countY = 0; countY < ord; countY++) {
    for ( int countX = 0; countX < ord; countX++) {
      cin >> mainMatrix[countY][countX];};
  };

  system("CLS");
  cout << "############## MATRIX DET ##############" << endl << endl;
  cout << endl << endl << " This is the input matrix:" << endl << endl << endl;
  printMatrix( mainMatrix );

  system("PAUSE");

  system("CLS");
  cout << "############## MATRIX DET ##############" << endl << endl;
  cout << " Working...!" << endl;
  double det = calcDet( mainMatrix, ord );
  system("CLS");

  cout << endl << endl << "############## MATRIX DET ##############" << endl << endl;
  cout << " Det =\t" << det << endl << endl;
  cout << " getMinimo() chiamata: " << getMinimoCount << " volte" << endl;
  cout << " calcDet() chiamata: " << calcDetCount << " volte" << endl << endl;
  
  return 0;
}

The concept is simple: you have a matrix of order n. While doing this by hand you'd prefer chosing a row that's particularly math friendly; since it's a computer doing the dirty work it really doesn't care about what number he's multiplying.

Every element a_IJ of a matrix has a minor. A minor is the determinant of the matrix without the I-th row and the J-th column. With this we can define the det of a matrix like so:

Sum (-1)^i+j * a_ij * M_ij

(where M_ij is the minimum of the element a_ij)

Once a matrix reach the order == 2 it just computes the determinant since is just a simple multiplication between 4 elements.

At first I had problems finding something that could be used as a matrix object and be passed around from one function to another. I tried bi-dimensional arrays, but they aren't dynamic and couldn't understand how I could pass an array object to a function. I was looking for some hidden matrix type or class to use but i had no luck what so ever.

I came up with "a vector of vectors" which worked but I'm not completely sure it's a really good idea. Not even counting that vector <vector<double> > looks awful.

Secondly, running some recursion tracking I found out that the time it takes ramps up so damn quickly:

• Order 3 = 4 calls
• Order 4 = 17 calls
• Order 5 = 86 calls
• ...
• Order 10 = 2.606.501 calls

Although, even in the latter case, it doesn't take too much: 3 to 4 seconds. Is there a way to reduce the steepness of the curve? This way it gets out of hand way too soon.

I don't have a huge programming experience so I know almost nothing on optimization or good practice either. Since I have some fresh code to work with I'd like to know any error I might be doing and way to optimize this algorithm.

How costly is to cast a type? The (type) x type of cast to be clear.

What could be a better way to input a matrix from the user?