How to use the petsc4py.PETSc.PC function in petsc4py

def create_solver():
    "Create the KSP solver"
    # create linear solver
    ksp = PETSc.KSP()
    ksp.create(PETSc.COMM_WORLD)

    # Use algebraic multigrid:
    pc = ksp.getPC()
    pc.setType(PETSc.PC.Type.GAMG)
    ksp.setFromOptions()

    ksp.setInitialGuessNonzero(True)

    return ksp

# Attach rigid deformations to A
    # Functions
    Z = [interpolate(Constant((1, 0)), V),
         interpolate(Constant((0, 1)), V),
         interpolate(Expression(('x[1]', '-x[0]'), degree=1), V)]
    # The basis
    Z = VectorSpaceBasis([z.vector() for z in Z])
    Z.orthonormalize()
    A.set_nullspace(Z)
    A.set_near_nullspace(Z)

    A = as_petsc(A)
    # Setup the preconditioner in petsc
    pc = PETSc.PC().create()
    pc.setType(PETSc.PC.Type.HYPRE)
    pc.setOperators(A)
    # Other options
    opts = PETSc.Options()
    opts.setValue('pc_hypre_boomeramg_cycle_type', 'V')
    opts.setValue('pc_hypre_boomeramg_relax_type_all',  'symmetric-SOR/Jacobi')
    opts.setValue('pc_hypre_boomeramg_coarsen_type', 'Falgout')  
    pc.setFromOptions()         

    # Wrap for cbc.block
    B00 = BlockPC(pc)
    # The Q norm via spectral
    Qi = Q.sub(0).collapse()
    B11 = inverse(VectorizedOperator(HsNorm(Qi, s=-0.5), Q))

    return block_diag_mat([B00, B11])

A = as_backend_type(assemble(b00))

    # Attach rigid deformations to A
    # Functions
    Z = [interpolate(Constant((1, 0)), V),
         interpolate(Constant((0, 1)), V),
         interpolate(Expression(('x[1]', '-x[0]'), degree=1), V)]
    # The basis
    Z = VectorSpaceBasis([z.vector() for z in Z])
    Z.orthonormalize()
    A.set_nullspace(Z)
    A.set_near_nullspace(Z)

    A = as_petsc(A)
    # Setup the preconditioner in petsc
    pc = PETSc.PC().create()
    pc.setType(PETSc.PC.Type.HYPRE)
    pc.setOperators(A)
    # Other options
    opts = PETSc.Options()
    opts.setValue('pc_hypre_boomeramg_cycle_type', 'V')
    opts.setValue('pc_hypre_boomeramg_relax_type_all',  'symmetric-SOR/Jacobi')
    opts.setValue('pc_hypre_boomeramg_coarsen_type', 'Falgout')  
    pc.setFromOptions()         

    # Wrap for cbc.block
    B00 = BlockPC(pc)
    # The Q norm via spectral
    Qi = Q.sub(0).collapse()
    B11 = inverse(VectorizedOperator(HsNorm(Qi, s=-0.5), Q))

    return block_diag_mat([B00, B11])

jac_mat = PETSc.Mat().createPython([(lsize, size), (lsize, size)],
                                               comm=system.comm)
            if trace: debug("petsc matrix creation DONE for %s" % voi)
            jac_mat.setPythonContext(self)
            jac_mat.setUp()

            if trace:  # pragma: no cover
                debug("creating KSP object for system", system.pathname)

            ksp = self.ksp[voi] = PETSc.KSP().create(comm=system.comm)
            if trace: debug("KSP creation DONE")

            ksp.setOperators(jac_mat)
            ksp.setType(self.options['ksp_type'])
            ksp.setGMRESRestart(1000)
            ksp.setPCSide(PETSc.PC.Side.RIGHT)
            ksp.setMonitor(Monitor(self))

            if trace:  # pragma: no cover
                debug("ksp.getPC()")
                debug("rhs_buf, sol_buf size: %d" % lsize)
            pc_mat = ksp.getPC()
            pc_mat.setType('python')
            pc_mat.setPythonContext(self)

        if trace:  # pragma: no cover
            debug("ksp setup done")

        if self.preconditioner:
            self.preconditioner.setup(system)

stop=L_range[0]+neqns,
                                            step=1,
                                            dtype="i")
        else: #assume blocked
            self.pressureDOF = numpy.arange(start=L_range[0],
                                            stop=L_range[0]+neqns,
                                            step=nc,
                                            dtype="i")
            velocityDOF = []
            for start in range(1,nc):
                velocityDOF.append(numpy.arange(start=L_range[0]+start,
                                                stop=L_range[0]+neqns,
                                                step=nc,
                                                dtype="i"))
            self.velocityDOF = numpy.vstack(velocityDOF).transpose().flatten()
        self.pc = p4pyPETSc.PC().create()
        self.pc.setType('fieldsplit')
        self.isp = p4pyPETSc.IS()
        self.isp.createGeneral(self.pressureDOF,comm=p4pyPETSc.COMM_WORLD)
        self.isv = p4pyPETSc.IS()
        self.isv.createGeneral(self.velocityDOF,comm=p4pyPETSc.COMM_WORLD)
        self.pc.setFieldSplitIS(('velocity',self.isv),('pressure',self.isp))

