nekStab / LightKrylov

Lightweight implementation of Krylov subspace techniques in Fortran.
BSD 3-Clause "New" or "Revised" License
18 stars 0 forks source link

Abstract opts type for linear solvers to have a standardized interface #25

Closed loiseaujc closed 1 year ago

loiseaujc commented 1 year ago

The interface for gmres reads

gmres(A, b, x, info, kdim, maxiter, tol, verbosity, transpose)

while that for bicgstab (before my modifications in #24) is

bicgstab(A, b, x, info, maxiter, tol, verbosity, transpose)

Likewise, the interface for cg is

cg(A, b, x, info, maxiter, tol, verbosity)

It would be nice to somehow standardise these interfaces, most notably if we include additional linear solvers later on which might require the definition of new input parameters. How about some abstract_opts_type that we extend for each different solvers? Taking gmres as an example, the header of the subroutine would look like

subroutine gmres(A, x, b, info, transpose, opts)
    !> Linear problem.
    class(abstract_linop)  , intent(in) :: A
    class(abstract_vector), intent(in) :: b
    !> Solution vector.
    class(abstract_vector), intent(inout) :: x
    !> Information flag.
    integer, intent(out) :: info
    !> Transpose flag.
    logical, optional, intent(in) :: transpose
    !> Options.
    type(gmres_opts), intent(in) :: opts

    ...

    return
end subroutine gmres

where gmres_opts is a type extending abstract_opts_type and containing the standard fields such as kdim, maxiter, tol and verbosity. The interface for all other linear solvers would then look the same, i.e.

solver(A, x, b, info, transpose, opts)

That would be particularly useful in the rational_arnoldi_factorization for instance where I pass a subroutine solve as input argument to solve the problem $\left( \mathbf{I} - \mathbf{A}/\sigma \right) \mathbf{x} = \mathbf{b}$. Having a standardized interface would allow us to pass any linear solver we want, even those that aren't yet implemented.

What do you think?

loiseaujc commented 1 year ago

This is probably the easiest issue among #27, #26 and this one. I'll start working on it later tonight. I think I can get it done in an hour or so. I will create a dedicated branch and send a PR as soon as I'm done.