josdejong/mathjs — pull request #12 by bcuccioli (2013-05-19 05:28:32+00:00)
Title: Compute determinant in polynomial-time.

The idea is to apply Gaussian elimination to the matrix. Each stage of Gaussian elimination changes the determinant in one of the following ways:

Swapping two rows negates the determinant
Scaling a row by a factor c also multiplies the determinant by c
Adding a multiple of one row to another has no effect
In a finite sequence of these operations, you end up with a "row-reduced" matrix, which is either the identity matrix or a matrix that is almost the identity matrix except for some number of zero rows. In the former case, the determinant is the cumulative number we've calculated. In the latter case, the matrix is singular so we just return 0.

An alternate way of thinking about this is Gaussian elimination factors your matrix into a product M = E_1 E_2 ... E_n I, where each E_i is an "elementary matrix" corresponding to one of the above 3 operations - swapping two rows, scaling a row by c, or adding a multiple of one row to another - and having determinant -1, c, or 1, respectively. The matrix I is either the identity or the identity missing some rows as described above, so that det(I) is either 1 or 0. Because det is a ring homomorphism (i.e. det(AB)=det(A)det(B)) this factorization gives an easy way of computing det(M) as the product det(E_1)...det(E_n)det(I).

The advantage of this method is twofold. The algorithm runs in O(n^3), as opposed to calculation by minors which scales as O(n!). It also doesn't use recursion, which might cause a stack overflow for matrices bigger than maybe 5x5.