Linear solver interfaces

[jlogan03]

I have a rough working prototype of serial, sparse-sparse BiCGSTAB here and it's got me thinking about how to structure interfaces to solvers.

It seems like there are a few specific actions that users need regularly that might be supported in different ways by different kinds of solvers -

• Just solve a system with minimal fiddling
• Prep a solver, keeping factorization, preconditioners, etc. to use for solving a system multiple times with different right-hand side
• Advance an iterative solver, then do something with the latest error estimate and solution
  • Convergence plots, checkpointing, etc
• Restart an iterative solver's internal state, preserving the latest solution
• Swap components of an initialized solver with minimum re-computation
  • Like keeping a symbolic factorization and using it to swap in a new A matrix for FEM with deformation
  • ...or swapping a new A or rhs into a partially-converged iterative solver to capture a nonlinear update

To that end, I'm thinking each solver should be associated with both

• A standalone function that provides the fuss-free solve, and
• A struct that keeps a workspace and implements methods for manipulating that workspace

All pretty vanilla so far.

Now, there's a pretty different pattern of usage between direct and iterative solvers, and I can't think of a way to stuff them into the same trait. So with that, we'd have two traits, say, DirectSolver and IterativeSolver, guiding the implementation of the solver structs, and possibly another two guiding the implementation of the convenience functions.

This leaves us with

• Two traits guiding implementation of direct and iterative solvers
• Optionally, another two traits or type aliases guiding the convenience functions
• For each solver,
  • A workspace struct implementing either DirectSolver or IterativeSolver
  • A convenience function, which either
    • is part of the impl on the struct
    • implements a trait or matches a provided type alias
    • is freeform and unconstrained in its construction

@mulimoen - I'm curious what your thoughts are on

• Whether now is a good time to implement such a structure, or whether new solvers should be added independently for now
• Whether there are important use-cases that I haven't captured
• Whether the convenience functions should be standalone (and maybe have their own traits or type aliases) or whether they should be part of the impl on the solver workspace structs
• Any other thoughts about the proposed structure

[mulimoen]

The implementation of BiCGSTAB would be great to have in this crate.

On generalisation into traits I would hold on a bit until we know what to expect on differences between the solvers. Maybe each solver should return residuals or other temporaries? Is the input similar between solvers and how does preconditioners fit in? It would be great to see what functionality is common between the solvers before settling on some traits.

[jlogan03]

Sounds good - I'll button that implementation up with my best guess at what an iterative solver interface might look like, but without encoding that interface in a trait just yet.