### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 3:34am UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/28 "2024-09-19T03:34:18Z")
> [@Richard\_Fateman](#):
>
> You are trying to grapple with the question of how to compare floating-point numbers, hoping that there a cure for “isapprox”. Lots of responses.
> Here’s mine.
> Don’t use it.
I don’t agree with a blanket statement like this. The most common application of approximate equality, by far, seems to be in software testing, where it is crucially important and appropriate. For example, if I’ve implemented a function to compute a Bessel function, I might want a regression test like:
```julia
@test mybesselj(3, some_x) ≈ known_result_for_x rtol=1e-13
```
to ensure that the code stays accurate to at least 13 significant digits in future revisions. Or if you write an eigensolver you might want to check the backwards accuracy `A ≈ X * Λ / Χ` to some tolerance for some `A` (note that the use of vector norms in `isapprox` is critical to do array comparisons like this reliably).
> [@Richard\_Fateman](#):
>
> Consider if what you intend by a particular use of isapprox(a,b) is
> -eps +a \< b \< eps+a or
> equivalently abs(a-b)\ where eps can be scaled,say by max(abs(a), abs(b),machine-epsilon) … or some other factor, depending on what you actually need
That’s pretty much [what `isapprox` does](https://docs.julialang.org/en/v1/base/math/#Base.isapprox), except that your pseudo-code `max(abs(a), abs(b),machine-epsilon)` seems to be falling into the classic mistake of treating the machine epsilon as an absolute tolerance rather than a relative tolerance.
(You seem to be advising people to roll their own approximate comparisons — which experience shows many people do badly — rather than just call `isapprox` and pass context-appropriate tolerances as needed? This advice mystifies me.)
But I absolutely agree that it’s best to think carefully about what tolerances you want to use an a given application, and to pass them explicitly to `isapprox`. The default `rtol` in `isapprox` is reasonable for a lot of simple `@test a ≈ b` checks where you just want to catch crude bugs (and aren’t comparing to zero): checking whether half the significant digits match will catch severe bugs but ignores most roundoff issues, but it is definitely better to put some thought into it.
> [@jlapeyre](#):
>
> It turns out that this design can be completely at odds with your intuition.
I don’t understand why you think the result is non-intuitive here. Why would you be doing 77-digit arithmetic if you don’t want the computer to consider the 35th digit to be significant by default (without a manually specified tolerance)?
But yes, `isapprox` cannot be a magical [DWIM](https://en.wikipedia.org/wiki/DWIM) function that always knows the tolerance you intend.
---
### Author:  [@Richard\_Fateman](https://discourse.julialang.org/u/Richard_Fateman)
#### Post date: [September 19, 2024, 5:11am UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/29 "2024-09-19T05:11:06Z")
If you are using isapprox to test for accuracy of software, that
might be fairly misleading. Why? Your program might be tolerable at random places, and therefore pass. But it might be really bad in edge cases. For instance, near zeroes of Bessel functions, where the calculations may
be more difficult and possibly more important.
If you are testing for nearness of floats as a termination condition for an iterative process, you should have some notion of what kind of convergence is desired. If isapprox(), whatever it does, floats your boat, so to speak,
it seems to me that nevertheless a subsequent reader of the code may find it more illuminating to see something explicit like abs(a-b)\< some\_epsilon.
The IEEE 754 spec also has denormalized numbers / gradual underflow,
which affects how small an abs(a-b) can be computed. Presumably smaller than machine-epsilon. In the past, higher-level languages tend to lack support for such notions, making life more difficult for error analysts. Perhaps newer languages are better.
---
### Author:  [@jlapeyre](https://discourse.julialang.org/u/jlapeyre)
#### Post date: [September 19, 2024, 11:52am UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/30 "2024-09-19T11:52:43Z")
> [@stevengj](#):
>
> I don’t understand why you think the result is non-intuitive here. Why would you be doing 77-digit arithmetic if you don’t want the computer to consider the 35th digit to be significant by default (without a manually specified tolerance)?
When I said “you” I wasn’t thinking of the most experienced users. I was thinking of a user whose first thought is that `≈` has the everyday meaning “about the same”. Based on experience, I strongly suspect that there are a lot of these users. But that’s really a question for pollsters, I suppose.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 12:31pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/31 "2024-09-19T12:31:37Z")
> [@Richard\_Fateman](#):
>
> If you are using isapprox to test for accuracy of software, that
> might be fairly misleading. Why? Your program might be tolerable at random places, and therefore pass. But it might be really bad in edge cases. For instance, near zeroes of Bessel functions, where the calculations may
> be more difficult and possibly more important.
When you test special function implementations, you don’t just pick random places to compare. You have to specially test corner case like zeros, or places where your internal approximation changes. (Indeed, for _any_ function you want to ensure your tests are designed to cover all branches and corner/extreme cases.) But the test is still going to look like `isapprox` with carefully chosen tolerances. And for most test cases what you want is a relative tolerance.
And, as I said, if you roll your own approximate comparison from scratch in every case (as you seem to be advocating), many people get key things wrong, e.g. treating machine precision as an absolute tolerance, or comparing arrays elementwise rather than in an overall norm. I don’t understand why you seem to think that every testing code should rewrite this basic function over and over again.
For example, for tests of array functions like FFTs where you [_can_ often use carefully chosen randomized tests](https://epubs.siam.org/doi/10.1137/S0097539796311168), it’s still important to test things with vector norms (as is done in `isapprox`) so that you use a tolerance based on the overall scale of the vector. (Elementwise checks, in contrast, tend to randomly fail if one element of the array is much smaller than the others.)
> If `isapprox()`, whatever it does, floats your boat, so to speak,
You’re writing as though `isapprox` does something mysterious and unknowable? How it does the comparison [is precisely documented](https://docs.julialang.org/en/v1/base/math/#Base.isapprox), with user-controllable tolerances.
I agree that if you have an iterative solver with a convergence criteria is as simple as `norm(Δx) ≤ tol*scale`, it’s perfectly reasonable to write that explicitly, especially if performance is critical (especially on vectors — `isapprox` might call a `norm` or allocate a temporary vector that you can avoid), though this is not a concern in software testing.
> [@Richard\_Fateman](#):
>
> The IEEE 754 spec also has denormalized numbers / gradual underflow,
> which affects how small an abs(a-b) can be computed. Presumably smaller than machine-epsilon. In the past, higher-level languages tend to lack support for such notions, making life more difficult for error analysts. Perhaps newer languages are better.
I’m not sure what higher-level languages “lack support” for subnormal values these days? It is [possible to flush them to zero](https://docs.julialang.org/en/v1/manual/performance-tips/#Treat-Subnormal-Numbers-as-Zeros), but do any languages do this by default? Or are you talking about symbolic/computer-algebra systems that default to arbitrary-precision arithmetic rather than hardware IEEE arithmetic? In any case, Julia certainly supports them fine (and does not flush them by default), as long as you don’t want to [trap the underflow exceptions](https://github.com/JuliaLang/julia/issues/27705).
But I’m not sure what your point is. Subnormals kick in nearly 300 orders of magnitude below machine epsilon in double precision, so for most calculations (i.e. calculations using most of the 600 orders of magnitude supported in double precision) the [threshold for subnormal numbers](https://docs.julialang.org/en/v1/base/base/#Base.floatmin) is irrelevant for approximate comparisons.
(And yes, in many mathematical functions you might want to specially test the corner case of behaving gracefully in underflow/overflow scenarios. It is even more common to want to make sure you don’t accidentally underflow/overflow, a classic example being the naive algorithm `sqrt(x^2+y^2)` for [`hypot(x,y)`](https://docs.julialang.org/en/v1/base/math/#Base.Math.hypot). But such tests still involve `isapprox` or equivalent, with due care about the tolerances.)
---
### Author:  [@dlakelan](https://discourse.julialang.org/u/dlakelan)
#### Post date: [September 19, 2024, 1:13pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/32 "2024-09-19T13:13:07Z")
Slightly off topic, but welcome Richard! For those who don’t know, he’s a major author of portions of the Maxima computer algebra system and professor emeritus at Berkeley specializing in computing and applied math. Glad to see you contributing to the Julia ecosystem.
---
### Author:  [@BioTurboNick](https://discourse.julialang.org/u/BioTurboNick)
#### Post date: [September 19, 2024, 1:59pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/33 "2024-09-19T13:59:14Z")
> [@stevengj](#):
>
> (You seem to be advising people to roll their own approximate comparisons — which experience shows many people do badly — rather than just call `isapprox` and pass context-appropriate tolerances as needed? This advice mystifies me.)
No one in this conversation is suggesting people can’t pass context-appropriate tolerances as needed. It’s a question of the best default tolerances. You keep asserting any other default is wrong. I urge you to consider that is a matter of perspective and usage, and that the people using this function with its default arguments probably aren’t the people who are using extremely high and extremely low magnitude numbers.
And with that, I’m out.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 2:01pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/34 "2024-09-19T14:01:22Z")
> [@stevengj](#):
>
> The most common application of approximate equality, by far, seems to be in software testing, where it is crucially important and appropriate.
I should also mention another scenario where I use `≈` (i.e. `isapprox` with default tolerances) all the time, both for my own purposes and for pedagogy: quick checks of formulas and conjectures.
For example, suppose I’ve worked out a nontrivial formula, like \mathrm{vec}(BCA^T) = (A \otimes B) \mathrm{vec}(C), and I want to make sure I didn’t make some simple algebra mistake. The first thing I do is something like `B = randn(10,20); C = randn(20,30); A = randn(40,30); vec(B*C*A') ≈ kron(A, B) * vec(C)` and look for a `true` result.
I find it incredibly convenient (and clearer) that I don’t have to write `norm(vec(B*C*A') - kron(A, B) * vec(C)) < 1e-8 * norm(kron(A, B) * vec(C))` over and over again for these kinds of things. And manually inspecting the output is a pain, especially for arrays. I don’t need a custom tolerance since I’m not trying to check for floating-point issues, only algebra goofs, and an 8-digit match is usually enough to rule out the latter with high probability, while being insensitive to the former. But I want a relative tolerance because my function might not be scaled to values near 1 (though in this example it is).
PS. `≈` is especially quick since the Mac keyboard has a shortcut option-x for it, though you can also tab-complete `\approx`.
---
### Author:  [@mcabbott](https://discourse.julialang.org/u/mcabbott)
#### Post date: [September 19, 2024, 2:10pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/35 "2024-09-19T14:10:34Z")
> [@stevengj](#):
>
> The problem here is that `x ≈ 0` is not just using `≈` wrong “as it is defined” — it fundamentally makes no sense in the context of floating-point arithmetic.
How crazy would it be for this to print a warning on first use? If called with the default tolerance. With an explanation of why this is a bad idea (or a link to one). It’s correct once you think about it, but certainly catches people out.
> [@stevengj](#):
>
> PS. `≈` is especially quick since the Mac keyboard has a shortcut option-x for it
TIL, thanks!
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 2:21pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/36 "2024-09-19T14:21:11Z")
> [@BioTurboNick](#):
>
> No one in this conversation is suggesting people can’t pass context-appropriate tolerances as needed.
I’m specifically responding to @Richard_Fateman’s suggestion that people shouldn’t use `isapprox` (with or without custom tolerances) and should rather just write ad-hoc comparisons from scratch in most cases.
> [@BioTurboNick](#):
>
> You keep asserting any other default is wrong. I urge you to consider that is a matter of perspective and usage, and that the people using this function with its default arguments probably aren’t the people who are using extremely high and extremely low magnitude numbers.
It is precisely _because_ `atol` depends on context/usage/perspective that we can’t assume a nonzero default.
`Float64` has about 600 orders of magnitude available, but even scaling by `1e±6` (which is not very extreme) will render an `atol=1e-8` useless or worse. Yes, I think that having a basic mathematical function assume that all numbers are scaled to unit magnitude is wrong. Floating-point numbers are not fixed-point numbers, though many people seem to confuse the two.
(Even if your numbers start out around 1, e.g. if you generate them with `randn`, they can quickly become much larger or much smaller depending on what functions you pass them through.)
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 2:24pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/37 "2024-09-19T14:24:08Z")
> [@mcabbott](#):
>
> How crazy would it be for this to print a warning on first use [of `x ≈ 0`]?
I think it could be reasonable for the `@test` macro to print a suggestion on an `x ≈ 0` test that fails, since it already prints a custom message for failed comparisons: [more informative @test x ≈ y error message · Issue #55812 · JuliaLang/julia · GitHub](https://github.com/JuliaLang/julia/issues/55812)
Probably not the `isapprox` function itself (which would be a breaking change, and in any case it’s not usual behavior for mathematical functions to print warnings like this).
---
### Author:  [@Jean\_Michel](https://discourse.julialang.org/u/Jean_Michel)
#### Post date: [September 19, 2024, 3:00pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/38 "2024-09-19T15:00:08Z")
> [@stevengj](#):
>
> For example, suppose I’ve worked out a nontrivial formula, like \mathrm{vec}(BCA^T) = (A \otimes B) \mathrm{vec}(C)vec(BCAT)=(A⊗B)vec(C)\mathrm{vec}(BCA^T) = (A \otimes B) \mathrm{vec}(C), and I want to make sure I didn’t make some simple algebra mistake. The first thing I do is something like `B = randn(10,20); C = randn(20,30); A = randn(40,30); vec(B*C*A') ≈ kron(A, B) * vec(C)` and look for a `true` result.
The first thing I would do is
`B = rand(Int,10,20); C = rand(Int,20,30); A = rand(Int,40,30);`
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 3:00pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/39 "2024-09-19T15:00:13Z")
> [@stevengj](#):
>
> I find it incredibly convenient (and clearer) that I don’t have to write `norm(vec(B*C*A') - kron(A, B) * vec(C)) < 1e-8 * norm(kron(A, B) * vec(C))` over and over again for these kinds of things.
Aside: one function that I _do_ find myself writing over and over is:
```julia
relerr(a, b) = norm(a - b) / max(norm(a), norm(b))
```
in order to inspect the relative error of some calculation. I wish this were built-in. I suppose I could put it in my `startup.jl`, but then it will bite me every time I go to a new machine or try to share code snippets with someone else.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 19, 2024, 3:03pm UTC](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347/40 "2024-09-19T15:03:00Z")
> [@Jean\_Michel](#):
>
> The first thing I would do is
> `B = randn(Int,10,20); C = randn(Int,20,30); A = randn(Int,40,30);`
Sure, in this _particular_ example the function maps integers to integers. (I think you mean `rand(Int, m, n)`. Luckily, this particular function is also insensitive to integer overflow.)
But that isn’t true of most functions.
[Previous page](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347.md?page=1)
[Next page](https://discourse.julialang.org/t/history-of-isapprox-in-languages/119347.md?page=3)