import numpy as np
from scipy.sparse.linalg import gmres, lgmres

data = np.load("data.npz")
M = data['M']
rhs = data['rhs']

global it
it = 0
norm_rhs = np.linalg.norm(rhs)
def callb(rk):
    global it
    res = np.linalg.norm(M.dot(rk) - rhs)
    print(str(it) + ", " + str(res))
    it += 1

x0, optinfo = lgmres(M, rhs, tol=1e-10, maxiter=55, callback=callb)
print("res:", np.linalg.norm(M.dot(x0) - rhs))