Bug 281001 - Improve accuracy of asinf(x) and acosf(x)
Summary: Improve accuracy of asinf(x) and acosf(x)
Status: Closed FIXED
Alias: None
Product: Base System
Classification: Unclassified
Component: bin (show other bugs)
Version: CURRENT
Hardware: Any Any
: --- Affects Only Me
Assignee: freebsd-numerics (Nobody)
URL:
Keywords:
Depends on:
Blocks:
 
Reported: 2024-08-22 22:23 UTC by Steve Kargl
Modified: 2024-09-08 08:03 UTC (History)
2 users (show)

See Also:

Attachments
patch to lib/msun/src/e_asinf.c (1.21 KB, patch)
2024-08-22 22:23 UTC, Steve Kargl 		no flags Details | Diff
new patch with changes for both asinf() and acosf(). (2.75 KB, patch)
2024-08-23 18:06 UTC, Steve Kargl 		no flags Details | Diff
Show Obsolete (1) View All Add an attachment (proposed patch, testcase, etc.)

Note You need to log in before you can comment on or make changes to this bug.
Description Steve Kargl freebsd_committer freebsd_triage 2024-08-22 22:23:07 UTC 
Created attachment 253025 [details]
patch to lib/msun/src/e_asinf.c

The attached patch improves the accuracy of asinf(x).  It contains a better rational
approximation.  For exhaustive testing in the interval, the current libm gives

% ./tlibm asin -fPED -x 0x1p-12f -X 1
Interval tested for asinf: [0.000244141,1]
       ulp <= 0.5:  97.916% 98564994 |  97.916% 98564994
0.5 <  ulp <  0.6:  2.038% 2051023 |  99.953% 100616017
0.6 <  ulp <  0.7:  0.047%   47254 | 100.000% 100663271
0.7 <  ulp <  0.8:  0.000%      25 | 100.000% 100663296
Max ulp: 0.729891 at 5.00732839e-01

which isn't too bad given that much of the computation is actually done in double
floating point.

With the new rational approximation, exhaustive testing yields

% ./tlibm asin -fPED -x 0x1p-12f -X 1
Interval tested for asinf: [0.000244141,1]
       ulp <= 0.5:  99.711% 100372643 |  99.711% 100372643
0.5 <  ulp <  0.6:  0.288%  290357 | 100.000% 100663000
0.6 <  ulp <  0.7:  0.000%     296 | 100.000% 100663296
Max ulp: 0.636344 at 5.09706438e-01

That is, it's better.
Comment 1 Steve Kargl freebsd_committer freebsd_triage 2024-08-23 01:26:31 UTC 
Just realized that same rational approximation in e_asinf(x) appears in e_cosf(x).  I'll test that one tomorrow.
Comment 2 Steve Kargl freebsd_committer freebsd_triage 2024-08-23 18:06:58 UTC 
Created attachment 253043 [details]
new patch with changes for both asinf() and acosf().

The new patch has updated the rational approximation in acosf() to match the new rational approximation in asinf().  This change provides only a small improve in accuracy in exhaustive testing over the relevant intervals.

Unpatched acosf() gives

% ./tlibm acos -fPED -x -1 -X -0x1p-12f
Interval tested for acosf: [-1,-0.000244141]
       ulp <= 0.5:  97.008% 97651921 |  97.008% 97651921
0.5 <  ulp <  0.6:   2.441%  2457242 |  99.450% 100109163
0.6 <  ulp <  0.7:   0.472%   475503 |  99.922% 100584666
0.7 <  ulp <  0.8:   0.071%    71309 |  99.993% 100655975
0.8 <  ulp <  0.9:   0.007%     7319 | 100.000% 100663294
0.9 <  ulp <  1.0:   0.000%        2 | 100.000% 100663296
Max ulp: 0.914007 at -5.01484931e-01

% ./tlibm acos -fPED -x 0x1p-12f -X 1
Interval tested for acosf: [0.000244141,1]
       ulp <= 0.5:  97.317% 97962530 |  97.317% 97962530
