Use Fipy in parallel with GMRES from Petsc

[aph4499]

Hello,
I am trying to run Fipy in parallel using Petsc GMRES to solve the linear system of equations. However I believe something is not correct in my code because the residual from GMRES iteration is empty at each sweep. I have placed a sample code below to illustrate what I am doing. If anyone could have a look and provide some feedback, that would be really helpful.
Thanks

# create grid
m = Grid1D(nx=len(z), Lx=profiles.z[0][-1])

# create variables
c = CellVariable(name="c", mesh=m, hasOld=True, value=profiles.c[t0])
...

# create equations
eq_c = TransientTerm(var=c) == ...
...
full_eq = eq_c & ...

soln = plt.recarray(...)

# initialize
soln.c[0] = c
...

# setup petsc with gmres solver and ilu pc
ksp = PETSc.KSP().create()
ksp.setType('gmres')
pc = ksp.getPC()
pc.setType('ilu')
ksp.setConvergenceHistory()
ksp.view()

# Solve
for t in range(nnt-1):

    c.updateOld()
    ...

    iter = 0
    res = 1.
    while ( res > 1e-7 and iter < 50 ):
        res = full_eq.sweep(dt=dt_soln)
        iter = iter + 1
    
    soln.c[t+1] = c
    ...

    residuals = ksp.getConvergenceHistory()
    print("residuals = ", residuals)

[guyer]

The ksp you define here is not used by full_eq.sweep(). The instance that FiPy creates when you, e.g., call full_eq.sweep(dt=dt_soln, solver=LinearGMRESSolver(precon='ilu')) is not exposed outside of the _solve_() method. I'm currently doing a major revamp on the solvers and I'll give some thought on how to make this information accessible.

[guyer]

In the meantime, you can modify your local copy of PETScKrylovSolver._solve_() to add this diagnostic information to its ksp.

[guyer]

Does this work for you?

$ PETSC_OPTIONS="-ksp_monitor" python examples/diffusion/steadyState/mesh20x20/tri2Dinput.py
  0 KSP Residual norm 4.827005221118e+00 
  1 KSP Residual norm 1.172699487022e+00 
  2 KSP Residual norm 7.080338319206e-01 
  3 KSP Residual norm 5.011871776826e-01 
  4 KSP Residual norm 3.748924659047e-01 
    :
    :
finished

[aph4499]

Yes it does work thanks. I also tried it in parallel on a test problem with a variety of Petsc options.
PETSC_OPTIONS="-options_view -ksp_monitor -ksp_view -view_sol -view_mat -ksp_type cg" srun -n 10 python ./poisson2d.py > output.dat