How can I calculate SVD in fenics


58
views
0
8 weeks ago by
Ali  
Hi All,


I am Trying to figure out how can I calculate the smallest singular value in fenics.
I have mentioned an example in below that I was trying to print out the smallest singular value of stiffness matrix A.
""" Solve the Inf-Sup in fenics.

TODO: consider to convert that to 3D form.
"""
import dolfin
from dolfin import *
import scipy
from scipy.linalg import sqrtm
from petsc4py import PETSc
from numpy import array2string
import numpy as np
import sys
import mshr


#Create mesh
mesh = UnitSquareMesh(12, 12)

#Define function spaces

BDM = VectorElement("BDM", mesh.ufl_cell(), 2, dim=2)
CG = VectorElement("CG", mesh.ufl_cell(), 1, dim=2)

element = MixedElement([CG, BDM])
V = FunctionSpace(mesh, element)
w = Function(V)
#Define DOFs
Gdof = V.dim()

# Define trial and test functions
Z ,Q = TrialFunctions(V)
yup ,R = TestFunctions(V)

# Define bilinear form and right hand side

a = (dot(Q[0], grad(yup[0])) + dot(Q[1], grad(yup[1])) \
+ dot(Z[0], yup[0]) + dot(Z[1], yup[1]) + dot(grad(Z[0]), grad(yup[0])) + dot(grad(Z[1]), grad(yup[1])) \
+ dot(Q[0], R[0]) + dot(Q[1], R[1]) + ((div(Q[0]) * div(R[0]))) + ((div(Q[1]) * div(R[1]))))*dx
#####################################################


print('Global_DOF',Gdof)

matA = assemble(a)
#print out stiffness matrix
stringA = str(matA.array())
outf = open("stiffmat.txt","w")
outf.write(array2string(matA.array(),max_line_width=sys.maxint))
outf.close()

exit()

#calculate SVD of stiffness matrix
S = SLEPc.SVD()
S.create()
S.setOperator(matA)
S.setType(S.Type.LANCZOS)
S.solve()

nconv = S.getConverged()
if nconv > 0:
    v, u = matA.getVecs()
    for i in range(nconv):
        sigma = S.getSingularTriplet(i, u, v)
        print sigma

But It is not working.



Thanks,
~Ali
Community: FEniCS Project
Please login to add an answer/comment or follow this question.

Similar posts:
Search »