Elementary Function Work

Elementary Function Work

[N.G. Timmons]

Hi all,

This was meant to be sent in February but the year really got away from 
me for a bit.

I sent this email to Carlos O'Donell:

==============================================================================
I was looking through the source for sin/cos/exp etc. and noticed that a
lot of the comments don't actually say what/why things are being done.
There also exists a union type called "mynumber" which is used to access
the high and low bits of a 64-bit floating point number.

I would like to comment the elementary functions with:
      - What algorithm is being used
      - Why is it being used
      - Any caveats on unusual operations so future changes can be made
safely
      - Actual testing information (the current statistics given at the
top of some functions doesn't say how they were tested). This is to
allow for verification and comparison if future changes are made.

I would like to rename 'mynumber' and a few other variables to
meaningful names so that the source code can be followed with greater
ease.

Would the community be open to me making these changes?
==============================================================================

They suggested that I post to this group to see if I could make some of 
these changes?

Thanks,

Nick Timmons
Digital Technology Group
University of Cambridge

Re: Elementary Function Work

[Joseph Myers]

On Mon, 7 Dec 2020, N.G. Timmons via Libc-alpha wrote:

> I was looking through the source for sin/cos/exp etc. and noticed that a
> lot of the comments don't actually say what/why things are being done.

exp is now one of the optimized Arm implementations, which have detailed 
comments.

sin and cos for double are from the much older IBM libm.  It's quite 
possible they could be replaced by new implementations along the lines of 
those for float and be made significantly faster in the process.

> There also exists a union type called "mynumber" which is used to access
> the high and low bits of a 64-bit floating point number.

Where this is used to initialize tables of floating-point values using 
specific representations, with #ifdef cases for the endianness, it would 
be better to replace uses by uses of double, with hex float constants 
(verifying that installed stripped shared libraries are unchanged in the 
process).

>      - Actual testing information (the current statistics given at the
> top of some functions doesn't say how they were tested). This is to
> allow for verification and comparison if future changes are made.

Typically testing information (beyond the glibc testsuite) would go in the 
commit message for a change, describing how that change was tested, not in 
the source code.  Comments need to be kept up to date to reflect the 
version of the source code they appear in, saying "tested with the glibc 
testsuite" is redundant while describing some other tests that were run 
either (a) imposes a load on everyone modifying that code in future to 
rerun exactly those tests, (b) leaves information in the source code that 
becomes out of date as soon as changes are made, or (c) has to be removed 
when changes are made as it may not describe how that changed version was 
tested.

-- 
Joseph S. Myers
joseph@codesourcery.com

Re: Elementary Function Work

[N.G. Timmons]

Hi,

Apologies, I should have been more clear in my original e-mail.

>> I was looking through the source for sin/cos/exp etc. and noticed that 
>> a
>> lot of the comments don't actually say what/why things are being done.

For this I was referring to the 32-bit float implementations, but in 
general I was thinking about describing the algorithm being used so that 
people can understand if they need to check the source for any reason.

I really like the idea of new implementations for double or float that 
would improve performance - if I could help with that I would be very 
grateful.

>> There also exists a union type called "mynumber" which is used to 
>> access
>> the high and low bits of a 64-bit floating point number.

For this one I was initially interested in renaming it to something 
meaningful because when I first looked through the code I was really 
confused when I first saw it and navigating the source (especially 
online) can be a bit tricky to find the definition.

> Typically testing information (beyond the glibc testsuite) would go in 
> the
> commit message for a change, describing how that change was tested, not 
> in
> the source code.  Comments need to be kept up to date to reflect the
> version of the source code they appear in, saying "tested with the 
> glibc
> testsuite" is redundant while describing some other tests that were run
> either (a) imposes a load on everyone modifying that code in future to
> rerun exactly those tests, (b) leaves information in the source code 
> that
> becomes out of date as soon as changes are made, or (c) has to be 
> removed
> when changes are made as it may not describe how that changed version 
> was
> tested.

This makes sense.
I would be interested in setting up a test suite to provide ULP error 
information, error distribution information and number of known 
incorrectly rounded results that could be run to generate a 
documentation page or automatically update the function comment with the 
information for easy accessibility.
This is of particular interest to me as I am studying for my PhD in 
approximate computation and have found that different libm 
implementations have different error distributions and mismatch counts 
which affects interoperability. Part of my research output would be 
raising awareness of this. So if this information would be valuable to 
this project, it would be great to be able to implement this helper 
tool.

Thanks
Nick



On 2020-12-07 19:19, Joseph Myers wrote:
> On Mon, 7 Dec 2020, N.G. Timmons via Libc-alpha wrote:
> 
>> I was looking through the source for sin/cos/exp etc. and noticed that 
>> a
>> lot of the comments don't actually say what/why things are being done.
> 
> exp is now one of the optimized Arm implementations, which have 
> detailed
> comments.
> 
> sin and cos for double are from the much older IBM libm.  It's quite
> possible they could be replaced by new implementations along the lines 
> of
> those for float and be made significantly faster in the process.
> 
>> There also exists a union type called "mynumber" which is used to 
>> access
>> the high and low bits of a 64-bit floating point number.
> 
> Where this is used to initialize tables of floating-point values using
> specific representations, with #ifdef cases for the endianness, it 
> would
> be better to replace uses by uses of double, with hex float constants
> (verifying that installed stripped shared libraries are unchanged in 
> the
> process).
> 
>>      - Actual testing information (the current statistics given at the
>> top of some functions doesn't say how they were tested). This is to
>> allow for verification and comparison if future changes are made.
> 
> Typically testing information (beyond the glibc testsuite) would go in 
> the
> commit message for a change, describing how that change was tested, not 
> in
> the source code.  Comments need to be kept up to date to reflect the
> version of the source code they appear in, saying "tested with the 
> glibc
> testsuite" is redundant while describing some other tests that were run
> either (a) imposes a load on everyone modifying that code in future to
> rerun exactly those tests, (b) leaves information in the source code 
> that
> becomes out of date as soon as changes are made, or (c) has to be 
> removed
> when changes are made as it may not describe how that changed version 
> was
> tested.

Re: Elementary Function Work

[Joseph Myers]

On Tue, 8 Dec 2020, N.G. Timmons wrote:

> I would be interested in setting up a test suite to provide ULP error
> information, error distribution information and number of known incorrectly
> rounded results that could be run to generate a documentation page or
> automatically update the function comment with the information for easy
> accessibility.

I suggest working with Paul Zimmermann since he's been working in that 
area, and it's naturally something separate from any particular libm 
implementation.  (Though different libm implementations may make different 
choices regarding the bounds on acceptable errors, requirements for 
correctness of floating-point exceptions raised, etc.)

https://homepages.loria.fr/PZimmermann/papers/#accuracy

There are obvious areas for extending his work, such as covering the 
ldbl-96 and ldbl-128ibm formats and <complex.h> functions and different 
rounding modes.

When function errors are found that are beyond the limits of those 
accepted for glibc, it's appropriate to file bug reports when such bug 
reports aren't already open for the functions in question (see his recent 
bug reports for lgamma and tgamma).  When function errors are found that 
are within the bounds of those accepted for glibc, but larger than those 
appearing in libm-test-ulps, it's appropriate to add the test inputs in 
question to auto-libm-test-in and regenerate the corresponding 
auto-libm-test-out-* file, so that the libm-test-ulps files (and thus the 
tables in the glibc manual generated from those files) more accurately 
reflect the largest known errors.

If a function gets reimplemented so as to reduce the largest errors (for 
most functions with errors larger than 1 or 2 ulps, I expect that could be 
done without making performance any worse), it then makes sense to remove 
entries for that function in libm-test-ulps files so they can be 
regenerated from scratch.

When large errors are found in other libm implementations, reporting bugs 
against those implementations is a good idea as well (depending on each 
implementation's policy on what errors count as bugs).

-- 
Joseph S. Myers
joseph@codesourcery.com