Cholesky factor of an assembled matrix

Linear Algebra
1

Apparently I messed up before and posted an unfinished question. I tried deleting it, don’t know how well that went. Anyway, here is the question:

I have an assembled matrix M from a form a, say

u = TrialFunction(V)
v = TestFunction(V)
a = f * u*v*dx

where f(\equiv f) is some function. If I know f>\alpha > 0 , then M is positive definite and there is a unique Cholesky factorization M = LL^T .

Does fenics have a way for me to get the matrix L ? (or L^T for that matter?)

There seem to be preconditioners based on Cholesky factorization, but I can’t seem to find a method which allow me to extract the factor?

2

So I have found a way to do this. Not the most elegant, but it appears to work. It requires the scikit.sparse. If M is the matrix to be factorized we proceed as follows.

Note that the will produce matrices L and P such that PMP^T = LL^T , where L is lower triangular and P is a permutation matrix to maximize sparsity, since that is what scikit.sparse supports.

from sksparse import cholesky
from dolfin import as_backend_type, Matrix, PETScMatrix
from scipy.sparse import csr_matrix
from petsc4py import PETSc

Mmat = as_backend_type(M).mat()
S = csr_matrix(Mmat.getValuesCSR()[::-1], shape=Mmat.size)
S = S.tocsc()
factor = cholesky(S)
L_ = factor.L().tocsr()
P_ = factor.P()
IL = L_.indptr
JL = L_.indices
DL = L_.data
shape = L_shape

L = PETSc.Mat().createAIJWithArrays(shape,(IL,JL,DL))

li = list(range(len(P_)+1))
P = PETSc.Mat().createAIJWithArrays(shape,(li,P_,[1 for _ in li[1:]]))