Re: correctly rounded mathematical functions

On Tue, Jan 4, 2022 at 3:44 AM Paul Zimmermann <[email protected]> wrote:
>
>        Dear Joel,
>
> thank you for your answer.
>
> > From: Joel Sherrill <[email protected]>
> > Date: Mon, 3 Jan 2022 14:41:20 -0600
> >
> > On Mon, Jan 3, 2022 at 6:58 AM Paul Zimmermann <[email protected]> wrote:
> > >
> > >       Dear Newlib developers,
> > >
> > > the current C working draft [1, p392] has reserved names for correctly
> > > rounded functions (cr_exp, cr_log, cr_sin, ...).
> > >
> > > We propose to provide such correctly rounded implementations
> > > for the three IEEE formats (binary32, binary64, binary128) and the
> > > "extended double" format (long double on x86_64).
> > >
> > > These implementations will be correctly rounded for all rounding modes,
> > > for example one could do the following to emulate interval arithmetic:
> > >
> > >    fesetround (FE_DOWNWARD);
> > >    y_lo = cr_exp (x_lo);
> > >    fesetround (FE_UPWARD);
> > >    y_hi = cr_exp (x_hi);
> > >
> > > Users who want a fast implementation will call the exp/log/sin/... functions,
> > > users who want a correctly rounded function and thus reproducible results
> > > (whatever the hardware, compiler or operating system) will use the
> > > cr_exp/cr_log/cr_sin/... functions. Our goal is nevertheless to get the
> > > best performance possible.
> > >
> > > Our objective is to provide open-source implementations that can be integrated
> > > in the major mathematical libraries (GNU libc, Intel Math Library, AMD Libm,
> > > Redhat Newlib, OpenLibm, Musl, llvm-libc, CUDA, ROCm).
> > >
> > > Are developers of Newlib interested by such functions?
> > > If so, we could discuss what would be the requirements for integration in
> > > Newlib in terms of license, table size, allowed operations.
> >
> > Speaking from the RTEMS perspective, we are very interested
> > in the addition of more POSIX and C Standard Library methods.
> > We have been having GSoC students add them where possible
> > for a few years.
>
> great!
>
> > The license has to be permissive and should have no advertising
> > although based on COPYING.NEWLIB, some must have a BSD
> > advertising clause still. Possibly those need review.
> >
> > https://sourceware.org/git/?p=newlib-cygwin.git;a=blob_plain;f=COPYING.NEWLIB;hb=HEAD
>
> since some Newlib files are already under LGPL v3+, for example
> include/cgen/basic-ops.h, I guess the current
> license we have (LGPL v3+) should be ok for Newlib.

I think those files are all under a linux/ or rdos/ subdirectory and do not
impact any other target unless they use iconvdata/. I grep'ed and didn't
see anything else. Jeff or Corrina should speak up and provide a
definitive answer.

> > As long as the method's footprint only impacts applications that use
> > it, there shouldn't be a huge concern on size but the target domain is
> > mostly smaller single process, no-MMU embedded systems. That is
> > ignoring Cygwin which has fewer constraints.  If you end up adding
> > megabytes, that's going to be bad in general.
>
> we'll try to have reasonable table sizes. And maybe for some functions
> we'll provide several implementations, with different table sizes.

Great! It's hard to comment on size without something to complain about. :)

> > > We have started to work on two functions (cbrt and acos), for which we
> > > provide presumably correctly rounded implementations (up to the knowledge
> > > of hard-to-round cases) [2].
> >
> > Great! What's the size profile on those?
>
> the binary32 cbrt code is quite small:
>
> -rw-r--r-- 1 zimmerma caramba   2766 Dec 23 14:27 cbrtf.c
>
> the binary128 code is larger:
>
> -rw-r--r-- 1 zimmerma caramba  38450 Dec 22 18:06 cbrtf128-v2.patch

How is binary128 implemented on target CPUs without native support?
Is this something GCC will include soft methods for?

> > What's the minimum compiler or C language version to build these?
>
> good question. So far we work with current gcc versions (for example 11.2.0)
> but did not exercise with older version. We don't use advanced C feature,
> except the "glibc patches" use some GNU libc macros (for example mul_split
> to split the product of two doubles into a+b, fast_two_sum, ...).

I didn't expect anything too fancy. Will they compile with -std=c99, c11? newlib
generally provides all methods but the header file guards protect the prototype
visibility.

glibc specific macros can be reimplemented. I assume that's going to be an
issue on other non-glibc libraries.  Hopefully that's not too hard to deal with.

--joel

> Best regards,
> Paul