Faster approach of LU decomposition for a symmetric sparse matrix than lu in Matlab?

Dear All,

I have a symmetric sparse matrix A. I want to obtain its LU decomposition. But I found lu(A) took me about 1 second. I am wondering if there is a faster approach than lu.

Thanks a lot in advance.

Bei

[Bruno Luong]

Make sure using LU with 4 or 5 output arguments, because you get filled LU without row/column permutations and this is strongly not recommended for sparse matrix

Answer 1

[Matt J]

Is the matrix positive definite? If so, maybe CHOL?

[Bruno Luong]

Note: LDL is suitable for indefinite Hermitian.

[Benson Gou]

Hi, Matt,

It is close to positive definite. A=H'*H, but H might not be full rank. For the case when H is not full rank, I think we can not use CHOL.

Thanks a lot. Bei

[Benson Gou]

Hi, Bruno,

Thanks for your suggestion. I checked LDL.m but it does not take sparse matrix. I think it may be time-consuming.

Thanks a lot again.

Bei

[Bruno Luong]

It takes sparse, but you have to call with 4 arguments, sparse matrix needs permutation for remains sparse, the exact same comment I put or LU. If you don't want to bother with permutation, work with FULL matrix, no point to make a fixation on SPARSE.

>> A=sprand(10,10,0.3)

A =

     (7,1)       0.2575
    (10,1)       0.4733
     (7,2)       0.8407
     (9,2)       0.1966
    (10,2)       0.3517
     (2,3)       0.5060
     (7,3)       0.2543
     (4,4)       0.9593
     (6,4)       0.1493
    (10,4)       0.8308
     (2,5)       0.6991
     (5,5)       0.5472
     (8,5)       0.2435
     (7,6)       0.8143
     (1,7)       0.7513
     (3,7)       0.8909
     (9,7)       0.2511
    (10,7)       0.5853
     (5,8)       0.1386
     (8,8)       0.9293
     (9,8)       0.6160
    (10,8)       0.5497
     (1,9)       0.2551
     (8,10)      0.3500
    (10,10)      0.9172

>> [L,D,P,S] = ldl(A)

L =

     (1,1)       1.0000
     (2,1)       0.3415
     (2,2)       1.0000
     (7,2)       0.2903
     (9,2)       1.1320
     (3,3)       1.0000
     (6,3)       0.2265
     (7,3)       0.0493
     (8,3)       1.0000
     (9,3)       0.0735
    (10,3)       0.5116
     (4,4)       1.0000
     (5,5)       1.0000
     (7,5)       0.3815
     (8,5)       1.0000
     (6,6)       1.0000
     (8,6)      -0.3815
     (9,6)      -0.2903
    (10,6)       0.5141
     (7,7)       1.0000
     (8,8)       1.0000
     (9,9)       1.0000
    (10,9)      -0.8834
    (10,10)      1.0000

D =

     (1,1)       1.0000
     (2,2)       0.8834
     (4,3)       1.0000
     (3,4)       1.0000
     (5,5)       1.0000
     (7,6)       1.0000
     (6,7)       1.0000
     (7,7)      -0.0345
     (8,8)      -1.0000
     (9,9)      -1.1320
    (10,10)      0.8834

P =

     (5,1)        1
     (8,2)        1
     (3,3)        1
     (7,4)        1
     (4,5)        1
     (1,6)        1
    (10,7)        1
     (6,8)        1
     (9,9)        1
     (2,10)       1

S =

     (1,1)       4.6979
     (2,2)       3.2506
     (3,3)      21.0084
     (4,4)       1.0210
     (5,5)       1.3518
     (6,6)       6.5604
     (7,7)       0.1872
     (8,8)       1.0374
     (9,9)       1.5648
    (10,10)      0.4498

>>

[Benson Gou]

Hi, Bruno,

Thanks for your great help.

Originally I am looking for the lower triangular matrix L and the diagonal matrix D in L'*D*L=A when A is not full rank. My purpose is to look for additional measurements which can make A full rank. To do this, I need to form a matrix W which is formed by rows from the inverse of L. Those rows in the inverse of L correspond to zero diagonals in D. My own ldl.m code is very slow which took more than 100 seconds. You suggested to use Matlab's ldl.m which is much faster than mine. But Matlab's ldl.m changes the order of rows in A. I want to find out the ldl.m code in Matlab so that I can modify it for my use.

Thanks a lot for your great help. You have a good weekend!

Bei