diff options
author | pelissip <pelissip@280ebfd0-de03-0410-8827-d642c229c3f4> | 2005-05-18 08:22:55 +0000 |
---|---|---|
committer | pelissip <pelissip@280ebfd0-de03-0410-8827-d642c229c3f4> | 2005-05-18 08:22:55 +0000 |
commit | 6ce84039a34830e8841936c370539b5e01386e79 (patch) | |
tree | 73ef169251ea14a30662f3ed220be13116dacf45 /cosh.c | |
parent | a2ac43cf79fb45c4dada0fc9880a184ac32b3878 (diff) | |
download | mpfr-6ce84039a34830e8841936c370539b5e01386e79.tar.gz |
Fix bug of overflow.
Fix bug of ternary value in case of overflow.
git-svn-id: svn://scm.gforge.inria.fr/svn/mpfr/trunk@3585 280ebfd0-de03-0410-8827-d642c229c3f4
Diffstat (limited to 'cosh.c')
-rw-r--r-- | cosh.c | 28 |
1 files changed, 17 insertions, 11 deletions
@@ -57,15 +57,12 @@ mpfr_cosh (mpfr_ptr y, mpfr_srcptr xt , mp_rnd_t rnd_mode) MPFR_SAVE_EXPO_MARK (expo); MPFR_TMP_INIT_ABS(x, xt); - /* General case */ { /* Declaration of the intermediary variable */ mpfr_t t, te; - /* Declaration of the size variable */ mp_prec_t Ny = MPFR_PREC(y); /* Precision of output variable */ - mp_prec_t Nt; /* Precision of the intermediary variable */ long int err; /* Precision of error */ MPFR_ZIV_DECL (loop); @@ -73,7 +70,7 @@ mpfr_cosh (mpfr_ptr y, mpfr_srcptr xt , mp_rnd_t rnd_mode) /* compute the precision of intermediary variable */ /* The optimal number of bits : see algorithms.tex */ Nt = Ny + 3 + MPFR_INT_CEIL_LOG2 (Ny); - + /* initialise of intermediary variables */ mpfr_init2 (t, Nt); mpfr_init2 (te, Nt); @@ -83,18 +80,28 @@ mpfr_cosh (mpfr_ptr y, mpfr_srcptr xt , mp_rnd_t rnd_mode) for (;;) { /* Compute cosh */ + mpfr_clear_flags (); mpfr_exp (te, x, GMP_RNDD); /* exp(x) */ + /* exp can overflow (but not underflow since x>0) */ + /* BUG/TODO/FIXME: exp can overflow but cosh may be representable! */ + if (MPFR_UNLIKELY (mpfr_overflow_p ())) + { + inexact = mpfr_overflow (y, rnd_mode, MPFR_SIGN_POS); + MPFR_SAVE_EXPO_UPDATE_FLAGS (expo, MPFR_FLAGS_OVERFLOW); + break; + } mpfr_ui_div (t, 1, te, GMP_RNDU); /* 1/exp(x) */ mpfr_add (t, te, t, GMP_RNDU); /* exp(x) + 1/exp(x)*/ mpfr_div_2ui (t, t, 1, GMP_RNDN); /* 1/2(exp(x) + 1/exp(x))*/ - - /* Estimation of the error */ + + /* Estimation of the error */ err = Nt - 3; - /* Check if we can round */ - if (MPFR_UNLIKELY (MPFR_IS_INF (t)) - || MPFR_LIKELY (MPFR_CAN_ROUND (t, err, Ny, rnd_mode))) - break; + if (MPFR_LIKELY (MPFR_CAN_ROUND (t, err, Ny, rnd_mode))) + { + inexact = mpfr_set (y, t, rnd_mode); + break; + } /* Actualisation of the precision */ MPFR_ZIV_NEXT (loop, Nt); @@ -102,7 +109,6 @@ mpfr_cosh (mpfr_ptr y, mpfr_srcptr xt , mp_rnd_t rnd_mode) mpfr_set_prec (te, Nt); } MPFR_ZIV_FREE (loop); - inexact = mpfr_set (y, t, rnd_mode); mpfr_clear (te); mpfr_clear (t); |