Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

As far as i know LU decomposition allow you to calculate matrix determinant following easy and cheap formula: Det[A] = Det[L] Det[U] = Det[U]

Trying this out in Mathematica 7 gives me correct result by absolute value, i.e. it ignores negative determinants and transforms them to positive ones.

Sample code:

matrixA = {{-1, 2, 3, 5}, {-7, -4, 5, 4}, {-89, 7, 8, -6}, {8, 6, -1, 4}};

Det[matrixA] gives out -2067

but

{lu, p, c} = LUDecomposition[matrixA]

u = lu SparseArray[{i_, j_} /; j >= i -> 1, {4, 4}]

Det[u] is 2067

Well, the question is obvious - how to get correct determinant in Mathematica using LU decomposition?

share|improve this question
add comment

2 Answers 2

up vote 8 down vote accepted

Well, this is because you forgot to take into account the permutation matrix that is output by the LU decomposition routine. You have to remember two facts: 1.) Gaussian elimination generally performs row interchanges for numerical stability reasons, and 2.) the determinant of a matrix changes sign if you interchange rows.

In Mathematica at least, the function that will rescue you is Signature[], which gives the signature of the permutation required to turn the permutation matrix (which Mathematica outputs as a scrambled list of the numbers from 1 to n, where n is the size of your matrix) into the identity matrix.

So, to use your example in Mathematica:

m={{-1, 2, 3, 5}, {-7, -4, 5, 4}, {-89, 7, 8, -6}, {8, 6, -1, 4}};

we compare

Det[m]
-2067

with the following:

{lu, p, c} = LUDecomposition[m]
{{{-1, 2, 3, 5}, {7, -18, -16, -31}, {-8, -11/9, 31/9, 55/9}, {89, 19/2, -963/31, 2067/62}}, {1, 2, 4, 3}, 1}

Now, LUDecomposition[] here outputs three things: the merged $\mathbf{L}$ and $\mathbf{U}$ matrices, the permutation, and the condition number. We can get the tentative determinant by multiplying together the diagonal elements of lu, thus:

Tr[lu, Times]
2067

Here is where Signature comes in:

Signature[{1, 2, 4, 3}]
-1

Note that to turn {1, 2, 4, 3} into {1, 2, 3, 4}, one needs one swap; namely, to swap the third and fourth elements. Since 1 is odd, the signature is -1.

Thus,

Signature[p]*Tr[lu, Times]
-2067

gives the correct answer.

share|improve this answer
    
Very clear and through answer - kudos! –  mth Oct 25 '10 at 11:04
add comment

From the documentation: "Therefore multiplying L and U recreates the original matrix permuted with the row pivot vector." So you need to multiply by Signature[p] to get the correct sign.

share|improve this answer
    
You beat me! No worries, +1! –  J. M. Oct 25 '10 at 11:01
1  
@J. M.: No wonder I was quicker, since your answer was a liiiiitle bit more ambitious! –  Hans Lundmark Oct 25 '10 at 11:22
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.