Determinant of an n*n Matrix

Started by
4 comments, last by Emergent 14 years, 10 months ago
Hi, I've essentially cloned the code on this page to produce the following determinant() method of my Matrix<_HEIGHT,_WIDTH> class:
template<int _HEIGHT, int _WIDTH>
float Matrix<_HEIGHT,_WIDTH>::determinant() const {
	if (_WIDTH==_HEIGHT) {
		if (_WIDTH==1) {
			return this->getElement(0,0);
		} else if (_WIDTH==2) {
			return this->getElement(0,0)*this->getElement(1,1) - this->getElement(1,0)*this->getElement(0,1);
		} else {
			float det=0;
			for (int j1=0; j1<_WIDTH; j1++) {
				Matrix<_HEIGHT-1,_WIDTH-1> matrix;
				for (int i=1; i<_WIDTH; i++) {
					int j2 = 0;
					for (int j=0; j<_WIDTH; j++) {
						if (j!=j1) {
							matrix.setElement(i-1,j2,this->getElement(i,j));
							j2++;
						}
					}
				}
				det += pow(-1.0f,j1) * this->getElement(0,j1) * matrix.determinant();
			}
			return det;
		}
	}
}



As you can see, this is a recursive function. On line 11, it creates an instance of a matrix of dimensions 1 less than the current matrix. Then on line 21, the determinant() method of this new matrix is called. A 3x3 matrix will require the determinants of three 2x2 matrices. Likewise, a 4x4 matrix will require the determinants of four 3x3 matrices. However, when I run this method for a 3x3 matrix, the size of the matrices for each iteration goes down to 0, then flickers between 1 and 0 a few times, then starts increasing until it gets to 1993 and then crashes. However, it shouldn't even create a 1x1 matrix, because a 2x2 matrix will end the recursion by line 6. I even put a breakpoint in my main() function, right at the start, and the problem still occurs. I'm sure this must be an issue with using templates (or just a stupid error I've made). Anybody know what the problem is? Thanks -sixfoottallrabbit
Advertisement
This might be a memory issue - check to see that you don't have a leak.

BTW, it is probably easier to do det(A) by LU decomposition. (b/c det(A) = det(LU) = det(L)det(U)). the determinant of triangular matrices is the product of the diagonal elements. Lots of free LU source code available. Plus you only have to do it once and then you can solve many matrix equations, and it is non recursive.
you really should never actually be computing the determinant like this. Wikipedia has an example
here. Via LU decomposition is much easier, and LU is pretty trivial to write.
I'll echo previous posters that LU beats Leibniz (your method) for this. It's not just easier; it's much faster: Leibniz has exponential complexity; Lu by Gaussian elimination is O(n^3).

If all you need is the absolute value of the determinant, QR decomposition is an even faster way; it's O(n^2). (It's actually O(n^3); LU is better; see subsequent posts).

[Edited by - Emergent on June 1, 2009 11:06:06 AM]
Quote:Original post by Emergent
If all you need is the absolute value of the determinant, QR decomposition is an even faster way; it's O(n^2).


??

Hmm... Doing a full QR decomp is O(n^3) (but more arithmetic than LU). Is there some way to extract the det w/o all the computations (this upsets my sensibilities)
Quote:Original post by gpu_fem
Quote:Original post by Emergent
If all you need is the absolute value of the determinant, QR decomposition is an even faster way; it's O(n^2).


??

Hmm... Doing a full QR decomp is O(n^3) (but more arithmetic than LU). Is there some way to extract the det w/o all the computations (this upsets my sensibilities)


Oops... I left out the "inner loop" that works on each element of a column. Good call: O(n^3) it is; LU beats it.

This topic is closed to new replies.

Advertisement