0
$\begingroup$

Use Mathematica to decompose a given matrix $Q$ as $LDU$, where $L$, $D$ is a triangular matrix, diagonal matrix, respectively, and $U$ is the transpose of $L$.

The built-in LUDecomposition function fails for non-singular matrices.

I searched on the internet and found the NCAlgebra suite contains such a function.

<<NCAlgebra`
{ldl,p,s,rank}=NCLDLDecomposition@mat;
{ll,dd,uu}=GetFullLDUMatrices[ldl,s]

But I checked that ll.dd.uu is not the same as mat. I am not sure whether this is a bug.

Are there any more stable alternative ways to compute matrix LDU decomposition?

Improve this question
$\endgroup$
8
  • 1
    $\begingroup$ Do you need factorization for matrices with exact numbers? Or can you settle with machine precision numbers? Do you really need the factors? Or do you just need a good solver for symmetric, indefinite matrices? Do you need it only for General dense matrices? Or do you need a efficient factorization for sparse matrices? $\endgroup$ Commented Aug 13 at 13:26
  • 1
    $\begingroup$ I presume you mean "exact" instead of "accurate"? I think you should expand a bit on your use case and give a usage example. It is still not clear to me if you want that factorization to work for symbolic matrices -- or which type of number you mean to employ. $\endgroup$ Commented Aug 13 at 15:00
  • 1
    $\begingroup$ Related: mathematica.stackexchange.com/q/83469/38178 $\endgroup$ Commented Aug 13 at 15:01
  • 1
    $\begingroup$ There is ResourceFunction["LDLDecomposition"] in the Wolfram Function Repository, $\endgroup$ Commented Aug 13 at 15:51
  • 1
    $\begingroup$ I would (maybe naively) expect that to indicate it cannot be factored in that way. If memory serves there is no opportunity for pivoting. Maybe I'm wrong about that though. $\endgroup$ Commented Aug 13 at 20:07

1 Answer 1

Reset to default
1
$\begingroup$

As documented here, the resulting decomposition is such that $L D L^T = P M P^T$ hence it is

ll.dd.uu - mat[[p,p]]

that will be zero. Without permutations such decompositions may not exist. Think

$$ m = \begin{bmatrix} 0 & 1\\ 1 & 1 \end{bmatrix} $$

Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.