0.5 <  ulp <  0.6:   2.340%  2355182 |  99.657% 100317712
0.6 <  ulp <  0.7:   0.314%   316134 |  99.971% 100633846
0.7 <  ulp <  0.8:   0.029%    29450 | 100.000% 100663296
Max ulp: 0.796035 at 4.99814630e-01


Patched acosf() gives

% ./tlibm acos -fPED -x -1 -X -0x1p-12f
Interval tested for acosf: [-1,-0.000244141]
       ulp <= 0.5:  97.010% 97653245 |  97.010% 97653245
0.5 <  ulp <  0.6:   2.442%  2458373 |  99.452% 100111618
0.6 <  ulp <  0.7:   0.473%   476012 |  99.925% 100587630
0.7 <  ulp <  0.8:   0.068%    68603 |  99.993% 100656233
0.8 <  ulp <  0.9:   0.007%     7063 | 100.000% 100663296
Max ulp: 0.896189 at -5.04511118e-01

% ./tlibm acos -fPED -x 0x1p-12f -X 1
Interval tested for acosf: [0.000244141,1]
       ulp <= 0.5:  97.650% 98298175 |  97.650% 98298175
0.5 <  ulp <  0.6:   2.028%  2041709 |  99.679% 100339884
0.6 <  ulp <  0.7:   0.292%   293555 |  99.970% 100633439
0.7 <  ulp <  0.8:   0.030%    29857 | 100.000% 100663296
Max ulp: 0.775875 at 4.91849005e-01
Comment 3 commit-hook freebsd_committer freebsd_triage 2024-08-29 19:46:45 UTC 
A commit in branch main references this bug:

URL: https://cgit.FreeBSD.org/src/commit/?id=41e016289f77deb88b0ef1ec3f7b2ab3515ac7c8

commit 41e016289f77deb88b0ef1ec3f7b2ab3515ac7c8
Author:     Steve Kargl <kargl@FreeBSD.org>
AuthorDate: 2024-08-29 19:44:48 +0000
Commit:     Dimitry Andric <dim@FreeBSD.org>
CommitDate: 2024-08-29 19:44:48 +0000

    Improve accuracy of asinf(3) and acosf(3)

    This uses a better rational approximation to improve the accuracy of
    both functions. For exhaustive testing of asinf(3) in the interval, the
    current libm gives:

        % ./tlibm asin -fPED -x 0x1p-12f -X 1
        Interval tested for asinf: [0.000244141,1]
               ulp <= 0.5:  97.916% 98564994 |  97.916% 98564994
        0.5 <  ulp <  0.6:  2.038% 2051023 |  99.953% 100616017
        0.6 <  ulp <  0.7:  0.047%   47254 | 100.000% 100663271
        0.7 <  ulp <  0.8:  0.000%      25 | 100.000% 100663296
        Max ulp: 0.729891 at 5.00732839e-01

    which isn't too bad given that much of the computation is actually done
    in double floating point.

    With the new rational approximation, exhaustive testing yields:

        % ./tlibm asin -fPED -x 0x1p-12f -X 1
        Interval tested for asinf: [0.000244141,1]
               ulp <= 0.5:  99.711% 100372643 |  99.711% 100372643
        0.5 <  ulp <  0.6:  0.288%  290357 | 100.000% 100663000
        0.6 <  ulp <  0.7:  0.000%     296 | 100.000% 100663296
        Max ulp: 0.636344 at 5.09706438e-01

    Similarly, for exhaustive testing of asinf(3) in the interval, the
    current libm gives:

        % ./tlibm acos -fPED -x -1 -X -0x1p-12f
        Interval tested for acosf: [-1,-0.000244141]
               ulp <= 0.5:  97.008% 97651921 |  97.008% 97651921
        0.5 <  ulp <  0.6:   2.441%  2457242 |  99.450% 100109163
        0.6 <  ulp <  0.7:   0.472%   475503 |  99.922% 100584666
        0.7 <  ulp <  0.8:   0.071%    71309 |  99.993% 100655975
        0.8 <  ulp <  0.9:   0.007%     7319 | 100.000% 100663294
        0.9 <  ulp <  1.0:   0.000%        2 | 100.000% 100663296
        Max ulp: 0.914007 at -5.01484931e-01

        % ./tlibm acos -fPED -x 0x1p-12f -X 1
        Interval tested for acosf: [0.000244141,1]
               ulp <= 0.5:  97.317% 97962530 |  97.317% 97962530
        0.5 <  ulp <  0.6:   2.340%  2355182 |  99.657% 100317712
        0.6 <  ulp <  0.7:   0.314%   316134 |  99.971% 100633846
        0.7 <  ulp <  0.8:   0.029%    29450 | 100.000% 100663296
        Max ulp: 0.796035 at 4.99814630e-01

    With the new rational approximation, exhaustive testing yields:

        % ./tlibm acos -fPED -x -1 -X -0x1p-12f
        Interval tested for acosf: [-1,-0.000244141]
               ulp <= 0.5:  97.010% 97653245 |  97.010% 97653245
        0.5 <  ulp <  0.6:   2.442%  2458373 |  99.452% 100111618
        0.6 <  ulp <  0.7:   0.473%   476012 |  99.925% 100587630
        0.7 <  ulp <  0.8:   0.068%    68603 |  99.993% 100656233
        0.8 <  ulp <  0.9:   0.007%     7063 | 100.000% 100663296
        Max ulp: 0.896189 at -5.04511118e-01

        % ./tlibm acos -fPED -x 0x1p-12f -X 1
        Interval tested for acosf: [0.000244141,1]
               ulp <= 0.5:  97.650% 98298175 |  97.650% 98298175
        0.5 <  ulp <  0.6:   2.028%  2041709 |  99.679% 100339884
        0.6 <  ulp <  0.7:   0.292%   293555 |  99.970% 100633439
        0.7 <  ulp <  0.8:   0.030%    29857 | 100.000% 100663296
        Max ulp: 0.775875 at 4.91849005e-01

    PR:             281001
    MFC after:      1 week

 lib/msun/src/e_acosf.c | 20 +++++++++++++-------
 lib/msun/src/e_asinf.c | 22 ++++++++++++++--------
 2 files changed, 27 insertions(+), 15 deletions(-)
