summaryrefslogtreecommitdiff
path: root/Build/source/libs/mpfr/mpfr-src/src/exp3.c
diff options
context:
space:
mode:
Diffstat (limited to 'Build/source/libs/mpfr/mpfr-src/src/exp3.c')
-rw-r--r--Build/source/libs/mpfr/mpfr-src/src/exp3.c28
1 files changed, 16 insertions, 12 deletions
diff --git a/Build/source/libs/mpfr/mpfr-src/src/exp3.c b/Build/source/libs/mpfr/mpfr-src/src/exp3.c
index 04f3e68519d..2940b471b49 100644
--- a/Build/source/libs/mpfr/mpfr-src/src/exp3.c
+++ b/Build/source/libs/mpfr/mpfr-src/src/exp3.c
@@ -42,7 +42,7 @@ static void
mpfr_exp_rational (mpfr_ptr y, mpz_ptr p, long r, int m,
mpz_t *Q, mpfr_prec_t *mult)
{
- unsigned long n, i, j;
+ mp_bitcnt_t n, h, i, j; /* unsigned type, which is >= unsigned long */
mpz_t *S, *ptoj;
mpfr_prec_t *log2_nb_terms;
mpfr_exp_t diff, expo;
@@ -58,8 +58,9 @@ mpfr_exp_rational (mpfr_ptr y, mpz_ptr p, long r, int m,
/* Normalize p */
MPFR_ASSERTD (mpz_cmp_ui (p, 0) != 0);
n = mpz_scan1 (p, 0); /* number of trailing zeros in p */
+ MPFR_ASSERTN (n <= LONG_MAX); /* This is a limitation. */
mpz_tdiv_q_2exp (p, p, n);
- r -= n; /* since |p/2^r| < 1 and p >= 1, r >= 1 */
+ r -= (long) n; /* since |p/2^r| < 1 and p >= 1, r >= 1 */
/* Set initial var */
mpz_set (ptoj[0], p);
@@ -75,6 +76,7 @@ mpfr_exp_rational (mpfr_ptr y, mpz_ptr p, long r, int m,
/* Main Loop */
n = 1UL << m;
+ MPFR_ASSERTN (n != 0); /* no overflow */
for (i = 1; (prec_i_have < precy) && (i < n); i++)
{
/* invariant: Q[0]*Q[1]*...*Q[k] equals i! */
@@ -112,14 +114,14 @@ mpfr_exp_rational (mpfr_ptr y, mpz_ptr p, long r, int m,
/* accumulate all products in S[0] and Q[0]. Warning: contrary to above,
here we do not have log2_nb_terms[k-1] = log2_nb_terms[k]+1. */
- l = 0; /* number of accumulated terms in the right part S[k]/Q[k] */
+ h = 0; /* number of accumulated terms in the right part S[k]/Q[k] */
while (k > 0)
{
j = log2_nb_terms[k-1];
mpz_mul (S[k], S[k], ptoj[j]);
mpz_mul (S[k-1], S[k-1], Q[k]);
- l += 1 << log2_nb_terms[k];
- mpz_mul_2exp (S[k-1], S[k-1], r * l);
+ h += (mp_bitcnt_t) 1 << log2_nb_terms[k];
+ mpz_mul_2exp (S[k-1], S[k-1], r * h);
mpz_add (S[k-1], S[k-1], S[k]);
mpz_mul (Q[k-1], Q[k-1], Q[k]);
k--;
@@ -144,7 +146,9 @@ mpfr_exp_rational (mpfr_ptr y, mpz_ptr p, long r, int m,
mpz_tdiv_q (S[0], S[0], Q[0]);
mpfr_set_z (y, S[0], MPFR_RNDD);
- MPFR_SET_EXP (y, MPFR_GET_EXP (y) + expo - r * (i - 1) );
+ /* TODO: Check/prove that the following expression doesn't overflow. */
+ expo = MPFR_GET_EXP (y) + expo - r * (i - 1);
+ MPFR_SET_EXP (y, expo);
}
#define shift (GMP_NUMB_BITS/2)
@@ -215,10 +219,10 @@ mpfr_exp_3 (mpfr_ptr y, mpfr_srcptr x, mpfr_rnd_t rnd_mode)
twopoweri = GMP_NUMB_BITS;
/* Allocate tables */
- P = (mpz_t*) (*__gmp_allocate_func) (3*(k+2)*sizeof(mpz_t));
+ P = (mpz_t*) mpfr_allocate_func (3*(k+2)*sizeof(mpz_t));
for (i = 0; i < 3*(k+2); i++)
mpz_init (P[i]);
- mult = (mpfr_prec_t*) (*__gmp_allocate_func) (2*(k+2)*sizeof(mpfr_prec_t));
+ mult = (mpfr_prec_t*) mpfr_allocate_func (2*(k+2)*sizeof(mpfr_prec_t));
/* Particular case for i==0 */
mpfr_extract (uk, x_copy, 0);
@@ -245,8 +249,8 @@ mpfr_exp_3 (mpfr_ptr y, mpfr_srcptr x, mpfr_rnd_t rnd_mode)
/* Clear tables */
for (i = 0; i < 3*(k+2); i++)
mpz_clear (P[i]);
- (*__gmp_free_func) (P, 3*(k+2)*sizeof(mpz_t));
- (*__gmp_free_func) (mult, 2*(k+2)*sizeof(mpfr_prec_t));
+ mpfr_free_func (P, 3*(k+2)*sizeof(mpz_t));
+ mpfr_free_func (mult, 2*(k+2)*sizeof(mpfr_prec_t));
if (shift_x > 0)
{
@@ -283,8 +287,8 @@ mpfr_exp_3 (mpfr_ptr y, mpfr_srcptr x, mpfr_rnd_t rnd_mode)
}
}
- if (mpfr_can_round (shift_x > 0 ? t : tmp, realprec, MPFR_RNDN, MPFR_RNDZ,
- MPFR_PREC(y) + (rnd_mode == MPFR_RNDN)))
+ if (MPFR_CAN_ROUND (shift_x > 0 ? t : tmp, realprec,
+ MPFR_PREC(y), rnd_mode))
{
inexact = mpfr_set (y, shift_x > 0 ? t : tmp, rnd_mode);
if (MPFR_UNLIKELY (scaled && MPFR_IS_PURE_FP (y)))