How to apply the colamd algorithm on a sparse matrix in python?

Question

I am trying to apply the Colamd (Column approximate minimum degree permutation) on a scipy sparse matrix.

A incomplete solution can be found from the scipy doc:

from scipy.sparse import csc_matrix, linalg as sla
A = csc_matrix([[1,2,0,4],[1,0,0,1],[1,0,2,1],[2,2,1,0.]])
lu = sla.splu(A,permc_spec = 'COLAMD')
print(lu.perm_c)

However, if A is a singular matrix like below, the 'Python Scipy.sparse RuntimeError: Factor is exactly singular' is raised.

A = csc_matrix([[0,1,0], [0,2,1], [0,3,0]])

To sum up, I am searching a colamd code/algorithm in python working like the colamd(S) matlab method works (ie with singular matrix). I cannot call matlab code from python.

Thank you

Answer 1

Not sure if this will give you the same exact solution (since it's not clear to me how stable AMD/COLAMD is); but I was able to get a solution by adding a small epsilon-identity matrix to make A positive definite. As you make epsilon small (e.g. machine epsilon), I expect to get some reasonable solution.

epsilon=1e-6
from scipy.sparse import csc_matrix, linalg as sla
A = csc_matrix([[0,1,0], [0,2,1], [0,3,0]])+epsilon*csc_matrix(np.eye(3))
lu = sla.splu(A,permc_spec = 'COLAMD')
print(lu.perm_c)

> [0 1 2]

(note for your example problem; pretty much any positive epsilon other than 1.0 gave this permutation)