In [21]:

    

from dolfin import *
import numpy as np
import matplotlib.pylab as plt
from scipy.sparse import csc_matrix, spdiags
mesh = UnitSquareMesh(2,2)
parameters['reorder_dofs_serial'] = False
V = VectorFunctionSpace(mesh,"CG",2)
Q = FunctionSpace(mesh,"CG",1)
W = V*Q
(u, p) = TrialFunctions(W)
(v, q) = TestFunctions(W)
print V.dim(), Q.dim()


    

    




DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.
DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.
DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.






    




50 9

In [22]:

    

def boundary(x, on_boundary):
    return on_boundary
bc = DirichletBC(W.sub(0),Expression(("0","0")), boundary)
parameters['linear_algebra_backend'] = 'uBLAS'
BQB = assemble(inner(u,v)*dx+ div(v)*p*dx+div(u)*q*dx)
bc.apply(BQB)
BQB=BQB.sparray()


    

    




DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.

In [23]:

    

Q = BQB[0:np.max(range(V.dim()))+1,0:np.max(range(V.dim()))+1]


    

In [24]:

    

Qdiag = Q.diagonal()


    

In [25]:

    

B = BQB[np.max(range(V.dim()))+1:np.max(range(W.dim()))+1,0:np.max(range(V.dim()))+1]
Bt = BQB[0:np.max(range(V.dim()))+1,np.max(range(V.dim()))+1:np.max(range(W.dim()))+1]


    

In [27]:

    

d = spdiags(1.0/Qdiag, 0, len(Qdiag), len(Qdiag))


    

In [38]:

    

L = B*d*Bt
L.shape


    

    
Out[38]:





(9, 9)

In [39]:

    

U, s, Vh =svd(L.todense())


    

In [42]:

    

import petsc4py
import sys
petsc4py.init(sys.argv)
from petsc4py import PETSc


    

In [44]:

    

PETSc.PC.setFactorShift


    

    




---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-44-69faea4440fe> in <module>()
----> 1 PETSc.PC.setFactorShift(1)

TypeError: descriptor 'setFactorShift' requires a 'petsc4py.PETSc.PC' object but received a 'int'