| Index: gcc/mpfr/const_log2.c
|
| diff --git a/gcc/mpfr/const_log2.c b/gcc/mpfr/const_log2.c
|
| deleted file mode 100644
|
| index 198ac6d677f0bd4c2b65b9c77c3425d14df0bd41..0000000000000000000000000000000000000000
|
| --- a/gcc/mpfr/const_log2.c
|
| +++ /dev/null
|
| @@ -1,192 +0,0 @@
|
| -/* mpfr_const_log2 -- compute natural logarithm of 2
|
| -
|
| -Copyright 1999, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 Free Software Foundation, Inc.
|
| -Contributed by the Arenaire and Cacao projects, INRIA.
|
| -
|
| -This file is part of the GNU MPFR Library.
|
| -
|
| -The GNU MPFR Library is free software; you can redistribute it and/or modify
|
| -it under the terms of the GNU Lesser General Public License as published by
|
| -the Free Software Foundation; either version 2.1 of the License, or (at your
|
| -option) any later version.
|
| -
|
| -The GNU MPFR Library is distributed in the hope that it will be useful, but
|
| -WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
|
| -or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
|
| -License for more details.
|
| -
|
| -You should have received a copy of the GNU Lesser General Public License
|
| -along with the GNU MPFR Library; see the file COPYING.LIB. If not, write to
|
| -the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston,
|
| -MA 02110-1301, USA. */
|
| -
|
| -#define MPFR_NEED_LONGLONG_H
|
| -#include "mpfr-impl.h"
|
| -
|
| -/* Declare the cache */
|
| -MPFR_DECL_INIT_CACHE(__gmpfr_cache_const_log2, mpfr_const_log2_internal);
|
| -
|
| -/* Set User interface */
|
| -#undef mpfr_const_log2
|
| -int
|
| -mpfr_const_log2 (mpfr_ptr x, mp_rnd_t rnd_mode) {
|
| - return mpfr_cache (x, __gmpfr_cache_const_log2, rnd_mode);
|
| -}
|
| -
|
| -/* Auxiliary function: Compute the terms from n1 to n2 (excluded)
|
| - 3/4*sum((-1)^n*n!^2/2^n/(2*n+1)!, n = n1..n2-1).
|
| - Numerator is T[0], denominator is Q[0],
|
| - Compute P[0] only when need_P is non-zero.
|
| - Need 1+ceil(log(n2-n1)/log(2)) cells in T[],P[],Q[].
|
| -*/
|
| -static void
|
| -S (mpz_t *T, mpz_t *P, mpz_t *Q, unsigned long n1, unsigned long n2, int need_P)
|
| -{
|
| - if (n2 == n1 + 1)
|
| - {
|
| - if (n1 == 0)
|
| - mpz_set_ui (P[0], 3);
|
| - else
|
| - {
|
| - mpz_set_ui (P[0], n1);
|
| - mpz_neg (P[0], P[0]);
|
| - }
|
| - if (n1 <= (ULONG_MAX / 4 - 1) / 2)
|
| - mpz_set_ui (Q[0], 4 * (2 * n1 + 1));
|
| - else /* to avoid overflow in 4 * (2 * n1 + 1) */
|
| - {
|
| - mpz_set_ui (Q[0], n1);
|
| - mpz_mul_2exp (Q[0], Q[0], 1);
|
| - mpz_add_ui (Q[0], Q[0], 1);
|
| - mpz_mul_2exp (Q[0], Q[0], 2);
|
| - }
|
| - mpz_set (T[0], P[0]);
|
| - }
|
| - else
|
| - {
|
| - unsigned long m = (n1 / 2) + (n2 / 2) + (n1 & 1UL & n2);
|
| - unsigned long v, w;
|
| -
|
| - S (T, P, Q, n1, m, 1);
|
| - S (T + 1, P + 1, Q + 1, m, n2, need_P);
|
| - mpz_mul (T[0], T[0], Q[1]);
|
| - mpz_mul (T[1], T[1], P[0]);
|
| - mpz_add (T[0], T[0], T[1]);
|
| - if (need_P)
|
| - mpz_mul (P[0], P[0], P[1]);
|
| - mpz_mul (Q[0], Q[0], Q[1]);
|
| -
|
| - /* remove common trailing zeroes if any */
|
| - v = mpz_scan1 (T[0], 0);
|
| - if (v > 0)
|
| - {
|
| - w = mpz_scan1 (Q[0], 0);
|
| - if (w < v)
|
| - v = w;
|
| - if (need_P)
|
| - {
|
| - w = mpz_scan1 (P[0], 0);
|
| - if (w < v)
|
| - v = w;
|
| - }
|
| - /* now v = min(val(T), val(Q), val(P)) */
|
| - if (v > 0)
|
| - {
|
| - mpz_div_2exp (T[0], T[0], v);
|
| - mpz_div_2exp (Q[0], Q[0], v);
|
| - if (need_P)
|
| - mpz_div_2exp (P[0], P[0], v);
|
| - }
|
| - }
|
| - }
|
| -}
|
| -
|
| -/* Don't need to save / restore exponent range: the cache does it */
|
| -int
|
| -mpfr_const_log2_internal (mpfr_ptr x, mp_rnd_t rnd_mode)
|
| -{
|
| - unsigned long n = MPFR_PREC (x);
|
| - mp_prec_t w; /* working precision */
|
| - unsigned long N;
|
| - mpz_t *T, *P, *Q;
|
| - mpfr_t t, q;
|
| - int inexact;
|
| - int ok = 1; /* ensures that the 1st try will give correct rounding */
|
| - unsigned long lgN, i;
|
| - MPFR_ZIV_DECL (loop);
|
| -
|
| - MPFR_LOG_FUNC (("rnd_mode=%d", rnd_mode), ("x[%#R]=%R inex=%d",x,x,inexact));
|
| -
|
| - mpfr_init2 (t, MPFR_PREC_MIN);
|
| - mpfr_init2 (q, MPFR_PREC_MIN);
|
| -
|
| - if (n < 1253)
|
| - w = n + 10; /* ensures correct rounding for the four rounding modes,
|
| - together with N = w / 3 + 1 (see below). */
|
| - else if (n < 2571)
|
| - w = n + 11; /* idem */
|
| - else if (n < 3983)
|
| - w = n + 12;
|
| - else if (n < 4854)
|
| - w = n + 13;
|
| - else if (n < 26248)
|
| - w = n + 14;
|
| - else
|
| - {
|
| - w = n + 15;
|
| - ok = 0;
|
| - }
|
| -
|
| - MPFR_ZIV_INIT (loop, w);
|
| - for (;;)
|
| - {
|
| - N = w / 3 + 1; /* Warning: do not change that (even increasing N!)
|
| - without checking correct rounding in the above
|
| - ranges for n. */
|
| -
|
| - /* the following are needed for error analysis (see algorithms.tex) */
|
| - MPFR_ASSERTD(w >= 3 && N >= 2);
|
| -
|
| - lgN = MPFR_INT_CEIL_LOG2 (N) + 1;
|
| - T = (mpz_t *) (*__gmp_allocate_func) (3 * lgN * sizeof (mpz_t));
|
| - P = T + lgN;
|
| - Q = T + 2*lgN;
|
| - for (i = 0; i < lgN; i++)
|
| - {
|
| - mpz_init (T[i]);
|
| - mpz_init (P[i]);
|
| - mpz_init (Q[i]);
|
| - }
|
| -
|
| - S (T, P, Q, 0, N, 0);
|
| -
|
| - mpfr_set_prec (t, w);
|
| - mpfr_set_prec (q, w);
|
| -
|
| - mpfr_set_z (t, T[0], GMP_RNDN);
|
| - mpfr_set_z (q, Q[0], GMP_RNDN);
|
| - mpfr_div (t, t, q, GMP_RNDN);
|
| -
|
| - for (i = 0; i < lgN; i++)
|
| - {
|
| - mpz_clear (T[i]);
|
| - mpz_clear (P[i]);
|
| - mpz_clear (Q[i]);
|
| - }
|
| - (*__gmp_free_func) (T, 3 * lgN * sizeof (mpz_t));
|
| -
|
| - if (MPFR_LIKELY (ok != 0
|
| - || mpfr_can_round (t, w - 2, GMP_RNDN, rnd_mode, n)))
|
| - break;
|
| -
|
| - MPFR_ZIV_NEXT (loop, w);
|
| - }
|
| - MPFR_ZIV_FREE (loop);
|
| -
|
| - inexact = mpfr_set (x, t, rnd_mode);
|
| -
|
| - mpfr_clear (t);
|
| - mpfr_clear (q);
|
| -
|
| - return inexact;
|
| -}
|
|
|