Error in indirect KKT solver: `tol` keyword in IterativeSolvers.jl deprecated about cosmo.jl HOT 6 CLOSED

Closed with PR above

from cosmo.jl.

Good question. It's not clear to me whether it should be abstol or reltol there. Note that a few lines below we have :

where the tolerance is divided by the initial residual. Given the documentation for cg! takes as a convergence test |r_k| ≤ max(reltol * |r_0|, abstol) (see here), it's not really clear to me what we are expecting to happen.

Another problem in the same place is that we probably shouldn't be normalizing by norm(y1) without some guard against a zero norm condition, i.e. something like max(norm(y1),1.) seems safer.

Not sure who wrote this originally - paging @migarstka and @nrontsis...

from cosmo.jl.

Ah good point. I mean that changing to reltol would give the same behavior as you had previously. However, it looks like use abstol should be used directly and the normalization should be removed since IterativeSolvers.jl multiplies by the norm of the initial residual?

Compare https://github.com/JuliaLinearAlgebra/IterativeSolvers.jl/blob/c91cfbc5e80c4451d4b0c7bc9433e5302915cccc/src/cg.jl#L141

with
https://github.com/JuliaLinearAlgebra/IterativeSolvers.jl/blob/49cc093d60babbf7a3edd0119958d93617771406/src/cg.jl#L132

(also the warning in v0.8)
https://github.com/JuliaLinearAlgebra/IterativeSolvers.jl/blob/152f89f39db9c693cf283a8dc7b6c75be0499200/src/cg.jl#L134

from cosmo.jl.

Having looked at the links above and some older versions of IterativeSolvers.jl (i.e. pre v0.8), I think the situation is:

• v0.7 only had a "tol" value, which was treated as a relative tolerance internally (see here).
• v0.8 added a warning if the caller passed "tol", preferring the name "reltol" and remapping it internally (as in the warning cited above).
• current version uses "reltol" in the same way, but with the addition of an "abstol" that acts as a hard floor on the tolerance.

Assuming that the above is correct then the purpose of the current code is basically to act as a hack to force the use of an absolute tolerance term even where it was no previously supported. That is compatible (I think) with the approach in this paper, where we want the linear solve errors to be absolutely summable.

To match the original behaviour, do we therefore need (...abstol=get_tolerance(S), reltol = 0.) ? That seems correct (and agrees withs @tjdiamandis above I think) but also feels like it is ill-advised. Maybe just letting "reltol" take its default would be preferable in practice.

from cosmo.jl.

Not sure who wrote this originally - paging @migarstka and @nrontsis...

$ git blame -L70 -- src/linear_solver/kktsolver_indirect.jl

8553294f src/kktsolver_indirect.jl (nrontsis     2019-09-02 17:55:25 +0100  70)         cg!(S.previous_solution, L, y1, tol=get_tolerance(S)/norm(y1))

One for you @nrontsis 😋

I think we should aim to support the latest stable release of IterativeSolvers.jl which seems to be v0.9+ and also add that constraint to Project.toml file. (I remember, we have set this up as an optional dependency that is only usable if the user deliberately installs IterativeSolvers. I don't know if you can set a version constraint on an optional dependency.)

I would go with:

cg!(S.previous_solution, L, y1, atol=get_tolerance(S)/max(norm(y1),1.))

Is anyone willing to create a PR for this? We should also check whether the basic tests here actually run.

from cosmo.jl.

I think this was a hack to get around an issue of IterativeSolvers.jl at the time.

I think we should go with:

(...abstol=get_tolerance(S), reltol = 0.)

that Paul suggested. Having this with IterativeSolvers.jl 0.9+, might actually not be a breaking change as compared to the behaviour I had when I committed the line under question (8553294f), except when |r_0| = 0. I can write a simple "test" that checks exact equivalence on toy data.

but also feels like it is ill-advised. Maybe just letting "reltol" take its default would be preferable in practice.

Perhaps reltol=0 might be okay-ish, given that COSMO.jl scales the problem data by default?

from cosmo.jl.