Scipy LinearOperator with FreeDofs

Hi there,

I’m working on optimizing an inverse of a large block sparse matrix : [ [A C], [B, D] ] I’ve found that I can use the smaller blocks, and now I just need to be able to do an inverse of the matrix W = A - C D-1 B. The matrix D is a mass matrix so I am able to invert it and use the freedofs from the space to give to the inverse. The other matrices are also derived from Bilinear forms, in particular A is derived from the curl curl form plus a weighted mass form.

To provide a method for inverting this matrix, I was thinking of using the Matrix Free method outlined in the section on interfacing with Numpy and Scipy in the documentation. It requires only a method for applying the matrix W, which is no problem.

The problem though is that I don’t know how to get this process to recognize the freedofs. In the documentation the example that is given is this:

[code]import scipy
import scipy.sparse.linalg

tmp1 = f.vec.CreateVector()
tmp2 = f.vec.CreateVector()
def matvec(v):
tmp1.FV().NumPy()[:] = v
tmp2.data = a.mat * tmp1
return tmp2.FV().NumPy()

A = scipy.sparse.linalg.LinearOperator( (a.mat.height,a.mat.width), matvec)

u.vec.FV().NumPy()[:], succ = scipy.sparse.linalg.cg(A, f.vec.FV().NumPy())[/code]

It mentions in particular that the matrix a.mat is derived without Dirichlet dofs.

My question then is how can I incorporated Dirichlet dofs into the above example, and hence into my own process? Is there a way to do this that someone knows?

Thanks,
Pieter Vandenberge
Portland State

Typically the freedofs enter into the preconditioner,
The range of the preconditioner defines the sub-space satisfying the (homogeneous) Dirichlet constraints.

The simplest possibility is to use a projector:

pre = Projector(mask=fes.FreeDofs(), range=True)

best,
Joachim