From 9453a7481a3e41ed18ce08d83bba9ad895c16d4d Mon Sep 17 00:00:00 2001 From: Paul Zimmermann Date: Thu, 13 Feb 2020 17:02:37 +0100 Subject: improve corner cases in asin, added test for acos --- src/asin.c | 6 ++++-- tests/tacos.c | 35 +++++++++++++++++++++++++++++++++++ 2 files changed, 39 insertions(+), 2 deletions(-) diff --git a/src/asin.c b/src/asin.c index 0d975e4..96749db 100644 --- a/src/asin.c +++ b/src/asin.c @@ -87,7 +87,7 @@ mpc_asin (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) mpfr_rnd_t rnd_re, rnd_im; mpc_t z1; int inex, inex_re, inex_im, loop = 0; - mpfr_exp_t saved_emin, saved_emax; + mpfr_exp_t saved_emin, saved_emax, err; /* special values */ if (mpfr_nan_p (mpc_realref (op)) || mpfr_nan_p (mpc_imagref (op))) @@ -206,11 +206,13 @@ mpc_asin (mpc_ptr rop, mpc_srcptr op, mpc_rnd_t rnd) rnd_im = MPC_RND_IM(rnd); p = p_re >= p_im ? p_re : p_im; mpc_init2 (z1, p); + err = 0; /* number of lost bits */ while (1) { - mpfr_exp_t ex, ey, err; + mpfr_exp_t ex, ey; loop ++; + p += err; /* add number of lost bits in previous loop */ p += (loop <= 2) ? mpc_ceil_log2 (p) + 3 : p / 2; mpfr_set_prec (mpc_realref(z1), p); mpfr_set_prec (mpc_imagref(z1), p); diff --git a/tests/tacos.c b/tests/tacos.c index 6ad1cc6..2f764b1 100644 --- a/tests/tacos.c +++ b/tests/tacos.c @@ -28,11 +28,46 @@ along with this program. If not, see http://www.gnu.org/licenses/ . #include "data_check.tpl" #include "tgeneric.tpl" +/* test with reduced exponent range */ +static void +test_reduced (void) +{ + mpfr_exp_t emin = mpfr_get_emin (); + mpfr_exp_t emax = mpfr_get_emax (); + mpc_t x, z, zr; + + mpfr_set_emin (-148); + mpfr_set_emax (128); + mpc_init2 (x, 24); + mpc_init2 (z, 24); + mpc_init2 (zr, 24); + mpfr_set_flt (mpc_realref (x), -0x0.01f28c10p0f, MPFR_RNDN); + mpfr_set_flt (mpc_imagref (x), -0x6.79cdd0p-68f, MPFR_RNDN); + mpc_acos (z, x, MPC_RNDNN); + /* square root is 0.00038652126908433 - 1.37634353022868e12*I */ + mpfr_set_flt (mpc_realref (zr), 0x1.941242p0f, MPFR_RNDN); + mpfr_set_flt (mpc_imagref (zr), 0x6.79da18p-68f, MPFR_RNDN); + if (mpc_cmp (z, zr)) + { + printf ("Incorrect acos with reduced exponent range:\n"); + mpfr_printf ("Expected (%Re,%Re)\n", mpc_realref (zr), mpc_imagref (zr)); + mpfr_printf ("Got (%Re,%Re)\n", mpc_realref (z), mpc_imagref (z)); + exit (1); + } + mpc_clear (x); + mpc_clear (z); + mpc_clear (zr); + mpfr_set_emin (emin); + mpfr_set_emax (emax); +} + int main (void) { test_start (); + test_reduced (); + data_check_template ("acos.dsc", "acos.dat"); tgeneric_template ("acos.dsc", 2, 512, 7, 7); -- cgit v1.2.1