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 `m_IJ` of a matrix has a minor. A minor is the determinant of a matrix without the Ith row and the Jth 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?