def __init__(self,L):
        L_sizes = L.getSizes()
        L_range = L.getOwnershipRange()
        print "L_sizes",L_sizes
        neqns = L_sizes[0][0]
        print "neqns",neqns
        self.saturationDOF = numpy.arange(L_range[0],L_range[0]+neqns/2,dtype="i")
        #print "saturation",self.saturationDOF
        self.pressureDOF = numpy.arange(L_range[0]+neqns/2,L_range[0]+neqns,dtype="i")
        #print "pressure",self.pressureDOF
        self.pc = p4pyPETSc.PC().create()
        self.pc.setType('fieldsplit')
        self.isp = p4pyPETSc.IS()
        self.isp.createGeneral(self.saturationDOF,comm=p4pyPETSc.COMM_WORLD)
        self.isv = p4pyPETSc.IS()
        self.isv.createGeneral(self.pressureDOF,comm=p4pyPETSc.COMM_WORLD)
        self.pc.setFieldSplitIS(self.isp)
        self.pc.setFieldSplitIS(self.isv)
    def setUp(self, global_ksp=None):

from petsc4py import PETSc
    from hsmg import HsNorm
    
    V, Q, Y = W
    
    # H1
    u, v = TrialFunction(V), TestFunction(V)
    b00 = inner(grad(u), grad(v))*dx + inner(u, v)*dx
    # NOTE: since interpolation is broken with MINI I don't interpolate
    # here the RM basis to attach the vectros to matrix
    A = assemble(b00)

    A = as_petsc(A)
    # Setup the preconditioner in petsc
    pc = PETSc.PC().create()
    pc.setType(PETSc.PC.Type.HYPRE)
    pc.setOperators(A)
    # Other options
    opts = PETSc.Options()
    opts.setValue('pc_hypre_boomeramg_cycle_type', 'V')
    opts.setValue('pc_hypre_boomeramg_relax_type_all',  'symmetric-SOR/Jacobi')
    opts.setValue('pc_hypre_boomeramg_coarsen_type', 'Falgout')  
    pc.setFromOptions()         
    # Wrap for cbc.block
    B00 = BlockPC(pc)

    p, q = TrialFunction(Q), TestFunction(Q)
    B11 = LumpedInvDiag(assemble(inner(p, q)*dx))
    
    # The Y norm via spectral
    Yi = Y.sub(0).collapse()
    B22 = inverse(VectorizedOperator(HsNorm(Yi, s=-0.5), Y))

sigma, tau = TrialFunction(V), TestFunction(V)
        
        a = inner(div(sigma), div(tau))*dx
        if not hdiv0:
            a += inner(sigma, tau)*dx

        f = Constant(np.zeros(V.ufl_element().value_shape()))
        L = inner(tau, f)*dx

        A, _ = assemble_system(a, L, bc)

        # AMS setup
        Q = FunctionSpace(mesh, 'CG', 1)
        G = DiscreteOperators.build_gradient(V, Q)

        pc = PETSc.PC().create(mesh.mpi_comm().tompi4py())
        pc.setType('hypre')
        pc.setHYPREType('ams')

        # Attach gradient
        pc.setHYPREDiscreteGradient(mat(G))

        # Constant nullspace (in case not mass and bcs)
        constants = [vec(interpolate(c, V).vector())
                     for c in (Constant((1, 0)), Constant((0, 1)))]

        pc.setHYPRESetEdgeConstantVectors(*constants)

        # NOTE: term mass term is accounted for automatically by Hypre
        # unless pc.setPoissonBetaMatrix(None)
        if hdiv0: pc.setHYPRESetBetaPoissonMatrix(None)

from numpy import hstack
    from petsc4py import PETSc
    from hsmg import HsNorm
    
    V, Q, Y = W
    
    # H1
    u, v = TrialFunction(V), TestFunction(V)
    b00 = inner(grad(u), grad(v))*dx + inner(u, v)*dx
    # NOTE: since interpolation is broken with MINI I don't interpolate
    # here the RM basis to attach the vectros to matrix
    A = assemble(b00)

    A = as_petsc(A)
    # Setup the preconditioner in petsc
    pc = PETSc.PC().create()
    pc.setType(PETSc.PC.Type.HYPRE)
    pc.setOperators(A)
    # Other options
    opts = PETSc.Options()
    opts.setValue('pc_hypre_boomeramg_cycle_type', 'V')
    opts.setValue('pc_hypre_boomeramg_relax_type_all',  'symmetric-SOR/Jacobi')
    opts.setValue('pc_hypre_boomeramg_coarsen_type', 'Falgout')  
    pc.setFromOptions()         
    # Wrap for cbc.block
    B00 = BlockPC(pc)

    p, q = TrialFunction(Q), TestFunction(Q)
    B11 = LumpedInvDiag(assemble(inner(p, q)*dx))
    
    # The Y norm via spectral
    Yi = Y.sub(0).collapse()