In [2]:

    

%load_ext autoreload
%autoreload 2


    

In [3]:

    

import numpy as np
import matplotlib.pyplot as plt
from scipy import io
import pyamg
import gen_matrix_2d as gm2
import cmocean
import scipy.sparse.linalg as sla
import scipy.sparse as sp
%matplotlib inline


    

In [4]:

    

#sla.spsolve(A,b)


    

In [7]:

    

dct = io.loadmat('../mat/conv_exp/conv_40x40.mat')
A = dct['A'].tocsr()
b = dct['b'].toarray()


    

In [8]:

    

n = A.shape[0]
M2 = sla.spilu(A)
M_x = lambda x: M2.solve(x)
M = sla.LinearOperator((n,n),M_x)


    

    




/home/oevans/local/anaconda3/lib/python3.6/site-packages/scipy-0.19.0-py3.6-linux-x86_64.egg/scipy/sparse/linalg/dsolve/linsolve.py:315: SparseEfficiencyWarning: splu requires CSC matrix format
  warn('splu requires CSC matrix format', SparseEfficiencyWarning)






    




-----------------------------------------
RuntimeErrorTraceback (most recent call last)
<ipython-input-8-6f6ed89cbd68> in <module>()
      1 n = A.shape[0]
----> 2 M2 = sla.spilu(A)
      3 M_x = lambda x: M2.solve(x)
      4 M = sla.LinearOperator((n,n),M_x)

/home/oevans/local/anaconda3/lib/python3.6/site-packages/scipy-0.19.0-py3.6-linux-x86_64.egg/scipy/sparse/linalg/dsolve/linsolve.py in spilu(A, drop_tol, fill_factor, drop_rule, permc_spec, diag_pivot_thresh, relax, panel_size, options)
    329         _options.update(options)
    330     return _superlu.gstrf(N, A.nnz, A.data, A.indices, A.indptr,
--> 331                           ilu=True, options=_options)
    332 
    333 

RuntimeError: Factor is exactly singular

In [ ]:

    

sla.lgmres(A,b,M)


    

In [ ]:

    

%%time
mg = pyamg.rootnode_solver(A)


    

In [ ]:

    

mg.psolve(np.ones_like(b))


    

In [ ]:

    




    

In [ ]:

    




    

In [ ]:

    

r = mg.solve(b)


    

In [ ]:

    

sla.lgmres(A,b,M=mg.aspreconditioner())


    

In [ ]:

    

mg.aspreconditioner()


    

In [ ]:

    

%%time
rad1 = mg.solve(b,accel='lgmres',tol=1e-4)


    

In [ ]:

    

%%time
rad2 = sla.lgmres(A,b)


    

In [ ]:

    

rad[np.isnan(rad)] = 0


    

In [ ]:

    

mg


    

In [ ]:

    

plt.plot(rad)


    

In [ ]:

    

%%time 
sla.lgmres(A,b)


    

In [ ]:

    

_


    

In [ ]:

    

_2


    

In [ ]:

    

rad1,info = Out[19]


    

In [ ]:

    

# Compare irradiance plots for two different methods
ns = 80
nth = 80
figsize=[10,8]
plt.figure(figsize=figsize)
gm2.plot_irrad(
    rad2[0],
    var_lengths=[ns,ns,nth],
    var_order=[2,1,0],
    title='80x80',
    cmap=cmocean.cm.thermal
)


    

In [ ]: