Preconditioned Conjugate Gradient and LinearOperator in python

Question

[Homework] I am going to solve the linear system Ax=b by the Preconditioned Conjugate Gradient method, and I use spilu function from scipy.sparse.linalg for the preconditioner. A is a sparse symmetric 162*162 matrix. Since the spilu gives an approximation to the inverse of A, say M approximates A, and so spilu(A) gives M^-1, which is the preconditioner. I find that we can directly gives the preconditioner in the python Conjugate Gradient function, but my code below does not work.

M_inverse=scipy.sparse.linalg.spilu(A)
M2=scipy.sparse.linalg.LinearOperator((162,162),M_inverse.solve)
x3=scipy.sparse.linalg.cg(A,b,M2)
TypeError                                 Traceback (most recent call last)
<ipython-input-84-86f8f91df8d2> in <module>()
----> 1 x3=scipy.sparse.linalg.cg(A,b,M2)

/Users/ruobinghan/anaconda/lib/python3.4/site-packages/scipy/sparse/linalg/isolve/iterative.py in cg(A, b, x0, tol, maxiter, xtype, M, callback)

/Users/ruobinghan/anaconda/lib/python3.4/site-packages/scipy/sparse/linalg/isolve/iterative.py in non_reentrant(func, *a, **kw)
     83     try:
     84         d['__entered'] = True
---> 85         return func(*a, **kw)
     86     finally:
     87         d['__entered'] = False

/Users/ruobinghan/anaconda/lib/python3.4/site-packages/scipy/sparse/linalg/isolve/iterative.py in cg(A, b, x0, tol, maxiter, xtype, M, callback)
    219 @non_reentrant
    220 def cg(A, b, x0=None, tol=1e-5, maxiter=None, xtype=None, M=None, callback=None):
--> 221     A,M,x,b,postprocess = make_system(A,M,x0,b,xtype)
    222 
    223     n = len(b)

/Users/ruobinghan/anaconda/lib/python3.4/site-packages/scipy/sparse/linalg/isolve/utils.py in make_system(A, M, x0, b, xtype)
    108         x = zeros(N, dtype=xtype)
    109     else:
--> 110         x = array(x0, dtype=xtype)
    111         if not (x.shape == (N,1) or x.shape == (N,)):
    112             raise ValueError('A and x have incompatible dimensions')

TypeError: float() argument must be a string or a number, not 'LinearOperator'

Also, the question hints I will need to use LinearOperator interface, I do not understand what is exactly LinearOperator doing and why we need it here.

Any suggestion would be appreciated! Thanks in advance!

It's strongly recommended to paste code and errors in formatted text rather than taking screenshots. — machine yearning, Sep 30 '15 at 12:58
@machineyearning,thanks, I use the LinearOperator, and now a new error occurred, isn't the page tells me to uss cg(A,b,M2)? Why now LinearOperator becomes a problem?http://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.cg.html#scipy.sparse.linalg.cg — Bob, Sep 30 '15 at 13:53
Test your code in the interactive interpreter and figure out exactly where it's breaking so you can provide a [Minimal, Complete, and Verifiable example](https://stackoverflow.com/help/mcve). Try the functions that appear to be breaking, with very simple values, so you can be sure the functions are working as you expect them to. If the simple examples you come up with still don't seem to work, then post a new question with your simple example code (you should be able to narrow it down to 10 lines or less) and the exact error message you're getting, which will have a very short stack trace. — machine yearning, Sep 30 '15 at 14:12

score 2 · Answer 1 · answered Mar 24 '16 at 15:24

I think the parameters are in wrong order,

x3=scipy.sparse.linalg.cg(A,b,M2)

In the Error message:

220 def cg(A, b, x0=None, tol=1e-5, maxiter=None, xtype=None, M=None, 
callback=None):
--> 221     A,M,x,b,postprocess = make_system(A,M,x0,b,xtype)

M2 is in the place of x0 - the initial guess of the solution but not the preconditioner. In my host, with correct order, class-LinearOperator is functioning well.

correct version

x3=scipy.sparse.linalg.cg(A,b,M=M2)

Please use "key word" arguments as often as possible.

Preconditioned Conjugate Gradient and LinearOperator in python

1 Answers1