Comment 4 Dimitry Andric freebsd_committer freebsd_triage 2024-08-29 19:50:53 UTC 
Btw I fixed one typo in the patch: approximatione -> approximation. I assumed you meant the singular.
Comment 5 commit-hook freebsd_committer freebsd_triage 2024-09-08 07:55:14 UTC 
A commit in branch stable/14 references this bug:

URL: https://cgit.FreeBSD.org/src/commit/?id=33c82f11c2674beb50304b14248ff46372def520

commit 33c82f11c2674beb50304b14248ff46372def520
Author:     Steve Kargl <kargl@FreeBSD.org>
AuthorDate: 2024-08-29 19:44:48 +0000
Commit:     Dimitry Andric <dim@FreeBSD.org>
CommitDate: 2024-09-08 07:37:52 +0000

    Improve accuracy of asinf(3) and acosf(3)

    This uses a better rational approximation to improve the accuracy of
    both functions. For exhaustive testing of asinf(3) in the interval, the
    current libm gives:

        % ./tlibm asin -fPED -x 0x1p-12f -X 1
        Interval tested for asinf: [0.000244141,1]
               ulp <= 0.5:  97.916% 98564994 |  97.916% 98564994
        0.5 <  ulp <  0.6:  2.038% 2051023 |  99.953% 100616017
        0.6 <  ulp <  0.7:  0.047%   47254 | 100.000% 100663271
        0.7 <  ulp <  0.8:  0.000%      25 | 100.000% 100663296
        Max ulp: 0.729891 at 5.00732839e-01

    which isn't too bad given that much of the computation is actually done
    in double floating point.

    With the new rational approximation, exhaustive testing yields:

        % ./tlibm asin -fPED -x 0x1p-12f -X 1
        Interval tested for asinf: [0.000244141,1]
               ulp <= 0.5:  99.711% 100372643 |  99.711% 100372643
        0.5 <  ulp <  0.6:  0.288%  290357 | 100.000% 100663000
        0.6 <  ulp <  0.7:  0.000%     296 | 100.000% 100663296
        Max ulp: 0.636344 at 5.09706438e-01

    Similarly, for exhaustive testing of asinf(3) in the interval, the
    current libm gives:

        % ./tlibm acos -fPED -x -1 -X -0x1p-12f
        Interval tested for acosf: [-1,-0.000244141]
               ulp <= 0.5:  97.008% 97651921 |  97.008% 97651921
        0.5 <  ulp <  0.6:   2.441%  2457242 |  99.450% 100109163
        0.6 <  ulp <  0.7:   0.472%   475503 |  99.922% 100584666
        0.7 <  ulp <  0.8:   0.071%    71309 |  99.993% 100655975
        0.8 <  ulp <  0.9:   0.007%     7319 | 100.000% 100663294
        0.9 <  ulp <  1.0:   0.000%        2 | 100.000% 100663296
        Max ulp: 0.914007 at -5.01484931e-01

        % ./tlibm acos -fPED -x 0x1p-12f -X 1
        Interval tested for acosf: [0.000244141,1]
               ulp <= 0.5:  97.317% 97962530 |  97.317% 97962530
        0.5 <  ulp <  0.6:   2.340%  2355182 |  99.657% 100317712
        0.6 <  ulp <  0.7:   0.314%   316134 |  99.971% 100633846
        0.7 <  ulp <  0.8:   0.029%    29450 | 100.000% 100663296
        Max ulp: 0.796035 at 4.99814630e-01

    With the new rational approximation, exhaustive testing yields:

        % ./tlibm acos -fPED -x -1 -X -0x1p-12f
        Interval tested for acosf: [-1,-0.000244141]
               ulp <= 0.5:  97.010% 97653245 |  97.010% 97653245
        0.5 <  ulp <  0.6:   2.442%  2458373 |  99.452% 100111618
        0.6 <  ulp <  0.7:   0.473%   476012 |  99.925% 100587630
        0.7 <  ulp <  0.8:   0.068%    68603 |  99.993% 100656233
        0.8 <  ulp <  0.9:   0.007%     7063 | 100.000% 100663296
        Max ulp: 0.896189 at -5.04511118e-01

        % ./tlibm acos -fPED -x 0x1p-12f -X 1
        Interval tested for acosf: [0.000244141,1]
               ulp <= 0.5:  97.650% 98298175 |  97.650% 98298175
        0.5 <  ulp <  0.6:   2.028%  2041709 |  99.679% 100339884
        0.6 <  ulp <  0.7:   0.292%   293555 |  99.970% 100633439
        0.7 <  ulp <  0.8:   0.030%    29857 | 100.000% 100663296
        Max ulp: 0.775875 at 4.91849005e-01

    PR:             281001
    MFC after:      1 week

    (cherry picked from commit 41e016289f77deb88b0ef1ec3f7b2ab3515ac7c8)

 lib/msun/src/e_acosf.c | 20 +++++++++++++-------
 lib/msun/src/e_asinf.c | 22 ++++++++++++++--------
 2 files changed, 27 insertions(+), 15 deletions(-)
