Re: tol, atol and rtol behavior of scipy.sparse.linalg solvers

 Interesting. I have different exposure from other fields for "composing"
the combined tolerance value. There is no right way since both can be used
to arrive at the same value with the same input but with different mental
model. And we don't need to go into max-plus algebra , it's too late here
for that :)

But it must be an intuitive way since you can define arbitrary metrics. So
if both of you feel that user should expect max then nevermind max() it is.

On Tue, May 2, 2023 at 2:45 AM Robert Kern <[email protected]> wrote:

> On Mon, May 1, 2023 at 8:17 PM Robert Kern <[email protected]> wrote:
>
>> On Mon, May 1, 2023 at 7:22 PM Stefan van der Walt <[email protected]>
>> wrote:
>>
>>> On Mon, May 1, 2023, at 11:49, Ilhan Polat wrote:
>>>
>>> If atol is actually set to something, then we return max(atol,
>>> tol*norm(b)). Why max is used I don't understand.
>>>
>>> Typically, atol, rtol pair is used as tol = atol + rtol*<some metric>.
>>> This is more or less what everybody expects from this pair and their naming
>>> is chosen to reflect this. But I'm a bit lost in all the issues I could
>>> read.
>>>
>>>
>>> The first formulation — max(atol, rtol*norm(b)) — is what I'd expect.
>>> I.e., you use atol when answers are small, or rtol for larger answers.  I'm
>>> a bit surprised to see them summed; is this standard?
>>>
>>
>> I concur. `max()` is the appropriate way to combine the two that I'd
>> expect to see in almost all circumstances. `+` does have a history in our
>> ecosystem as we do use `+` in `np.isclose()`, though, but I think that is a
>> bit of a special case of comparing two numbers rather than the size of one
>> number, and I've never been really convinced by the arguments for it.
>> `math.isclose()` works like `max()` (though it is implemented with explicit
>> `||` operations).
>>
>> Julia's iterative solvers use `max()`.
>>
>>
>> https://github.com/JuliaLinearAlgebra/IterativeSolvers.jl/blob/master/src/cg.jl#L141
>>
>> MATLAB's numerical integration routines use `max()`.
>>
>>   https://www.mathworks.com/help/matlab/ref/integral.html#btc_m8o-4
>>
>
> That said, `+` is not at all uncommon outside of our own ecosystem either.
> Julia's Krylov solvers use it:
>
>   https://docs.juliahub.com/Krylov/0fcC3/0.5.5/solvers/
>
> Here's some numerical analysis course that uses it:
>
>   http://web.mit.edu/10.001/Web/Tips/Converge.htm
>
> I've never really seen any *comparison* between the two outside of some
> half-remembered arguments over `np.isclose()` early in development. I think
> `max()` is the most justified and interpretable. It's just a way to code
> the ORing of the two conditions. If either the absolute tolerance is met OR
> the relative tolerance is met, then stop. With `+`, it's a weird mix of
> both that's close to the `max()` result in most sane circumstances, but
> diverges in some cases. When the user chooses really bad choices for the
> two tolerances, neither will work particularly well, but `max()` is easier
> to reason about why it doesn't work well in that case.
>
> I strongly suspect `+` gets its lineage either from C codes or pre-F77
> codes where there was no builtin `max()` to use.
>
> --
> Robert Kern
>

_______________________________________________
SciPy-Dev mailing list -- [email protected]
To unsubscribe send an email to [email protected]
https://mail.python.org/mailman3/lists/scipy-dev.python.org/
Member address: [email protected]