public inbox for glibc-cvs@sourceware.org help / color / mirror / Atom feed
From: Adhemerval Zanella <azanella@sourceware.org> To: glibc-cvs@sourceware.org Subject: [glibc] math: Fix inaccuracy of j0f for x >= 2^127 when sin(x)+cos(x) is tiny Date: Fri, 7 Aug 2020 19:33:34 +0000 (GMT) [thread overview] Message-ID: <20200807193334.26B703857C55@sourceware.org> (raw) https://sourceware.org/git/gitweb.cgi?p=glibc.git;h=b7dd366dbeeff3016d0554c08a5f5c8cfe4c31d1 commit b7dd366dbeeff3016d0554c08a5f5c8cfe4c31d1 Author: Paul Zimmermann <Paul.Zimmermann@inria.fr> Date: Fri Aug 7 16:14:53 2020 -0300 math: Fix inaccuracy of j0f for x >= 2^127 when sin(x)+cos(x) is tiny Checked on x86_64-linux-gnu and i686-linux-gnu. Diff: --- math/auto-libm-test-in | 2 ++ sysdeps/ieee754/flt-32/e_j0f.c | 17 ++++++++++++++++- 2 files changed, 18 insertions(+), 1 deletion(-) diff --git a/math/auto-libm-test-in b/math/auto-libm-test-in index 4414e54d93..5d488a8711 100644 --- a/math/auto-libm-test-in +++ b/math/auto-libm-test-in @@ -5748,6 +5748,8 @@ j0 0x1p16382 j0 0x1p16383 # the next value generates larger error bounds on x86_64 (binary32) j0 0x2.602774p+0 xfail-rounding:ibm128-libgcc +# the next value exercises the flt-32 code path for x >= 2^127 +j0 0x8.2f4ecp+124 j1 -1.0 j1 0.0 diff --git a/sysdeps/ieee754/flt-32/e_j0f.c b/sysdeps/ieee754/flt-32/e_j0f.c index c89b9f2688..5d29611eb7 100644 --- a/sysdeps/ieee754/flt-32/e_j0f.c +++ b/sysdeps/ieee754/flt-32/e_j0f.c @@ -55,7 +55,22 @@ __ieee754_j0f(float x) z = -__cosf(x+x); if ((s*c)<zero) cc = z/ss; else ss = z/cc; - } + } else { + /* We subtract (exactly) a value x0 such that + cos(x0)+sin(x0) is very near to 0, and use the identity + sin(x-x0) = sin(x)*cos(x0)-cos(x)*sin(x0) to get + sin(x) + cos(x) with extra accuracy. */ + float x0 = 0xe.d4108p+124f; + float y = x - x0; /* exact */ + /* sin(y) = sin(x)*cos(x0)-cos(x)*sin(x0) */ + z = __sinf (y); + float eps = 0x1.5f263ep-24f; + /* cos(x0) ~ -sin(x0) + eps */ + z += eps * __cosf (x); + /* now z ~ (sin(x)-cos(x))*cos(x0) */ + float cosx0 = -0xb.504f3p-4f; + cc = z / cosx0; + } /* * j0(x) = 1/sqrt(pi) * (P(0,x)*cc - Q(0,x)*ss) / sqrt(x) * y0(x) = 1/sqrt(pi) * (P(0,x)*ss + Q(0,x)*cc) / sqrt(x)
reply other threads:[~2020-08-07 19:33 UTC|newest] Thread overview: [no followups] expand[flat|nested] mbox.gz Atom feed
Reply instructions: You may reply publicly to this message via plain-text email using any one of the following methods: * Save the following mbox file, import it into your mail client, and reply-to-all from there: mbox Avoid top-posting and favor interleaved quoting: https://en.wikipedia.org/wiki/Posting_style#Interleaved_style * Reply using the --to, --cc, and --in-reply-to switches of git-send-email(1): git send-email \ --in-reply-to=20200807193334.26B703857C55@sourceware.org \ --to=azanella@sourceware.org \ --cc=glibc-cvs@sourceware.org \ /path/to/YOUR_REPLY https://kernel.org/pub/software/scm/git/docs/git-send-email.html * If your mail client supports setting the In-Reply-To header via mailto: links, try the mailto: linkBe sure your reply has a Subject: header at the top and a blank line before the message body.
This is a public inbox, see mirroring instructions for how to clone and mirror all data and code used for this inbox; as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).