diff options
author | zimmerma <zimmerma@280ebfd0-de03-0410-8827-d642c229c3f4> | 2000-12-21 17:08:38 +0000 |
---|---|---|
committer | zimmerma <zimmerma@280ebfd0-de03-0410-8827-d642c229c3f4> | 2000-12-21 17:08:38 +0000 |
commit | c74d024d1bfb5cffbd9cf27664f824bb4087dda8 (patch) | |
tree | 5785612341cfeeaab51d555c975ec973b480292a | |
parent | 783f93b38a3d8264f81fb1b303c11c9e9cd9c27e (diff) | |
download | mpfr-c74d024d1bfb5cffbd9cf27664f824bb4087dda8.tar.gz |
k2r -> ansi style
removed #include <math.h> by defining auxiliary functions
fixed several tiny remaining bugs with NaN/Inf
git-svn-id: svn://scm.gforge.inria.fr/svn/mpfr/trunk@925 280ebfd0-de03-0410-8827-d642c229c3f4
-rw-r--r-- | add.c | 8 | ||||
-rw-r--r-- | add_ui.c | 6 | ||||
-rw-r--r-- | agm.c | 67 | ||||
-rw-r--r-- | cmp.c | 33 | ||||
-rw-r--r-- | cmp_ui.c | 10 | ||||
-rw-r--r-- | div.c | 33 | ||||
-rw-r--r-- | div_ui.c | 6 | ||||
-rw-r--r-- | exp.c | 24 | ||||
-rw-r--r-- | exp2.c | 37 | ||||
-rw-r--r-- | exp3.c | 32 | ||||
-rw-r--r-- | get_str.c | 15 | ||||
-rw-r--r-- | inp_str.c | 3 | ||||
-rw-r--r-- | log.c | 30 | ||||
-rw-r--r-- | log2.c | 15 | ||||
-rw-r--r-- | out_str.c | 1 | ||||
-rw-r--r-- | pi.c | 11 | ||||
-rw-r--r-- | pow.c | 42 | ||||
-rw-r--r-- | print_raw.c | 5 | ||||
-rw-r--r-- | reldiff.c | 14 | ||||
-rw-r--r-- | set_d.c | 15 | ||||
-rw-r--r-- | sin_cos.c | 27 | ||||
-rw-r--r-- | sqrt.c | 21 | ||||
-rw-r--r-- | sqrt_ui.c | 1 | ||||
-rw-r--r-- | sub.c | 6 | ||||
-rw-r--r-- | sub_ui.c | 13 | ||||
-rw-r--r-- | tests/tcmp.c | 17 | ||||
-rw-r--r-- | tests/tcmp2.c | 15 | ||||
-rw-r--r-- | tests/tcmp_ui.c | 1 | ||||
-rw-r--r-- | tests/tdump.c | 1 | ||||
-rw-r--r-- | tests/teq.c | 22 | ||||
-rw-r--r-- | tests/tmul_ui.c | 3 | ||||
-rw-r--r-- | tests/trandom.c | 18 | ||||
-rw-r--r-- | tests/ttrunc.c | 7 | ||||
-rw-r--r-- | trunc.c | 5 | ||||
-rw-r--r-- | ui_div.c | 5 | ||||
-rw-r--r-- | ui_sub.c | 2 |
36 files changed, 358 insertions, 213 deletions
@@ -383,9 +383,9 @@ mpfr_add1 (a, b, c, rnd_mode, diff_exp) void #if __STDC__ -mpfr_add(mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mp_rnd_t rnd_mode) +mpfr_add (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mp_rnd_t rnd_mode) #else -mpfr_add(a, b, c, rnd_mode) +mpfr_add (a, b, c, rnd_mode) mpfr_ptr a; mpfr_srcptr b; mpfr_srcptr c; @@ -397,6 +397,8 @@ mpfr_add(a, b, c, rnd_mode) if (MPFR_IS_NAN(b) || MPFR_IS_NAN(c)) { MPFR_SET_NAN(a); return; } + + MPFR_CLEAR_NAN(a); if (MPFR_IS_INF(b)) { @@ -425,6 +427,8 @@ mpfr_add(a, b, c, rnd_mode) return; } + MPFR_CLEAR_INF(a); /* clear Inf flag */ + if (!MPFR_NOTZERO(b)) { mpfr_set(a, c, rnd_mode); return; } if (!MPFR_NOTZERO(c)) { mpfr_set(a, b, rnd_mode); return; } @@ -28,9 +28,9 @@ MA 02111-1307, USA. */ void #if __STDC__ -mpfr_add_ui(mpfr_ptr y, mpfr_srcptr x, unsigned long int u, mp_rnd_t rnd_mode) +mpfr_add_ui (mpfr_ptr y, mpfr_srcptr x, unsigned long int u, mp_rnd_t rnd_mode) #else -mpfr_add_ui(y, x, u, rnd_mode) +mpfr_add_ui (y, x, u, rnd_mode) mpfr_ptr y; mpfr_srcptr x; unsigned long int u; @@ -49,4 +49,6 @@ mpfr_add_ui(y, x, u, rnd_mode) mpfr_add(y, x, uu, rnd_mode); } + else + mpfr_set (y, x, rnd_mode); } @@ -20,17 +20,55 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "mpfr.h" #include "mpfr-impl.h" +/* returns ceil(log(d)/log(2)) */ +long _mpfr_ceil_log2 (double d) +{ + long exp; + union ieee_double_extract x; + + x.d = d; + exp = x.s.exp - 1023; + x.s.exp = 1023; /* value for 1 <= d < 2 */ + if (x.d != 1.0) exp++; + return exp; +} + +/* returns floor(log(d)/log(2)) */ +long _mpfr_floor_log2 (double d) +{ + long exp; + union ieee_double_extract x; + + x.d = d; + exp = x.s.exp - 1023; + x.s.exp = 1023; /* value for 1 <= d < 2 */ + return exp; +} + +/* returns y >= 2^d */ +double _mpfr_ceil_exp2 (double d) +{ + long exp; + union ieee_double_extract x; + + exp = (long) d; + if (d != (double) exp) exp++; + /* now exp = ceil(d) */ + x.d = 1.0; + x.s.exp = 1023 + exp; + return x.d; +} + void #ifdef __STDC__ -mpfr_agm(mpfr_ptr r, mpfr_srcptr op2, mpfr_srcptr op1, mp_rnd_t rnd_mode) +mpfr_agm (mpfr_ptr r, mpfr_srcptr op2, mpfr_srcptr op1, mp_rnd_t rnd_mode) #else -mpfr_agm(r, op2, op1, rnd_mode) +mpfr_agm (r, op2, op1, rnd_mode) mpfr_ptr r; mpfr_srcptr op2; mpfr_srcptr op1; @@ -45,16 +83,21 @@ mpfr_agm(r, op2, op1, rnd_mode) TMP_DECL(marker1); TMP_DECL(marker2); - /* TODO : CHECK FOR INFINITY. */ - /* If a or b is NaN, the result is NaN */ if (MPFR_IS_NAN(op1) || MPFR_IS_NAN(op2)) { MPFR_SET_NAN(r); return; } - - /* If a or b is negative, the result is NaN */ + /* If a or b is negative (including -Infinity), the result is NaN */ if ((MPFR_SIGN(op1) < 0) || (MPFR_SIGN(op2) < 0)) { MPFR_SET_NAN(r); return; } + + MPFR_CLEAR_NAN(r); + + /* If a or b is +Infinity, the result is +Infinity */ + if (MPFR_IS_INF(op1) || MPFR_IS_INF(op2)) + { MPFR_SET_INF(r); MPFR_SET_SAME_SIGN(r, op1); return; } + + MPFR_CLEAR_INF(r); /* If a or b is 0, the result is 0 */ if ((MPFR_NOTZERO(op1) && MPFR_NOTZERO(op2)) == 0) @@ -66,7 +109,6 @@ mpfr_agm(r, op2, op1, rnd_mode) q = MPFR_PREC(r); p = q + 15; - /* Initialisations */ go_on=1; @@ -105,16 +147,19 @@ mpfr_agm(r, op2, op1, rnd_mode) eq=0; - err=ceil((3.0/2.0*log((double)p)/LOG2+1.0)*exp(-(double)p*LOG2)+3.0*exp(-2.0*(double)p*uo*LOG2/(vo-uo))); + err=1 + (int) ((3.0/2.0*_mpfr_ceil_log2(p)+1.0)*_mpfr_ceil_exp2(-(double)p) + +3.0*_mpfr_ceil_exp2(-2.0*(double)p*uo/(vo-uo))); if(p-err-3<=q) { p=q+err+4; - err=ceil((3.0/2.0*log((double)p)/LOG2+1.0)*exp(-(double)p*LOG2)+3.0*exp(-2.0*(double)p*uo*LOG2/(vo-uo))); + err= 1 + + (int) ((3.0/2.0*_mpfr_ceil_log2(p)+1.0)*_mpfr_ceil_exp2(-(double)p) + +3.0*_mpfr_ceil_exp2(-2.0*(double)p*uo/(vo-uo))); } /* Calculus of un and vn */ while (eq<=p-2) { mpfr_mul(tmp,u,v,GMP_RNDN); - mpfr_sqrt(tmpu,tmp,GMP_RNDN); + mpfr_sqrt (tmpu, tmp, GMP_RNDN); mpfr_add(tmp,u,v,GMP_RNDN); mpfr_div_2exp(tmpv,tmp,1,GMP_RNDN); mpfr_set(u,tmpu,GMP_RNDN); @@ -42,9 +42,9 @@ of b and c, i.e. one plus the number of bits shifts to align b and c /* compares b and sign(s)*c */ int #if __STDC__ -mpfr_cmp3 ( mpfr_srcptr b, mpfr_srcptr c, long int s) +mpfr_cmp3 (mpfr_srcptr b, mpfr_srcptr c, long int s) #else -mpfr_cmp3(b, c, s) +mpfr_cmp3 (b, c, s) mpfr_srcptr b; mpfr_srcptr c; long int s; @@ -54,9 +54,17 @@ mpfr_cmp3(b, c, s) unsigned long bn, cn; mp_limb_t *bp, *cp; + if (MPFR_IS_NAN(b) || MPFR_IS_NAN(c)) return 1; + + if (MPFR_IS_INF(b)) { + if (MPFR_IS_INF(c) && (MPFR_SIGN(b) * s * MPFR_SIGN(c) > 0)) + return 0; + else + return MPFR_SIGN(b); + } + if (!MPFR_NOTZERO(b)) { if (!MPFR_NOTZERO(c)) return 0; else return -(s*MPFR_SIGN(c)); - /*TODO: bug ou feature ? s pas pris en compte... */ } else if (!MPFR_NOTZERO(c)) return MPFR_SIGN(b); @@ -64,9 +72,6 @@ mpfr_cmp3(b, c, s) if (s<0) return(MPFR_SIGN(b)); /* now signs are equal */ - if (MPFR_IS_INF(b)) - return MPFR_SIGN(b) * !MPFR_IS_INF(c); - diff_exp = MPFR_EXP(b)-MPFR_EXP(c); s = (MPFR_SIGN(b) > 0) ? 1 : -1; @@ -191,7 +196,7 @@ mpfr_cmp2(b, c) z = 0; /* number of possibly cancelled bits - 1 */ /* thus result is either k if low(b) >= low(c) or k+z+1 if low(b) < low(c) */ - if (d > BITS_PER_MP_LIMB) { return k; } + if (d > BITS_PER_MP_LIMB) return k; while (bn >= 0) /* the next limb of b to be considered is b[bn] */ { @@ -221,16 +226,20 @@ mpfr_cmp2(b, c) else z += BITS_PER_MP_LIMB; } - if (cn >= 0) - count_leading_zeros(cc, ~(cp[cn--] << (BITS_PER_MP_LIMB - d))); - else { cc = 0; } + /* warning: count_leading_zeros doesn't work with zero */ + if ((cn >= 0) && (u = ~(cp[cn--] << (BITS_PER_MP_LIMB - d)))) + count_leading_zeros(cc, u); + else + cc = 0; - k += cc + d; /* here d=0 or d=1: if d=1, we have one more cancelled bit */ + /* here d=0 or d=1: if d=1, we have one more cancelled bit if we don't + shift cp[cn] */ + k += cc; if (cc < d) return k; while (cn >= 0 && !~cp[cn]) { z += BITS_PER_MP_LIMB; cn--; } /* now either cn<0 or cp[cn] is not 111...111 */ - if (cn >= 0) { count_leading_zeros(cc, ~cp[cn]); return (k + z + cc); } + if (cn >= 0) { count_leading_zeros(cc, ~cp[cn]); return (k + d+ z + cc); } return k; /* We **need** that the nonsignificant limbs of c are set to zero there */ @@ -33,7 +33,7 @@ MA 02111-1307, USA. */ int #if __STDC__ -mpfr_cmp_ui_2exp ( mpfr_srcptr b, unsigned long int i, int f ) +mpfr_cmp_ui_2exp (mpfr_srcptr b, unsigned long int i, int f ) #else mpfr_cmp_ui_2exp (b, i, f) mpfr_srcptr b; @@ -43,6 +43,8 @@ mpfr_cmp_ui_2exp (b, i, f) { int e, k, bn; mp_limb_t c, *bp; + if (MPFR_IS_NAN(b)) return 1; + if (MPFR_SIGN(b) < 0) return -1; /* now b>=0 */ else if (MPFR_IS_INF(b)) return 1; @@ -79,9 +81,9 @@ mpfr_cmp_ui_2exp (b, i, f) int #if __STDC__ -mpfr_cmp_si_2exp ( mpfr_srcptr b, long int i, int f ) +mpfr_cmp_si_2exp (mpfr_srcptr b, long int i, int f ) #else -mpfr_cmp_si_2exp(b, i, f) +mpfr_cmp_si_2exp (b, i, f) mpfr_srcptr b; long int i; int f; @@ -89,6 +91,8 @@ mpfr_cmp_si_2exp(b, i, f) { int e, k, bn, si; mp_limb_t c, *bp; + if (MPFR_IS_NAN(b)) return 1; + si = (i<0) ? -1 : 1; /* sign of i */ if (MPFR_SIGN(b) * i < 0 || MPFR_IS_INF(b)) return MPFR_SIGN(b); /* both signs differ */ @@ -19,7 +19,6 @@ along with the MPFR Library; see the file COPYING.LIB. If not, write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ -#include <math.h> #include <stdio.h> #include <stdlib.h> #include "gmp.h" @@ -55,30 +54,32 @@ mpfr_div (r, u, v, rnd_mode) TMP_DECL (marker); if (MPFR_IS_NAN(u) || MPFR_IS_NAN(v)) { MPFR_SET_NAN(r); return; } + + MPFR_CLEAR_NAN(r); + if (MPFR_IS_INF(u)) { if (MPFR_IS_INF(v)) - { - MPFR_SET_NAN(r); return; - } + MPFR_SET_NAN(r); else { MPFR_SET_INF(r); if (MPFR_SIGN(r) != MPFR_SIGN(u) * MPFR_SIGN(v)) - { MPFR_CHANGE_SIGN(r); } - return; - } - } - else - { - if (MPFR_IS_INF(v)) - { - MPFR_SET_ZERO(r); - if (MPFR_SIGN(r) != MPFR_SIGN(u) * MPFR_SIGN(v)) - { MPFR_CHANGE_SIGN(r); } - return; + MPFR_CHANGE_SIGN(r); } + return; } + else + if (MPFR_IS_INF(v)) + { + MPFR_CLEAR_INF(r); + MPFR_SET_ZERO(r); + if (MPFR_SIGN(r) != MPFR_SIGN(u) * MPFR_SIGN(v)) + MPFR_CHANGE_SIGN(r); + return; + } + + MPFR_CLEAR_INF(r); /* clear Inf flag */ usize = (MPFR_PREC(u) - 1)/BITS_PER_MP_LIMB + 1; vsize = (MPFR_PREC(v) - 1)/BITS_PER_MP_LIMB + 1; @@ -45,10 +45,12 @@ mpfr_div_ui(y, x, u, rnd_mode) if (MPFR_IS_NAN(x)) { MPFR_SET_NAN(y); return 1; } + MPFR_CLEAR_NAN(y); /* clear NaN flag */ + if (MPFR_IS_INF(x)) { MPFR_SET_INF(y); - if (MPFR_SIGN(y) * MPFR_SIGN(x) * u < 0) + if (MPFR_SIGN(y) * MPFR_SIGN(x) < 0) /* consider u=0 as +0 */ MPFR_CHANGE_SIGN(y); return 0; /* TODO: semantique de la division par un zero entier ? signe ? */ @@ -64,6 +66,8 @@ mpfr_div_ui(y, x, u, rnd_mode) } } + MPFR_CLEAR_INF(y); + TMP_MARK(marker); xn = (MPFR_PREC(x)-1)/BITS_PER_MP_LIMB + 1; yn = (MPFR_PREC(y)-1)/BITS_PER_MP_LIMB + 1; @@ -20,7 +20,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "mpfr.h" @@ -38,9 +37,9 @@ extern int mpfr_exp3 _PROTO((mpfr_ptr, mpfr_srcptr, mp_rnd_t)); */ int #if __STDC__ -mpfr_exp(mpfr_ptr y, mpfr_srcptr x, mp_rnd_t rnd_mode) +mpfr_exp (mpfr_ptr y, mpfr_srcptr x, mp_rnd_t rnd_mode) #else -mpfr_exp(y, x, rnd_mode) +mpfr_exp (y, x, rnd_mode) mpfr_ptr y; mpfr_srcptr x; mp_rnd_t rnd_mode; @@ -49,14 +48,25 @@ mpfr_exp(y, x, rnd_mode) int expx, precy; if (MPFR_IS_NAN(x)) { MPFR_SET_NAN(y); return 1; } + + MPFR_CLEAR_NAN(y); + if (MPFR_IS_INF(x)) { - if (MPFR_SIGN(x) > 0) - { MPFR_SET_INF(y); if (MPFR_SIGN(y) == -1) { MPFR_CHANGE_SIGN(y); } } - else - { MPFR_SET_ZERO(y); if (MPFR_SIGN(y) == -1) { MPFR_CHANGE_SIGN(y); } } + if (MPFR_SIGN(x) > 0) { + MPFR_SET_INF(y); + if (MPFR_SIGN(y) < 0) MPFR_CHANGE_SIGN(y); + } + else { + MPFR_CLEAR_INF(y); + MPFR_SET_ZERO(y); + if (MPFR_SIGN(y) < 0) MPFR_CHANGE_SIGN(y); + } /* TODO: conflits entre infinis et zeros ? */ + return 1; } + + MPFR_CLEAR_INF(y); if (!MPFR_NOTZERO(x)) { mpfr_set_ui(y, 1, GMP_RNDN); return 0; } @@ -21,7 +21,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "mpfr.h" @@ -33,6 +32,31 @@ mp_exp_t mpz_normalize (mpz_t, mpz_t, int); int mpz_normalize2 (mpz_t, mpz_t, int, int); int mpfr_exp2 (mpfr_ptr, mpfr_srcptr, mp_rnd_t); +/* returns floor(sqrt(n)) */ +unsigned long _mpfr_isqrt (unsigned long n) +{ + unsigned long s; + + s = 1; + do { + s = (s + n / s) / 2; + } while (!(s*s <= n && n <= s*(s+2))); + return s; +} + +/* returns floor(n^(1/3)) */ +unsigned long _mpfr_cuberoot (unsigned long n) +{ + double s, is; + + s = 1.0; + do { + s = (2*s*s*s + (double) n) / (3*s*s); + is = (double) ((int) s); + } while (!(is*is*is <= (double) n && (double) n < (is+1)*(is+1)*(is+1))); + return (unsigned long) is; +} + #define SWITCH 100 /* number of bits to switch from O(n^(1/2)*M(n)) method to O(n^(1/3)*M(n)) method */ @@ -139,15 +163,14 @@ mpfr_exp2 (y, x, rnd_mode) mpfr_sub_one_ulp(y); return 1; } - n = (int) floor(mpfr_get_d(x)/LOG2); + n = (int) (mpfr_get_d(x) / LOG2); /* for the O(n^(1/2)*M(n)) method, the Taylor series computation of n/K terms costs about n/(2K) multiplications when computed in fixed point */ - K = (int) (precy<SWITCH) ? sqrt( (double) (precy + 1)/ 2.0 ) - : pow( 4.0 * (double) precy, 1.0/3.0); + K = (precy<SWITCH) ? _mpfr_isqrt((precy + 1) / 2) : _mpfr_cuberoot (4*precy); l = (precy-1)/K + 1; - err = K + (int) ceil(log(2.0*(double)l+18.0)/LOG2); + err = K + (int) _mpfr_ceil_log2 (2.0 * (double) l + 18.0); /* add K extra bits, i.e. failure probability <= 1/2^K = O(1/precy) */ q = precy + err + K + 3; mpfr_init2(r, q); mpfr_init2(s, q); mpfr_init2(t, q); @@ -161,7 +184,7 @@ mpfr_exp2 (y, x, rnd_mode) /* if n<0, we have to get an upper bound of log(2) in order to get an upper bound of r = x-n*log(2) */ - mpfr_const_log2(s, (n>=0) ? GMP_RNDZ : GMP_RNDU); + mpfr_const_log2 (s, (n>=0) ? GMP_RNDZ : GMP_RNDU); #ifdef DEBUG printf("n=%d log(2)=",n); mpfr_print_raw(s); putchar('\n'); #endif @@ -323,7 +346,7 @@ mpfr_exp2_aux2(s, r, q, exps) /* estimate value of l */ l = q / (-MPFR_EXP(r)); - m = (int) sqrt((double) l); + m = (int) _mpfr_isqrt (l); /* we access R[2], thus we need m >= 2 */ if (m < 2) m = 2; TMP_MARK(marker); @@ -20,7 +20,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "mpfr.h" @@ -28,25 +27,11 @@ MA 02111-1307, USA. */ /* #define DEBUG */ -int mylog2 (int); int mpfr_exp_rational (mpfr_ptr, mpz_srcptr, int, int); int mpfr_exp3 (mpfr_ptr, mpfr_srcptr, mp_rnd_t); int #if __STDC__ -mylog2 (int x) -#else -mylog2 (x) - int x; -#endif -{ - int i = 0; - for ( ; x != 1; x >>= 1, i++) ; - return i; -} - -int -#if __STDC__ mpfr_exp_rational (mpfr_ptr y, mpz_srcptr p, int r, int m) #else mpfr_exp_rational (y, p, r, m) @@ -103,7 +88,7 @@ mpfr_exp_rational (y, p, r, m) l = 0; accu = 0; while (k > 0){ - mpz_mul(S[k], S[k], ptoj[mylog2(nb_terms[k])]); + mpz_mul(S[k], S[k], ptoj[_mpfr_ceil_log2((double) nb_terms[k])]); mpz_mul(S[k-1], S[k-1], P[k]); accu += nb_terms[k]; mpz_mul_2exp(S[k-1], S[k-1], r* accu); @@ -189,13 +174,9 @@ mpfr_exp3 (y, x, rnd_mode) /* Decomposer x */ /* on commence par ecrire x = 1.xxxxxxxxxxxxx ----- k bits -- */ - prec_x = (int) ceil(log - ((double) (MPFR_PREC(x)) / (double) BITS_PER_MP_LIMB) - /LOG2); + prec_x = _mpfr_ceil_log2 ((double) (MPFR_PREC(x)) / BITS_PER_MP_LIMB); if (prec_x < 0) prec_x = 0; - logn = (int) ceil(log - ((double) prec_x+MPFR_PREC(y)) - /LOG2); + logn = _mpfr_ceil_log2 ((double) prec_x + MPFR_PREC(y)); if (logn < 2) logn = 2; ttt = MPFR_EXP(x); mpfr_init2(x_copy,MPFR_PREC(x)); @@ -211,10 +192,9 @@ mpfr_exp3 (y, x, rnd_mode) mpz_init (uk); while (!good){ Prec = realprec+shift+2+shift_x; - k = (int) ceil(log - ((double) (Prec) / (double) BITS_PER_MP_LIMB) - /LOG2); - /* Maintenant, il faut extraire : */ + k = _mpfr_ceil_log2 ((double) Prec / BITS_PER_MP_LIMB); + + /* now we have to extract */ mpfr_init2(t, Prec); mpfr_init2(tmp, Prec); mpfr_set_ui(tmp,1,GMP_RNDN); @@ -21,7 +21,6 @@ MA 02111-1307, USA. */ /* #define DEBUG */ -#include <math.h> #include <stdio.h> #include <stdlib.h> #include <string.h> @@ -99,17 +98,22 @@ char *mpfr_get_str(str, expptr, base, n, op, rnd_mode) /* first determines the exponent */ e = MPFR_EXP(op); - d = fabs(mpfr_get_d2(op, 0)); + d = ABS(mpfr_get_d2(op, 0)); /* the absolute value of op is between 1/2*2^e and 2^e */ /* the output exponent f is such that base^(f-1) <= |op| < base^f i.e. f = 1 + floor(log(|op|)/log(base)) = 1 + floor((log(|m|)+e*log(2))/log(base)) */ - f = 1 + (int) floor((log(d)+(double)e*LOG2)/log((double)base)); + /* f = 1 + (int) floor((log(d)/LOG2+(double)e)*LOG2/log((double)base)); */ + d = ((double) e + _mpfr_floor_log2((double) d)) + * __mp_bases[base].chars_per_bit_exactly; + /* warning: (int) d rounds towards 0 */ + f = (int) d; /* f equals floor(d) if d >= 0 and ceil(d) if d < 0 */ + if (f <= d) f++; if (n==0) { /* performs exact rounding, i.e. returns y such that for GMP_RNDU for example, we have: x*2^(e-p) <= y*base^(f-n) */ - n = (int) ((double)MPFR_PREC(op)*LOG2/log((double)base)); + n = (int) (__mp_bases[base].chars_per_bit_exactly * MPFR_PREC(op)); if (n==0) n=1; } #ifdef DEBUG @@ -118,7 +122,8 @@ char *mpfr_get_str(str, expptr, base, n, op, rnd_mode) /* now the first n digits of the mantissa are obtained from rnd(op*base^(n-f)) */ if (pow2) prec = n*pow2; - else prec = (long) ceil((double)n*log((double)base)/LOG2); + else + prec = 1 + (long) ((double) n / __mp_bases[base].chars_per_bit_exactly); #ifdef DEBUG printf("prec=%d\n", prec); #endif @@ -23,8 +23,9 @@ MA 02111-1307, USA. */ #include <stdio.h> #include <ctype.h> #include "gmp.h" -#include "mpfr.h" #include "gmp-impl.h" +#include "mpfr.h" +#include "mpfr-impl.h" size_t #if __STDC__ @@ -20,7 +20,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "mpfr.h" @@ -44,9 +43,9 @@ MA 02111-1307, USA. */ int #if __STDC__ -mpfr_log(mpfr_ptr r, mpfr_srcptr a, mp_rnd_t rnd_mode) +mpfr_log (mpfr_ptr r, mpfr_srcptr a, mp_rnd_t rnd_mode) #else -mpfr_log(r, a, rnd_mode) +mpfr_log (r, a, rnd_mode) mpfr_ptr r; mpfr_srcptr a; mp_rnd_t rnd_mode; @@ -60,27 +59,30 @@ mpfr_log(r, a, rnd_mode) TMP_DECL(marker); /* If a is NaN, the result is NaN */ - if (MPFR_IS_NAN(a)) - { MPFR_CLEAR_FLAGS(r); MPFR_SET_NAN(r); return 1; } - - if (MPFR_IS_ZERO(a)) + if (MPFR_IS_NAN(a)) { + MPFR_SET_NAN(r); + return 1; + } + + MPFR_CLEAR_NAN(r); + + /* check for infinity before zero */ + if (MPFR_IS_INF(a)) { - MPFR_CLEAR_FLAGS(r); MPFR_SET_INF(r); - if (MPFR_SIGN(r) != -1) { MPFR_CHANGE_SIGN(r); } + if (MPFR_SIGN(r) > 0) MPFR_CHANGE_SIGN(r); return 1; } - if (MPFR_IS_INF(a)) + if (MPFR_IS_ZERO(a)) { - MPFR_CLEAR_FLAGS(r); MPFR_SET_INF(r); - if (MPFR_SIGN(r) != 1) { MPFR_CHANGE_SIGN(r); } + if (MPFR_SIGN(r) < 0) MPFR_CHANGE_SIGN(r); return 1; } /* Now we can clear the flags without damage even if r == a */ - MPFR_CLEAR_FLAGS(r); + MPFR_CLEAR_INF(r); /* If a is 1, the result is 0 */ if (mpfr_cmp_ui_2exp(a,1,0)==0){ @@ -108,7 +110,7 @@ mpfr_log(r, a, rnd_mode) printf("p=%d\n", p); #endif /* Calculus of m (depends on p) */ - m=(int) ceil(((double) p)/2.0) -MPFR_EXP(a)+1; + m = (p + 1) / 2 - MPFR_EXP(a) + 1; /* All the mpfr_t needed have a precision of p */ TMP_MARK(marker); @@ -20,7 +20,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "longlong.h" @@ -72,14 +71,10 @@ mpfr_const_aux_log2 (mylog, rnd_mode) int prec_x; mpz_init(cst); - logn = (int) ceil(log - ((double) MPFR_PREC(mylog)) - /LOG2); + logn = _mpfr_ceil_log2 ((double) MPFR_PREC(mylog)); prec_x = prec_i_want + logn; while (!good){ - prec = (int) ceil(log - ((double) (prec_x)) - /LOG2); + prec = _mpfr_ceil_log2 ((double) prec_x); mpfr_init2(tmp1, prec_x); mpfr_init2(result, prec_x); mpfr_init2(tmp2, prec_x); @@ -132,9 +127,9 @@ mpfr_const_aux_log2 (mylog, rnd_mode) */ void #if __STDC__ -mpfr_const_log2(mpfr_ptr x, mp_rnd_t rnd_mode) +mpfr_const_log2 (mpfr_ptr x, mp_rnd_t rnd_mode) #else -mpfr_const_log2(x, rnd_mode) mpfr_ptr x; mp_rnd_t rnd_mode; +mpfr_const_log2 (x, rnd_mode) mpfr_ptr x; mp_rnd_t rnd_mode; #endif { int N, oldN, k, precx; mpz_t s, t, u; @@ -156,7 +151,7 @@ mpfr_const_log2(x, rnd_mode) mpfr_ptr x; mp_rnd_t rnd_mode; N=2; do { oldN = N; - N = precx + (int)ceil(log((double)N)/LOG2); + N = precx + _mpfr_ceil_log2 ((double) N); } while (N != oldN); mpz_init_set_ui(s,0); mpz_init(u); @@ -20,7 +20,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include <string.h> #include "gmp.h" #include "gmp-impl.h" @@ -20,7 +20,6 @@ the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ #include <stdio.h> -#include <math.h> #include "gmp.h" #include "gmp-impl.h" #include "longlong.h" @@ -59,15 +58,11 @@ mpfr_pi_machin3 (mylog, rnd_mode) mpz_t cst; MPFR_CLEAR_FLAGS(mylog); - logn = (int) ceil(log - ((double) MPFR_PREC(mylog)) - /LOG2); + logn = _mpfr_ceil_log2 ((double) MPFR_PREC(mylog)); prec_x = prec_i_want + logn + 5; mpz_init(cst); while (!good){ - prec = (int) ceil(log - ((double) (prec_x)) - /LOG2); + prec = _mpfr_ceil_log2 ((double) prec_x); mpfr_init2(tmp1, prec_x); mpfr_init2(tmp2, prec_x); @@ -189,7 +184,7 @@ mpfr_const_pi(x, rnd_mode) N=1; do { oldN = N; - N = (prec+3)/4 + (int)ceil(log((double)N+1.0)/LOG2); + N = (prec+3)/4 + _mpfr_ceil_log2((double) N + 1.0); } while (N != oldN); mpz_init(pi); mpz_init(num); mpz_init(den); mpz_init(d3); mpz_init(d2); mpz_init(tmp); @@ -20,7 +20,6 @@ along with the MPFR Library; see the file COPYING.LIB. If not, write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ -#include <math.h> #include <stdio.h> #include "gmp.h" #include "mpfr.h" @@ -38,29 +37,48 @@ mpfr_pow_ui (x, y, n, rnd) mp_rnd_t rnd; #endif { - long int i; double err; + long int i; + unsigned long m; + double err; + mpfr_t ycopy; + + if (MPFR_IS_NAN(y)) { MPFR_SET_NAN(x); return 0; } + + MPFR_CLEAR_NAN(x); if (MPFR_IS_INF(y)) { - MPFR_CLEAR_FLAGS(x); if (n == 0) { MPFR_SET_NAN(x); return 0; } - else if ((MPFR_SIGN(y) < 0 && (MPFR_SIGN(x) > 0 || n % 2 == 0)) || - (MPFR_SIGN(y) > 0 && (MPFR_SIGN(x) < 0 && n % 2 == 1))) + /* Inf^n = Inf, (-Inf)^n = Inf for n even, -Inf for n odd */ + else if ((MPFR_SIGN(y) < 0 && (MPFR_SIGN(x) > 0 || n % 2 == 1)) || + (MPFR_SIGN(y) > 0 && (MPFR_SIGN(x) < 0))) MPFR_CHANGE_SIGN(x); MPFR_SET_INF(x); return 0; } - if (MPFR_IS_NAN(y)) { MPFR_CLEAR_FLAGS(x); MPFR_SET_NAN(x); return 0; } + + MPFR_CLEAR_INF(x); - MPFR_CLEAR_FLAGS(x); if (n==0) { mpfr_set_ui(x, 1, rnd); return 0; } - mpfr_set(x, y, rnd); err = 1.0; - for (i=0;(1<<i)<=n;i++); + + if (x == y) { + mpfr_init2 (ycopy, MPFR_PREC(y)); + mpfr_set (ycopy, y, GMP_RNDN); /* exact */ + } + + mpfr_set(x, y, rnd); + err = 1.0; + for (m=n, i=0; m; i++, m>>=1); /* now 2^(i-1) <= n < 2^i */ for (i-=2; i>=0; i--) { mpfr_mul(x, x, x, rnd); err = 2.0*err+1.0; - if (n & (1<<i)) { mpfr_mul(x, x, y, rnd); err += 1.0; } + if (n & (1<<i)) { + mpfr_mul(x, x, (x == y) ? ycopy : y, rnd); + err += 1.0; + } } - return (int) ceil(log(err)/LOG2); + + if (x == y) mpfr_clear (ycopy); + return _mpfr_ceil_log2 (err); } /* sets x to y^n, and returns ceil(log2(max ulp error)) */ @@ -97,5 +115,5 @@ mpfr_ui_pow_ui (x, y, n, rnd) err = err + 1.0; } } - return (int) ceil(log(err)/LOG2); + return _mpfr_ceil_log2 (err); } diff --git a/print_raw.c b/print_raw.c index 38c1a4f59..69e682757 100644 --- a/print_raw.c +++ b/print_raw.c @@ -81,7 +81,10 @@ mpfr_print_raw(x) else if (MPFR_IS_INF(x)) { if (MPFR_SIGN(x) == 1) { printf("Inf"); } else printf("-Inf"); } - else if (!MPFR_NOTZERO(x)) printf("0"); + else if (!MPFR_NOTZERO(x)) { + if (MPFR_SIGN(x) < 0) printf("-"); + printf("0"); + } else { /* 3 char for sign + 0 + binary point + MPFR_ABSSIZE(x) * BITS_PER_MP_LIMB for mantissa @@ -21,8 +21,9 @@ MA 02111-1307, USA. */ #include <stdio.h> #include "gmp.h" -#include "mpfr.h" #include "gmp-impl.h" +#include "mpfr.h" +#include "mpfr-impl.h" /* reldiff(b, c) = abs(b-c)/b */ void @@ -36,6 +37,8 @@ mpfr_reldiff(a, b, c, rnd_mode) mp_rnd_t rnd_mode; #endif { + mpfr_t b_copy; + if (MPFR_IS_NAN(b) || MPFR_IS_NAN(c)) { MPFR_CLEAR_FLAGS(a); MPFR_SET_NAN(a); return; } if (MPFR_IS_INF(b)) @@ -58,9 +61,16 @@ mpfr_reldiff(a, b, c, rnd_mode) /* TODO: faire preciser la SEMANTIQUE DE CE FOUTOIR. */ mpfr_set_ui(a, MPFR_SIGN(c), rnd_mode); else { + if (a == b) { + mpfr_init2 (b_copy, MPFR_PREC(b)); + mpfr_set (b_copy, b, GMP_RNDN); + } + mpfr_sub(a, b, c, rnd_mode); mpfr_abs(a, a, rnd_mode); /* for compatibility with MPF */ - mpfr_div(a, a, b, rnd_mode); + mpfr_div(a, a, (a == b) ? b_copy : b, rnd_mode); + + if (a == b) mpfr_clear (b_copy); } } @@ -242,9 +242,9 @@ __mpfr_scale2 (d, exp) void #if __STDC__ -mpfr_set_d(mpfr_ptr r, double d, mp_rnd_t rnd_mode) +mpfr_set_d (mpfr_ptr r, double d, mp_rnd_t rnd_mode) #else -mpfr_set_d(r, d, rnd_mode) +mpfr_set_d (r, d, rnd_mode) mpfr_ptr r; double d; mp_rnd_t rnd_mode; @@ -253,7 +253,16 @@ mpfr_set_d(r, d, rnd_mode) int signd, sizer; unsigned int cnt; MPFR_CLEAR_FLAGS(r); - if (d == 0) { MPFR_SET_ZERO(r); return; } + if (d == 0) { + union ieee_double_extract x; + MPFR_SET_ZERO(r); + /* set correct sign */ + x.d = d; + if (((x.s.sig==1) && (MPFR_SIGN(r)>0)) + || ((x.s.sig==0) && (MPFR_SIGN(r)<0))) + MPFR_CHANGE_SIGN(r); + return; + } else if (DOUBLE_ISNAN(d)) { MPFR_SET_NAN(r); return; } else if (DOUBLE_ISINF(d)) { @@ -19,7 +19,6 @@ along with the MPFR Library; see the file COPYING.LIB. If not, write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ -#include <math.h> #include <stdio.h> #include <stdlib.h> #include "gmp.h" @@ -57,10 +56,10 @@ int mpfr_sin_cos (mpfr_ptr sinus, mpfr_ptr cosinus, mpfr_srcptr x, mp_rnd_t rnd_mode) #else mpfr_sin_cos (sinus, cosinus, x, rnd_mode) -mpfr_ptr sinus; -mpfr_ptr cosinus; -mpfr_srcptr x; -mp_rnd_t rnd_mode; + mpfr_ptr sinus; + mpfr_ptr cosinus; + mpfr_srcptr x; + mp_rnd_t rnd_mode; #endif { mpfr_t t_sin, t_cos; @@ -85,6 +84,11 @@ mp_rnd_t rnd_mode; int tmp_factor; int tmpi; + if (sinus == cosinus) { + fprintf (stderr, "Error in mpfr_sin_cos: 1st and 2nd operands must be different\n"); + exit (1); + } + if (MPFR_IS_NAN(x) || MPFR_IS_INF(x)) { MPFR_SET_NAN(sinus); MPFR_SET_NAN(cosinus); @@ -97,9 +101,7 @@ mp_rnd_t rnd_mode; return 0; /* exact results */ } - prec_x = (int) ceil(log - ((double) (MPFR_PREC(x)) / (double) BITS_PER_MP_LIMB) - /LOG2); + prec_x = _mpfr_ceil_log2 ((double) MPFR_PREC(x) / BITS_PER_MP_LIMB); ttt = MPFR_EXP(x); mpfr_init2(x_copy, MPFR_PREC(x)); mpfr_set(x_copy,x,GMP_RNDD); @@ -113,17 +115,16 @@ mp_rnd_t rnd_mode; } target_prec = MPFR_PREC(sinus); if (MPFR_PREC(cosinus) > target_prec) target_prec = MPFR_PREC(cosinus); - logn = (int) ceil(log((double)target_prec)/LOG2); + logn = _mpfr_ceil_log2 ((double) target_prec); if (logn < 2) logn = 2; factor = logn; realprec = target_prec + logn; mpz_init (uk); while (!good) { Prec = realprec + 2*shift + 2 + shift_x + factor; - k = (int) ceil(log - ((double) (Prec) / (double) BITS_PER_MP_LIMB) - /LOG2); - /* Maintenant, il faut extraire : */ + k = _mpfr_ceil_log2 ((double) Prec / BITS_PER_MP_LIMB); + + /* now we extract */ mpfr_init2(t_sin, Prec); mpfr_init2(t_cos, Prec); mpfr_init2(tmp, Prec); @@ -19,7 +19,6 @@ along with the MPFR Library; see the file COPYING.LIB. If not, write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ -#include <math.h> #include <stdio.h> #include <stdlib.h> #include "gmp.h" @@ -50,19 +49,33 @@ mpfr_sqrt (r, u, rnd_mode) char can_round = 0; TMP_DECL (marker); TMP_DECL(marker0); - if (MPFR_IS_NAN(u) || MPFR_SIGN(u) == -1) { + if (MPFR_IS_NAN(u)) { MPFR_SET_NAN(r); return 1; } - if (MPFR_SIGN(r) != 1) { MPFR_CHANGE_SIGN(r); } + if (MPFR_SIGN(u) < 0) { + if (MPFR_IS_INF(u) || MPFR_NOTZERO(u)) { + MPFR_SET_NAN(r); + return 1; + } + else { /* sqrt(-0) = -0 */ + MPFR_SET_ZERO(r); + if (MPFR_SIGN(r) > 0) MPFR_CHANGE_SIGN(r); + return 0; + } + } + + MPFR_CLEAR_NAN(r); + + if (MPFR_SIGN(r) < 0) MPFR_CHANGE_SIGN(r); if (MPFR_IS_INF(u)) { MPFR_SET_INF(r); return 1; } - MPFR_CLEAR_FLAGS (r); + MPFR_CLEAR_INF(r); prec = MPFR_PREC(r); @@ -19,7 +19,6 @@ along with the MPFR Library; see the file COPYING.LIB. If not, write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ -#include <math.h> #include <stdio.h> #include <stdlib.h> #include "gmp.h" @@ -525,6 +525,8 @@ mpfr_sub (a, b, c, rnd_mode) MPFR_SET_NAN(a); return; } + + MPFR_CLEAR_NAN(a); if (MPFR_IS_INF(b)) { @@ -543,18 +545,20 @@ mpfr_sub (a, b, c, rnd_mode) MPFR_SET_INF(a); if (MPFR_SIGN(b) != MPFR_SIGN(a)) { MPFR_CHANGE_SIGN(a); } } + return; } else if (MPFR_IS_INF(c)) { MPFR_SET_INF(a); if (MPFR_SIGN(c) == MPFR_SIGN(a)) { MPFR_CHANGE_SIGN(a); } + return; } if (!MPFR_NOTZERO(b)) { mpfr_neg(a, c, rnd_mode); return; } if (!MPFR_NOTZERO(c)) { mpfr_set(a, b, rnd_mode); return; } - MPFR_CLEAR_FLAGS (a); + MPFR_CLEAR_INF(a); diff_exp = MPFR_EXP(b)-MPFR_EXP(c); if (MPFR_SIGN(b) == MPFR_SIGN(c)) { @@ -41,17 +41,6 @@ mpfr_sub_ui (y, x, u, rnd_mode) mp_limb_t up[1]; unsigned long cnt; - if (MPFR_IS_NAN(x)) { - MPFR_SET_NAN(y); - return; - } - - if (MPFR_IS_INF(x)) - { - MPFR_SET_INF(y); - if (MPFR_SIGN(x) != MPFR_SIGN(y)) { MPFR_CHANGE_SIGN(y); } return; - } - if (u) { /* if u=0, do nothing */ MPFR_INIT1(up, uu, BITS_PER_MP_LIMB, 1); count_leading_zeros(cnt, (mp_limb_t) u); @@ -60,4 +49,6 @@ mpfr_sub_ui (y, x, u, rnd_mode) mpfr_sub (y, x, uu, rnd_mode); } + else + mpfr_set (y, x, rnd_mode); } diff --git a/tests/tcmp.c b/tests/tcmp.c index 46c06e8b3..3d033314d 100644 --- a/tests/tcmp.c +++ b/tests/tcmp.c @@ -23,7 +23,6 @@ MA 02111-1307, USA. */ #include <stdlib.h> #include <math.h> #include "gmp.h" -#include "longlong.h" #include "mpfr.h" #include "mpfr-test.h" @@ -101,6 +100,22 @@ int main() exit(1); } + mpfr_set_d(xx, Infp, GMP_RNDN); + mpfr_set_d(yy, Infp, GMP_RNDN); + if (mpfr_cmp(xx, yy) != 0) { + fprintf(stderr, + "Error in mpfr_cmp(Infp, Infp), gives %d\n", mpfr_cmp(xx, yy)); + exit(1); + } + + mpfr_set_d(xx, Infm, GMP_RNDN); + mpfr_set_d(yy, Infm, GMP_RNDN); + if (mpfr_cmp(xx, yy) != 0) { + fprintf(stderr, + "Error in mpfr_cmp(Infm, Infm), gives %d\n", mpfr_cmp(xx, yy)); + exit(1); + } + mpfr_set_d(xx, Infm, GMP_RNDN); mpfr_set_d(yy, 2346.09234, GMP_RNDN); if (mpfr_cmp(xx, yy) >= 0) { diff --git a/tests/tcmp2.c b/tests/tcmp2.c index 316ab1cf4..57ebe1a38 100644 --- a/tests/tcmp2.c +++ b/tests/tcmp2.c @@ -1,6 +1,6 @@ /* Test file for mpfr_cmp2. -Copyright (C) 1999 Free Software Foundation. +Copyright (C) 1999-2000 Free Software Foundation. This file is part of the MPFR Library. @@ -23,18 +23,15 @@ MA 02111-1307, USA. */ #include <stdlib.h> #include <math.h> #include "gmp.h" -#include "longlong.h" #include "mpfr.h" #include "mpfr-impl.h" -#ifdef __mips -#include <sys/fpu.h> -#endif -extern int isnan(); +void tcmp2 (double, double, int); -void tcmp2(x, y, i) double x, y; int i; +void tcmp2 (double x, double y, int i) { - mpfr_t xx,yy; int j; + mpfr_t xx, yy; + int j; if (i==-1) { if (x==y) i=53; @@ -79,7 +76,7 @@ int main() tcmp2(1.0, 1.0-x, i); x /= 2.0; } - for (j=0;j<1000000;j++) { + for (j=0; j<100000; j++) { x = drand48(); y = drand48(); if (x<y) { z=x; x=y; y=z; } diff --git a/tests/tcmp_ui.c b/tests/tcmp_ui.c index b1b8de517..ce0eb2ca5 100644 --- a/tests/tcmp_ui.c +++ b/tests/tcmp_ui.c @@ -23,7 +23,6 @@ MA 02111-1307, USA. */ #include <stdlib.h> #include <math.h> #include "gmp.h" -#include "longlong.h" #include "mpfr.h" int main() diff --git a/tests/tdump.c b/tests/tdump.c index 42ac7c4b5..17b000496 100644 --- a/tests/tdump.c +++ b/tests/tdump.c @@ -31,6 +31,7 @@ int main() mpfr_init2(z, 100); mpfr_set_ui(z, 0, GMP_RNDN); mpfr_dump(z, GMP_RNDD); + printf(" ^--- 0.e1 printed above is ok\n"); mpfr_clear(z); return 0; } diff --git a/tests/teq.c b/tests/teq.c index ce862653e..b5a3ff6f0 100644 --- a/tests/teq.c +++ b/tests/teq.c @@ -1,6 +1,6 @@ -/* Test file for mpfr_cmp2. +/* Test file for mpfr_eq. -Copyright (C) 1999 Free Software Foundation. +Copyright (C) 1999-2000 Free Software Foundation. This file is part of the MPFR Library. @@ -23,19 +23,19 @@ MA 02111-1307, USA. */ #include <stdlib.h> #include <math.h> #include "gmp.h" -#include "gmp-impl.h" -#include "longlong.h" #include "mpfr.h" #include "mpfr-impl.h" -void teq(mpfr_srcptr x) +void teq (mpfr_t); + +void teq (mpfr_t x) { mpfr_t y; long k, px, mx; mpfr_init2(y, MPFR_PREC(x)); - mx = (MPFR_PREC(x) - 1)/BITS_PER_MP_LIMB; - px = BITS_PER_MP_LIMB - 2; + mx = (MPFR_PREC(x) - 1)/mp_bits_per_limb; + px = mp_bits_per_limb - 2; for (k = 2; k < MPFR_PREC(x); k++) { @@ -55,7 +55,7 @@ void teq(mpfr_srcptr x) exit(-1); } - if (px) { --px; } else { --mx; px = BITS_PER_MP_LIMB - 1; } + if (px) { --px; } else { --mx; px = mp_bits_per_limb - 1; } } mpfr_clear(y); } @@ -64,11 +64,11 @@ int main() { int j; mpfr_t x; - mpfr_init2(x, 1000); + mpfr_init2 (x, 1000); for (j=0;j<1000;j++) { - mpfr_random(x); - teq(x); + mpfr_random (x); + teq (x); } mpfr_clear (x); return 0; diff --git a/tests/tmul_ui.c b/tests/tmul_ui.c index 4a65aab36..e226b35ba 100644 --- a/tests/tmul_ui.c +++ b/tests/tmul_ui.c @@ -22,7 +22,6 @@ MA 02111-1307, USA. */ #include <stdio.h> #include <stdlib.h> #include "gmp.h" -#include "gmp-impl.h" #include "mpfr.h" #include "mpfr-impl.h" @@ -36,7 +35,7 @@ main(int argc, char **argv) /* checks that result is normalized */ mpfr_set_d(y, 6.93147180559945286227e-01, GMP_RNDZ); mpfr_mul_ui(x, y, 1, GMP_RNDZ); - if (MPFR_MANT(x)[MPFR_PREC(x)/BITS_PER_MP_LIMB] >> (BITS_PER_MP_LIMB-1) == 0) { + if (MPFR_MANT(x)[MPFR_PREC(x)/mp_bits_per_limb] >> (mp_bits_per_limb-1) == 0) { fprintf(stderr, "Error in mpfr_mul_ui: result not normalized\n"); exit(1); } diff --git a/tests/trandom.c b/tests/trandom.c index 1233c45b5..7229bd46f 100644 --- a/tests/trandom.c +++ b/tests/trandom.c @@ -1,6 +1,6 @@ /* Test file for the various mpfr_random fonctions. -Copyright (C) 1999 Free Software Foundation. +Copyright (C) 1999-2000 Free Software Foundation. This file is part of the MPFR Library. @@ -23,12 +23,14 @@ MA 02111-1307, USA. */ #include <stdlib.h> #include <math.h> #include "gmp.h" -#include "gmp-impl.h" -#include "longlong.h" #include "mpfr.h" #include "mpfr-impl.h" -void test_random(unsigned long nbtests, unsigned long prec, int verbose) +void test_random (unsigned long, unsigned long, int); +void test_random2 (unsigned long, unsigned long, int); +void test_urandomb (unsigned long, unsigned long, int); + +void test_random (unsigned long nbtests, unsigned long prec, int verbose) { mpfr_t x; int *tab, size_tab, k; @@ -73,7 +75,7 @@ void test_random(unsigned long nbtests, unsigned long prec, int verbose) return; } -void test_random2(unsigned long nbtests, unsigned long prec, int verbose) +void test_random2 (unsigned long nbtests, unsigned long prec, int verbose) { mpfr_t x; int *tab, size_tab, k; @@ -86,7 +88,7 @@ void test_random2(unsigned long nbtests, unsigned long prec, int verbose) for (k = 0; k < size_tab; ++k) tab[k] = 0; for (k = 0; k < nbtests; k++) { - mpfr_random2(x, ABSSIZE(x), 0); + mpfr_random2 (x, MPFR_ABSSIZE(x), 0); d = mpfr_get_d(x); av += d; var += d*d; if (d < 1) tab[(int)(size_tab * d)]++; @@ -115,7 +117,7 @@ void test_random2(unsigned long nbtests, unsigned long prec, int verbose) return; } -void test_urandomb(unsigned long nbtests, unsigned long prec, int verbose) +void test_urandomb (unsigned long nbtests, unsigned long prec, int verbose) { mpfr_t x; int *tab, size_tab, k; @@ -162,7 +164,7 @@ void test_urandomb(unsigned long nbtests, unsigned long prec, int verbose) return; } -int main(int argc, char **argv) +int main (int argc, char **argv) { unsigned long nbtests, prec; int verbose = 0; diff --git a/tests/ttrunc.c b/tests/ttrunc.c index d30ae287d..41987170d 100644 --- a/tests/ttrunc.c +++ b/tests/ttrunc.c @@ -1,6 +1,6 @@ -/* Test file for mpfr_cmp2. +/* Test file for mpfr_trunc, mpfr_ceil, mpfr_floor. -Copyright (C) 1999 Free Software Foundation. +Copyright (C) 1999-2000 Free Software Foundation. This file is part of the MPFR Library. @@ -21,10 +21,7 @@ MA 02111-1307, USA. */ #include <stdio.h> #include <stdlib.h> -#include <math.h> #include "gmp.h" -#include "gmp-impl.h" -#include "longlong.h" #include "mpfr.h" #include "mpfr-impl.h" @@ -88,12 +88,15 @@ FUNC_NAME (r, u) return; } + MPFR_CLEAR_NAN(r); + if (MPFR_IS_INF(u)) { MPFR_SET_INF(r); + if (MPFR_SIGN(r) != MPFR_SIGN(u)) MPFR_CHANGE_SIGN(r); return; } - MPFR_CLEAR_FLAGS(r); + MPFR_CLEAR_INF(r); if (!MPFR_NOTZERO(u)) { MPFR_SET_ZERO(r); @@ -47,15 +47,18 @@ mpfr_ui_div (y, u, x, rnd_mode) return; } - MPFR_CLEAR_FLAGS(y); + MPFR_CLEAR_NAN(y); if (MPFR_IS_INF(x)) { + MPFR_CLEAR_INF(y); MPFR_SET_ZERO(y); if (MPFR_SIGN(x) != MPFR_SIGN(y)) MPFR_CHANGE_SIGN(y); return; } + MPFR_CLEAR_INF(y); + if (u) { MPFR_INIT1(up, uu, BITS_PER_MP_LIMB, 1); count_leading_zeros(cnt, (mp_limb_t) u); @@ -47,6 +47,8 @@ mpfr_ui_sub (y, u, x, rnd_mode) return; } + MPFR_CLEAR_NAN(y); + if (MPFR_IS_INF(x)) { MPFR_SET_INF(y); |