Search code examples
juliasparse-matrixmatrix-factorizationsuitesparse

Julia: converting CHOLMOD factor to sparse matrix and back again


I have a CHOLMOD factorization of a sparse matrix H, and I want to edit the sparse representation of the upper, lower, and block diagonal factors. How can I do this? When I run the below, the last line doesn't work.

H = sprand(10,10,0.5)
fac = ldltfact(H; shift=0.0)
fD = fac[:D]
D = Base.SparseArrays.CHOLMOD.Sparse(fD)

And is there any way to go in the reverse direction from a sparse matrix to a CHOLMOD.factor?


Solution

  • Extracting the relevant factorization matrices of ldltfact can be a little tedious. The following example shows an example similar to the one in the question with a final test that the extracted matrices recover the original factorized one:

    srand(1)
    pre = sprand(10,10,0.5)
    H = pre + pre' + speye(10,10)
    
    fac = ldltfact(H; shift=0.0)
    P = sparse(1:size(H,1),fac[:p],ones(size(H,1)))
    LD = sparse(fac[:LD]) # this matrix contains both D and L embedded in it
    
    L = copy(LD)
    for i=1:size(L,1)
      L[i,i] = 1.0
    end
    
    D = sparse(1:size(L,1),1:size(L,1),diag(LD))
    
    PHP = P*H*P'
    LDL = L*D*L'
    
    using Base.Test
    @test PHP ≈ LDL
    

    The expected output (and actual on Julia v0.6.3):

    julia> @test PHP ≈ LDL
    Test Passed
    

    Hope this helps.