[-- Attachment #1: Type: text/plain, Size: 231 bytes --] Hi, the attached patch documents the fact that the "Known Maximum Errors" given in the reference manual might not be maximal. See thread starting at https://sourceware.org/pipermail/libc-alpha/2023-December/153279.html. Paul [-- Attachment #2: 0001-document-the-fact-that-Known-Maximum-Errors-might-no.patch --] [-- Type: application/octet-stream, Size: 1042 bytes --] From f08ba79e35de677decef10982a8ecf5230b8456a Mon Sep 17 00:00:00 2001 From: Paul Zimmermann <Paul.Zimmermann@inria.fr> Date: Fri, 16 Feb 2024 09:41:14 +0100 Subject: [PATCH] document the fact that "Known Maximum Errors" might not be maximal --- manual/math.texi | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/manual/math.texi b/manual/math.texi index 2f6ee253b9..3abd83c835 100644 --- a/manual/math.texi +++ b/manual/math.texi @@ -1276,7 +1276,9 @@ extensions. @cindex ulps This section lists the known errors of the functions in the math -library. Errors are measured in ``units of the last place''. This is a +library. However, for some functions like the Bessel functions, +larger errors are known, but are documented elsewhere (bugzilla). +Errors are measured in ``units of the last place''. This is a measure for the relative error. For a number @math{z} with the representation @math{d.d@dots{}d@mul{}2^e} (we assume IEEE floating-point numbers with base 2) the ULP is represented by -- 2.34.1

On 2024-02-16 09:48:07 +0100, Paul Zimmermann wrote: > the attached patch documents the fact that the "Known Maximum Errors" > given in the reference manual might not be maximal. This section lists the known errors of the functions in the math -library. Errors are measured in ``units of the last place''. This is a +library. However, for some functions like the Bessel functions, +larger errors are known, but are documented elsewhere (bugzilla). +Errors are measured in ``units of the last place''. This is a Instead of saying that they are "documented elsewhere", shouldn't the manual be updated? BTW, the manual contains fmaf - - - - - fma - - - - - fmal - - - - - fmaf128 - - - - - fma_ldoublef - - - - - fma_ldouble - - - - - fma_ldoublel - - - - - fma_ldoublef128 - - - - - However, correct rounding is required for these functions. And I can't see any reference to f128, ldoublef, etc. in the manual except in this table. -- Vincent Lefèvre <vincent@vinc17.net> - Web: <https://www.vinc17.net/> 100% accessible validated (X)HTML - Blog: <https://www.vinc17.net/blog/> Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)

Hi Vincent, > Instead of saying that they are "documented elsewhere", shouldn't the > manual be updated? this is one solution I proposed in https://sourceware.org/pipermail/libc-alpha/2023-December/153279.html, and Carlos answered: > (a) There are known defects where ULPs may reach values that are not useful > for talking about the library in general. > > (b) There is value in being clear about the worst case known ULPs for an > implementation of a given algorithm. > > If a test is marked as XFAIL then it is clearly (a) and listing that worst > case ULPs in the manual may not be useful. thus the proposed patch documents the existence of the "XFAIL" entries. Paul

On 2024-02-16 11:02:40 +0100, Paul Zimmermann wrote: > Hi Vincent, > > > Instead of saying that they are "documented elsewhere", shouldn't the > > manual be updated? > > this is one solution I proposed in > https://sourceware.org/pipermail/libc-alpha/2023-December/153279.html, > and Carlos answered: > > > (a) There are known defects where ULPs may reach values that are not useful > > for talking about the library in general. > > > > (b) There is value in being clear about the worst case known ULPs for an > > implementation of a given algorithm. > > > > If a test is marked as XFAIL then it is clearly (a) and listing that worst > > case ULPs in the manual may not be useful. > > thus the proposed patch documents the existence of the "XFAIL" entries. If I understand correctly, Carlos means (unintended) bugs, in which case, values may even be completely wrong (as already seen). I think that saying "documented elsewhere" is misleading. And note that bugs are not specific to the math functions. Any buggy library function will not behave as described in the manual. You should rather say that the given bounds obviously do not take bugs into account, if you think that this is worth recalling. -- Vincent Lefèvre <vincent@vinc17.net> - Web: <https://www.vinc17.net/> 100% accessible validated (X)HTML - Blog: <https://www.vinc17.net/blog/> Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)

```
Vincent,
> If I understand correctly, Carlos means (unintended) bugs, in which
> case, values may even be completely wrong (as already seen). I think
> that saying "documented elsewhere" is misleading. And note that bugs
> are not specific to the math functions. Any buggy library function
> will not behave as described in the manual. You should rather say
> that the given bounds obviously do not take bugs into account, if
> you think that this is worth recalling.
of course *unknown* bugs might exist, but we are speaking of *known* defects here.
The issue is that the manual says "*Known* Maximum Errors in Math Functions".
What do you suggest to avoid users deducing from the manual that the ulp error on say y0
is bounded by 3 ulps on x86_64, whereas it can be as large as 5.93e15 ulps?
Checking y0 with glibc-2.39/build and rndn
GNU libc version: 2.39
GNU libc release: stable
y0 0 -1 0x1.c982eb8d417eap-1 [5920543797734652] [5.93e+15] 5.92055e+15 5920543797734652
libm gives -0x1.8p-55
mpfr gives -0x1.af74bfa0f1304p-56
Paul
```

On 2024-02-16 11:54:03 +0100, Paul Zimmermann wrote: > Vincent, > > > If I understand correctly, Carlos means (unintended) bugs, in which > > case, values may even be completely wrong (as already seen). I think > > that saying "documented elsewhere" is misleading. And note that bugs > > are not specific to the math functions. Any buggy library function > > will not behave as described in the manual. You should rather say > > that the given bounds obviously do not take bugs into account, if > > you think that this is worth recalling. > > of course *unknown* bugs might exist, but we are speaking of *known* > defects here. Well, known bugs in non-math functions are never documented either. For instance, concerning regexps, there is bug 10844[*], which has been known since 2009 (more than 14 years ago) and is not documented. [*] https://sourceware.org/bugzilla/show_bug.cgi?id=10844 > The issue is that the manual says "*Known* Maximum Errors in Math Functions". > > What do you suggest to avoid users deducing from the manual that the ulp error on say y0 > is bounded by 3 ulps on x86_64, whereas it can be as large as 5.93e15 ulps? > > Checking y0 with glibc-2.39/build and rndn > GNU libc version: 2.39 > GNU libc release: stable > y0 0 -1 0x1.c982eb8d417eap-1 [5920543797734652] [5.93e+15] 5.92055e+15 5920543797734652 > libm gives -0x1.8p-55 > mpfr gives -0x1.af74bfa0f1304p-56 If I understand correctly, the *intent* is to give an error bounded by 3 ulp, i.e. an error up to 3 ulp is not something that the user should expect to be "fixed". Perhaps the glibc manual should provide a way to find known bugs more easily. That said, typing "y0" in bugzilla immediately gives https://sourceware.org/bugzilla/show_bug.cgi?id=16492 (for some other functions, which are also common words, this is not that obvious). -- Vincent Lefèvre <vincent@vinc17.net> - Web: <https://www.vinc17.net/> 100% accessible validated (X)HTML - Blog: <https://www.vinc17.net/blog/> Work: CR INRIA - computer arithmetic / AriC project (LIP, ENS-Lyon)

On Fri, 16 Feb 2024, Vincent Lefevre wrote: > BTW, the manual contains > > fmaf - - - - - > fma - - - - - > fmal - - - - - > fmaf128 - - - - - > fma_ldoublef - - - - - > fma_ldouble - - - - - > fma_ldoublel - - - - - > fma_ldoublef128 - - - - - > > However, correct rounding is required for these functions. IBM long double doesn't have such a thing as correct rounding, so we allow errors of up to 3ulp (the documented accuracy for IBM long double division) for functions expected to be correctly rounded for otehr formats. > And I can't see any reference to f128, ldoublef, etc. in the manual > except in this table. See @deftypefunx _FloatN fmafN (_Float@var{N} @var{x}, _Float@var{N} @var{y}, _Float@var{N} @var{z}) @deftypefunx _FloatNx fmafNx (_Float@var{N}x @var{x}, _Float@var{N}x @var{y}, _Float@var{N}x @var{z}) in the manual for fmaf128 (which should always be correctly rounded, since _Float128 is always IEEE binary128). The "_ldouble" suffix is how libm-test-ulps indicates the argument type for narrowing functions that are parametrized by two types (that is, it's an implementation detail that ideally wouldn't appear in the manual at all). IBM long double arguments (thus, float or double returns - functions {f,d}{add,sub,mul,div,sqrt,fma}l) is the only case where the narrowing functions should ever have results that are not correctly rounded. -- Joseph S. Myers josmyers@redhat.com

```
On Fri, 16 Feb 2024, Vincent Lefevre wrote:
> If I understand correctly, the *intent* is to give an error bounded
> by 3 ulp, i.e. an error up to 3 ulp is not something that the user
> should expect to be "fixed".
The intent is that the functions should be both "accurate" (where more
details of the precise goals are given in the manual - see the first
several paragraphs of the "Errors in Math Functions" section) and "fast".
3 ulp isn't a specific intent. The testsuite always treats any error
above 9 ulp (16 ulp for IBM long double) as a bug. In practice I think
it's likely that most or all of the functions with errors above 1 ulp for
round-to-nearest could be improved (other than for IBM long double) to
have errors of at most 1 ulp without making them any slower, using better
implementations. (The definition of ulp used by the testsuite and the
table in the manual compares to the correctly rounded result, *not* to the
mathematically precise result. So fractional ulp values are not possible,
except for the special case where the result returned by the function has
a lower exponent than the correctly rounded result.)
(The main cases that would be genuinely hard to get down to low errors for
all arguments are the cpow, jn and yn functions.)
--
Joseph S. Myers
josmyers@redhat.com
```