Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[libc][math] Implement double precision expm1 function correctly rounded for all rounding modes. #67048

Merged
merged 3 commits into from
Sep 28, 2023

Conversation

lntue
Copy link
Contributor

@lntue lntue commented Sep 21, 2023

Implementing expm1 function for double precision based on exp function algorithm:

  • Reduced x = log2(e) * (hi + mid1 + mid2) + lo, where:
    • hi is an integer
    • mid1 * 2^-6 is an integer
    • mid2 * 2^-12 is an integer
    • |lo| < 2^-13 + 2^-30
  • Then exp(x) - 1 = 2^hi * 2^mid1 * 2^mid2 * exp(lo) - 1 ~ 2^hi * (2^mid1 * 2^mid2 * (1 + lo * P(lo)) - 2^(-hi) )
  • We evaluate fast pass with P(lo) is a degree-3 Taylor polynomial of (e^lo - 1) / lo in double precision
  • If the Ziv accuracy test fails, we use degree-6 Taylor polynomial of (e^lo - 1) / lo in double double precision
  • If the Ziv accuracy test still fails, we re-evaluate everything in 128-bit precision.

@llvmbot llvmbot added the libc label Sep 21, 2023
@lntue
Copy link
Contributor Author

lntue commented Sep 21, 2023

@zimmermann6

Copy link
Contributor

@michaelrj-google michaelrj-google left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Overall LGTM from the libc side, with minor nits.


#include <errno.h>

// #define DEBUGDEBUG
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Did you mean to leave the debug code in here?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I found myself retyping these debug codes quite many times, and I couldn't find a better place to share this type of debugging codes yet. So I think I should leave it at least in one function so that a simple copy paste during development of other functions would be enough.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

SGTM

0x3FD79289C6E6A5C0,
0x3FD05DE80A173EA0, // 0x1.05de80a173eap-2
0xbf1eb7a4cb841fcc, // -0x1.eb7a4cb841fccp-14
0xbf19a61fb925970d,
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

should this one have a comment on its value as well?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fixed.

@github-actions
Copy link

github-actions bot commented Sep 27, 2023

✅ With the latest revision this PR passed the C/C++ code formatter.

@zimmermann6
Copy link

looks good to me. All core-math tests do pass for all 4 rounding modes, and efficiency is quite good:

zimmerma@biscotte:~/svn/core-math$ LIBM=/localdisk/zimmerma/llvm-project/build/projects/libc/lib/libllvmlibc.a CORE_MATH_PERF_MODE=rdtsc ./perf.sh expm1
GNU libc version: 2.37
GNU libc release: stable
[####################] 100 %
Ntrial = 20 ; Min = 44.681 + 0.398 clc/call; Median-Min = 0.318 clc/call; Max = 45.533 clc/call;
[####################] 100 %
Ntrial = 20 ; Min = 44.593 + 0.350 clc/call; Median-Min = 0.287 clc/call; Max = 46.071 clc/call;
[####################] 100 %
Ntrial = 20 ; Min = 29.969 + 0.284 clc/call; Median-Min = 0.317 clc/call; Max = 30.416 clc/call;
zimmerma@biscotte:~/svn/core-math$ LIBM=/localdisk/zimmerma/llvm-project/build/projects/libc/lib/libllvmlibc.a CORE_MATH_PERF_MODE=rdtsc PERF_ARGS=--latency ./perf.sh expm1
GNU libc version: 2.37
GNU libc release: stable
[####################] 100 %
Ntrial = 20 ; Min = 73.321 + 0.314 clc/call; Median-Min = 0.302 clc/call; Max = 73.831 clc/call;
[####################] 100 %
Ntrial = 20 ; Min = 120.724 + 0.314 clc/call; Median-Min = 0.312 clc/call; Max = 121.278 clc/call;
[####################] 100 %
Ntrial = 20 ; Min = 70.920 + 0.449 clc/call; Median-Min = 0.418 clc/call; Max = 71.737 clc/call;

(first timings are for core-math, last for llvm)

@lntue lntue merged commit da28593 into llvm:main Sep 28, 2023
3 checks passed
legrosbuffle pushed a commit to legrosbuffle/llvm-project that referenced this pull request Sep 29, 2023
…ded for all rounding modes. (llvm#67048)

Implementing expm1 function for double precision based on exp function
algorithm:

- Reduced x = log2(e) * (hi + mid1 + mid2) + lo, where:
  * hi is an integer
  * mid1 * 2^-6 is an integer
  * mid2 * 2^-12 is an integer
  * |lo| < 2^-13 + 2^-30
- Then exp(x) - 1 = 2^hi * 2^mid1 * 2^mid2 * exp(lo) - 1 ~ 2^hi *
(2^mid1 * 2^mid2 * (1 + lo * P(lo)) - 2^(-hi) )
- We evaluate fast pass with P(lo) is a degree-3 Taylor polynomial of
(e^lo - 1) / lo in double precision
- If the Ziv accuracy test fails, we use degree-6 Taylor polynomial of
(e^lo - 1) / lo in double double precision
- If the Ziv accuracy test still fails, we re-evaluate everything in
128-bit precision.
@lntue lntue deleted the expm1 branch October 19, 2023 14:31
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

4 participants