This is a read only copy of the old FEniCS QA forum. Please visit the new QA forum to ask questions

cbc.block matrices to petsc matrices

+1 vote

How do I convert a block.block_mat.block_mat created using cbc.block to a petsc matrix?

I would like to use the SLEPcEigenSolver but I created block matrices using cbc.block. Alternatively, is it possible to create a block matrix when I already have the blocks using just PETSc matrices?

Thanks!

asked Jul 27, 2016 by sarah FEniCS Novice (430 points)

1 Answer

+1 vote
 
Best answer

Hi, here I did something similar to what you are asking. Blocks are extracted from cbc.block_mat and PETSc.MatNest object is created that is then used by the eigenvalue solver from slepc4py. Alternatively, for many eigenvalue solvers only the action of the operator is required. This is illustrated here (A is a matrix, W is a list of PETSc.Vec)

class EnergyNorm(object):
    def __init__(self, A, W):
        self.A = A
        self.W = W

    def mult(self, mat, x, y):
        alphas = [wi.dot(x) for wi in self.W]
        self.A.mult(x, y)
        for alpha, wi in zip(alphas, self.W):
            y.axpy(alpha, wi)
energy_norm = EnergyNorm(A, W)

# need dummy vector for size
vec = assemble(inner(Constant((1, 1, 1)), v)*dx)
vec = as_backend_type(vec).vec()

from petsc4py import PETSc

E = PETSc.Mat().createPython([vec.getSizes(), vec.getSize()],
                             comm=mpi_comm_world())
E.setPythonContext(energy_norm)
E.setUp()

from slepc4py import SLEPc
# Setup the eigensolver
solver = SLEPc.EPS().create()
solver.setOperators(E, B.mat())
solver.setType(solver.Type.KRYLOVSCHUR)
solver.setDimensions(3, PETSc.DECIDE)
answered Jul 27, 2016 by MiroK FEniCS Expert (80,920 points)
selected Aug 16, 2016 by johannr

Thank you for your help!

...