Using L-infinity norm as convergence criteria when using FFT-preconditioned Conjugate Gradient Solver

[xkykai]

When we use the conjugate gradient solver and set some value of reltol and abstol, implicitly I suspect people imagine that the residuals of the solution are roughly uniformly distributed within the entire immersed domain. However, when the FFT solver is used a preconditioner, what happens is that the flow divergence occurs near the immersed boundaries, therefore the residuals of the CG solver actually concentrates near the immersed boundary.

See below for an example of a converged time step for an immersed boundary simulation:

Thus when we use a L2 norm as we are doing right now, the L2 norm would actually have much smaller values than the maximum of the residuals due to us averaging over a region with many much smaller numbers.

Perhaps it is a useful feature/default to set the convergence criteria using the L-infinity norm (ie: the maximum value) of the residuals instead? This way we guarantee no values of residuals are orders of magnitude larger than our tolerance.

@amontoison @glwagner @simone-silvestri interested in your thoughts!

[amontoison]

If you are interested by the residuals, CR or MINRES is a more relevant alternative than CG.
CR and MINRES minimizes the residual norm (||r_k||_2) at each iteration while CG minimizes the energy norm (||e_k||_A).

[xkykai]

I am specifically referring to the choice of norm in this line:

https://github.com/CliMA/Oceananigans.jl/blob/90086fd02c501dc3e41e170c2e87ff1661811b8f/src/Solvers/conjugate_gradient_solver.jl#L170C5-L170C42

where in our CG iterations we use the L2 norm of the residual as the convergence criteria.

Your suggestions of using CR and MINRES would be easy to implement once we get Krylov.jl in place! It would be incredibly exciting

[glwagner]

I suspect it's about a day of work to get that in