Comment 6 commit-hook freebsd_committer freebsd_triage 2024-09-08 07:56:18 UTC 
A commit in branch stable/13 references this bug:

URL: https://cgit.FreeBSD.org/src/commit/?id=75bca9e7553c70d97b36df30ab05bbbbacb8a499

commit 75bca9e7553c70d97b36df30ab05bbbbacb8a499
Author:     Steve Kargl <kargl@FreeBSD.org>
AuthorDate: 2024-08-29 19:44:48 +0000
Commit:     Dimitry Andric <dim@FreeBSD.org>
CommitDate: 2024-09-08 07:37:57 +0000

    Improve accuracy of asinf(3) and acosf(3)

    This uses a better rational approximation to improve the accuracy of
    both functions. For exhaustive testing of asinf(3) in the interval, the
    current libm gives:

        % ./tlibm asin -fPED -x 0x1p-12f -X 1
        Interval tested for asinf: [0.000244141,1]
               ulp <= 0.5:  97.916% 98564994 |  97.916% 98564994
        0.5 <  ulp <  0.6:  2.038% 2051023 |  99.953% 100616017
        0.6 <  ulp <  0.7:  0.047%   47254 | 100.000% 100663271
        0.7 <  ulp <  0.8:  0.000%      25 | 100.000% 100663296
        Max ulp: 0.729891 at 5.00732839e-01

    which isn't too bad given that much of the computation is actually done
    in double floating point.

    With the new rational approximation, exhaustive testing yields:

        % ./tlibm asin -fPED -x 0x1p-12f -X 1
        Interval tested for asinf: [0.000244141,1]
               ulp <= 0.5:  99.711% 100372643 |  99.711% 100372643
        0.5 <  ulp <  0.6:  0.288%  290357 | 100.000% 100663000
        0.6 <  ulp <  0.7:  0.000%     296 | 100.000% 100663296
        Max ulp: 0.636344 at 5.09706438e-01

    Similarly, for exhaustive testing of asinf(3) in the interval, the
    current libm gives:

        % ./tlibm acos -fPED -x -1 -X -0x1p-12f
        Interval tested for acosf: [-1,-0.000244141]
               ulp <= 0.5:  97.008% 97651921 |  97.008% 97651921
        0.5 <  ulp <  0.6:   2.441%  2457242 |  99.450% 100109163
        0.6 <  ulp <  0.7:   0.472%   475503 |  99.922% 100584666
        0.7 <  ulp <  0.8:   0.071%    71309 |  99.993% 100655975
        0.8 <  ulp <  0.9:   0.007%     7319 | 100.000% 100663294
        0.9 <  ulp <  1.0:   0.000%        2 | 100.000% 100663296
        Max ulp: 0.914007 at -5.01484931e-01

        % ./tlibm acos -fPED -x 0x1p-12f -X 1
        Interval tested for acosf: [0.000244141,1]
               ulp <= 0.5:  97.317% 97962530 |  97.317% 97962530
        0.5 <  ulp <  0.6:   2.340%  2355182 |  99.657% 100317712
        0.6 <  ulp <  0.7:   0.314%   316134 |  99.971% 100633846
        0.7 <  ulp <  0.8:   0.029%    29450 | 100.000% 100663296
        Max ulp: 0.796035 at 4.99814630e-01

    With the new rational approximation, exhaustive testing yields:

        % ./tlibm acos -fPED -x -1 -X -0x1p-12f
        Interval tested for acosf: [-1,-0.000244141]
               ulp <= 0.5:  97.010% 97653245 |  97.010% 97653245
        0.5 <  ulp <  0.6:   2.442%  2458373 |  99.452% 100111618
        0.6 <  ulp <  0.7:   0.473%   476012 |  99.925% 100587630
        0.7 <  ulp <  0.8:   0.068%    68603 |  99.993% 100656233
        0.8 <  ulp <  0.9:   0.007%     7063 | 100.000% 100663296
        Max ulp: 0.896189 at -5.04511118e-01

        % ./tlibm acos -fPED -x 0x1p-12f -X 1
        Interval tested for acosf: [0.000244141,1]
               ulp <= 0.5:  97.650% 98298175 |  97.650% 98298175
        0.5 <  ulp <  0.6:   2.028%  2041709 |  99.679% 100339884
        0.6 <  ulp <  0.7:   0.292%   293555 |  99.970% 100633439
        0.7 <  ulp <  0.8:   0.030%    29857 | 100.000% 100663296
        Max ulp: 0.775875 at 4.91849005e-01

    PR:             281001
    MFC after:      1 week

    (cherry picked from commit 41e016289f77deb88b0ef1ec3f7b2ab3515ac7c8)

 lib/msun/src/e_acosf.c | 20 +++++++++++++-------
 lib/msun/src/e_asinf.c | 22 ++++++++++++++--------
 2 files changed, 27 insertions(+), 15 deletions(-)