LU factorisation of dense and sparse matrix

General Usage
,
1

Hi guys,

I’ve been trying convert MATLAB code to Julia, and I observed LU factorization for sparse and dense matrix of the same matrix give different number of non-zero fill-ins.

Consider the following example:

using LinearAlgebra
using SparseArrays
using Test

#LU factorisation of sparse matrix
A_sparse = sprand(Float64, 10, 10, 0.3)
F = lu(A_sparse)
sL, sU = F.L, F.U

#LU factorisation of dense matrix
L, U = lu(Matrix(A_sparse))

#Testing
nnz(sL) .== nnz(sparse(L))
nnz(sU) .== nnz(sparse(U))

@test Matrix(sL) ≈ L atol=1e-5
@test Matrix(sU) ≈ U atol=1e-5

All tests fail for me. Any ideas why this might be?

2

My guess would be that sparse/dense routines use different permutations, where the sparse version tries to minimise fill-in while the dense prioritises numeric stability.

3

Definitely. The main trick of sparse LU is to find a fill-reducing permutation (such as minimum degree).

However, it’s possible to specify a custom permutation to the sparse-LU routine (`lu` support for custom permutation (like in `cholesky`) · Issue #116 · JuliaSparse/SparseArrays.jl · GitHub), so in principle you should be able to force it to use the same permutation as the dense LU routine to get them to produce the same fill-in. Not sure why this would be worth the trouble, however.