Crash when specifying "linear_solver = iterative"

[Marie Rognes]

The adjoint computation below fails. If "linear_solver" is "gmres", "cg" or "lu" all is well.

from dolfin import *
from dolfin_adjoint import *

mesh = UnitSquareMesh(2, 2)
V = FunctionSpace(mesh, "CG", 1)
f = Constant(1.0)
u = project(f, V)
v = TestFunction(V)

F = u*u*v*dx - f*v*dx

info_green("Running forward")
pde = NonlinearVariationalProblem(F, u, J=derivative(F, u))
solver = NonlinearVariationalSolver(pde)
solver.parameters["linear_solver"] = "iterative"
solver.solve()

J = Functional(u*u*dx*dt[FINISH_TIME])

info_green("Computing adjoint")
z = compute_adjoint(J)
for x in z:
    print x

[Patrick Farrell (pefarrell)]

I'm going to pin this one on dolfin. Because there is no identity representation in UFL, we have to use KrylovSolver/LUSolver even if solve(a == L), solve(F == 0), or the OO interface is used in the forward model. This is caused by the inconsistent treatment of "iterative" across the LA interfaces in dolfin.

There are at least three opinions about what "iterative" means in dolfin. In dolfin/fem/LinearVariationalSolver.cpp:

    if (solver_type == "iterative" || solver_type == "krylov")
    {
      // Adjust iterative solver type
      if (symmetric)
        solver_type = "cg";
      else
        solver_type = "gmres";
    }

In dolfin/la/LinearSolver.cpp:

  else if (method == "iterative")
    method = "gmres";

(no mention of "cg")

and dolfin/la/KrylovSolver.cpp crashes if you pass it method="iterative" (which is what is happening here):

  File "/home/pef/src/dolfin-adjoint/dolfin_adjoint/adjlinalg.py", line 399, in wrap_solve
    solver = dolfin.KrylovSolver(method, pc)
*** -------------------------------------------------------------------------
*** Error: Unable to solve linear system using Krylov iteration.
*** Reason: Unknown Krylov method "iterative". Use list_krylov_solver_methods() to list available Krylov methods.

I could patch over this in dolfin-adjoint, but the Right Thing To Do is fix it in dolfin (preferably by taking the LinearVariationalSolver logic and porting it to LinearSolver, KrylovSolver and anywhere else it's needed, like the backend-specific objects).

[Marie Rognes (meg-simula)]

Agree that the "iterative" keyword is the design flaw. It is easy to get around this by specifying the linear solver directly instead ("cg" or "gmres" for instance), and that reads much clearer too.