The left_preconditioner option in newton_solve!() does not work as expected, with current implementations of preconditioners using right_preconditioner instead.
In a cut-down version of the nonlinear_solvers.jl module, this bug was tracked down to the combinations of these lines
https://github.com/mabarnes/moment_kinetics/blob/e9e1505d0d93f462398bd1d4119db00834a289a3/moment_kinetics/src/nonlinear_solvers.jl#L1362-L1365
https://github.com/mabarnes/moment_kinetics/blob/e9e1505d0d93f462398bd1d4119db00834a289a3/moment_kinetics/src/nonlinear_solvers.jl#L1408
which applies the left-preconditioner to rhs_delta - residual0 and residual0 , respectively. This double application to residual0 is incorrect, since the second call left_preconditioner(residual0) alters the underlying data stored in residual0 permanently, even though we subsequently calculate left_preconditioner(rhs_delta - residual0) in nonlinear_solvers_approximate_Jacobian_vector_product!(). In the cut-down code, the bug was removed by using a buffer array to set the first Krylov vector, e.g.,
# To start with we use 'v' as a buffer to make a copy of residual0 to which we can apply the left-preconditioner.
@. v = residual0
left_preconditioner(v)
# Now we actually set 'w' as the first Krylov vector, and normalise it.
@. w = -v
The
left_preconditioneroption innewton_solve!()does not work as expected, with current implementations of preconditioners usingright_preconditionerinstead.In a cut-down version of the
nonlinear_solvers.jlmodule, this bug was tracked down to the combinations of these lineshttps://github.com/mabarnes/moment_kinetics/blob/e9e1505d0d93f462398bd1d4119db00834a289a3/moment_kinetics/src/nonlinear_solvers.jl#L1362-L1365
https://github.com/mabarnes/moment_kinetics/blob/e9e1505d0d93f462398bd1d4119db00834a289a3/moment_kinetics/src/nonlinear_solvers.jl#L1408
which applies the left-preconditioner to
rhs_delta - residual0andresidual0, respectively. This double application toresidual0is incorrect, since the second callleft_preconditioner(residual0)alters the underlying data stored inresidual0permanently, even though we subsequently calculateleft_preconditioner(rhs_delta - residual0)innonlinear_solvers_approximate_Jacobian_vector_product!(). In the cut-down code, the bug was removed by using a buffer array to set the first Krylov vector, e.g.,