3
$\begingroup$

If A is a symmetric matrix, is there a fast algorithm for LU factorization? I know this algorithm for non-symmetric matrix.

    For k = 1,..,n        For i = k + 1,...,n          mult := a_{ik}/a_{kk}          a_{jk} := mult          For j = k+1,...,n              a_{ij} := a_{ij} - mult * a_{kj}          Endfor      Endfor    Endfor 
  • 0
    Use Lapack. $$ $ $2012-11-26

2 Answers 2

1

See Wikipedia article, try using Cholesky decomposition.

Added Later: The matrix has to be positive definite. So this is just a partial answer to the question.

  • 0
    @J.M. Yes, that's correct. Sorry. I need to add a condition that the matrix is positive definite.2012-08-24
1

The cholesky-decomposition can be made "robust" - just keep track of zeros and negative signs in the diagonal. Here is a code-snippet in my (proprietary,sorry) MatMate-language, which should be easily translatable into a C- or Basic - or something-else routine:

TMP = X           // X is the symmetric matrix L   = Null(X)     // Null-matrix of same size as X;                   //         shall become the cholesky-factor sg  = L[*,1]      // vector which keeps track of the signs in the diagonal  c  = 0           // index for the current column/row     // repeat the following up to the number of rows/colums of X c=c+1             d = w(TMP[c,c])   // the value of the diagonal element TMP[c,c]             // if zero or lower/equal required precision finish procedure s = sign(d) sg[c] = s         // store the sign of the diagonal element L[*,c] = TMP[*,c]/sqrt(s*d)  // fill the c'th column TMP = TMP - s * L[*,c] * L[*,c]'  // reduce tmp;                                   //   the apostroph means transposition 

Check result:

chk = L * mkdiag(sg) * L'  // check: re-compose cholesky-factors                                    // with signed-diagonal err = X - chk   // should be zero 

We can even include, that the algorithm terminates when d becomes zero: the cholesky-factor L reflects then exactly the rank of X