Annotation of OpenXM_contrib/gmp/mpfr/tests/tcmp2.c, Revision 1.1.1.2
1.1 maekawa 1: /* Test file for mpfr_cmp2.
2:
1.1.1.2 ! ohara 3: Copyright 1999, 2000, 2001, 2002 Free Software Foundation.
1.1 maekawa 4:
5: This file is part of the MPFR Library.
6:
7: The MPFR Library is free software; you can redistribute it and/or modify
1.1.1.2 ! ohara 8: it under the terms of the GNU Lesser General Public License as published by
! 9: the Free Software Foundation; either version 2.1 of the License, or (at your
1.1 maekawa 10: option) any later version.
11:
12: The MPFR Library is distributed in the hope that it will be useful, but
13: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.1.1.2 ! ohara 14: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
1.1 maekawa 15: License for more details.
16:
1.1.1.2 ! ohara 17: You should have received a copy of the GNU Lesser General Public License
1.1 maekawa 18: along with the MPFR Library; see the file COPYING.LIB. If not, write to
19: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
20: MA 02111-1307, USA. */
21:
22: #include <stdio.h>
23: #include <stdlib.h>
24: #include <math.h>
25: #include "gmp.h"
26: #include "mpfr.h"
27: #include "mpfr-impl.h"
1.1.1.2 ! ohara 28: #include "mpfr-test.h"
1.1 maekawa 29:
1.1.1.2 ! ohara 30: void tcmp2 _PROTO ((double, double, int));
! 31: void special _PROTO ((void));
! 32: void worst_cases _PROTO ((void));
! 33: void set_bit _PROTO ((mpfr_t, unsigned int, int));
! 34:
! 35: /* set bit n of x to b, where bit 0 is the most significant one */
! 36: void
! 37: set_bit (mpfr_t x, unsigned int n, int b)
! 38: {
! 39: unsigned l;
! 40: mp_size_t xn;
! 41:
! 42: xn = (MPFR_PREC(x) - 1) / mp_bits_per_limb;
! 43: l = n / mp_bits_per_limb;
! 44: n %= mp_bits_per_limb;
! 45: n = mp_bits_per_limb - 1 - n;
! 46: if (b)
! 47: MPFR_MANT(x)[xn - l] |= (mp_limb_t) 1 << n;
! 48: else
! 49: MPFR_MANT(x)[xn - l] &= ~((mp_limb_t) 1 << n);
! 50: }
! 51:
! 52: /* check that for x = 1.u 1 v 0^k low(x)
! 53: y = 1.u 0 v 1^k low(y)
! 54: mpfr_cmp2 (x, y) returns 1 + |u| + |v| + k for low(x) >= low(y),
! 55: and 1 + |u| + |v| + k + 1 otherwise */
! 56: void
! 57: worst_cases (void)
! 58: {
! 59: mpfr_t x, y;
! 60: unsigned int i, j, k, b, expected;
! 61: mp_prec_t l;
! 62:
! 63: mpfr_init2 (x, 200);
! 64: mpfr_init2 (y, 200);
! 65:
! 66: mpfr_set_ui (y, 1, GMP_RNDN);
! 67: for (i=1; i<MPFR_PREC(x); i++)
! 68: {
! 69: mpfr_set_ui (x, 1, GMP_RNDN);
! 70: mpfr_div_2exp (y, y, 1, GMP_RNDN); /* y = 1/2^i */
! 71:
! 72: l = 0;
! 73: if (mpfr_cmp2 (x, y, &l) <= 0 || l != 1)
! 74: {
! 75: fprintf (stderr, "Error in mpfr_cmp2:\nx=");
! 76: mpfr_out_str (stderr, 2, 0, x, GMP_RNDN);
! 77: fprintf (stderr, "\ny=");
! 78: mpfr_out_str (stderr, 2, 0, y, GMP_RNDN);
! 79: fprintf (stderr, "\ngot %lu instead of %u\n", l, 1);
! 80: exit(1);
! 81: }
! 82:
! 83: mpfr_add (x, x, y, GMP_RNDN); /* x = 1 + 1/2^i */
! 84: l = 0;
! 85: if (mpfr_cmp2 (x, y, &l) <= 0 || l != 0)
! 86: {
! 87: fprintf (stderr, "Error in mpfr_cmp2:\nx=");
! 88: mpfr_out_str (stderr, 2, 0, x, GMP_RNDN);
! 89: fprintf (stderr, "\ny=");
! 90: mpfr_out_str (stderr, 2, 0, y, GMP_RNDN);
! 91: fprintf (stderr, "\ngot %lu instead of %u\n", l, 0);
! 92: exit(1);
! 93: }
! 94: }
! 95:
! 96: for (i=0; i<64; i++) /* |u| = i */
! 97: {
! 98: mpfr_random (x);
! 99: mpfr_set (y, x, GMP_RNDN);
! 100: set_bit (x, i + 1, 1);
! 101: set_bit (y, i + 1, 0);
! 102: for (j=0; j<64; j++) /* |v| = j */
! 103: {
! 104: b = random () % 2;
! 105: set_bit (x, i + j + 2, b);
! 106: set_bit (y, i + j + 2, b);
! 107:
! 108: for (k=0; k<64; k++)
! 109: {
! 110:
! 111: if (k) set_bit (x, i + j + k + 1, 0);
! 112: if (k) set_bit (y, i + j + k + 1, 1);
! 113:
! 114: set_bit (x, i + j + k + 2, 1);
! 115: set_bit (y, i + j + k + 2, 0);
! 116: l = 0; mpfr_cmp2 (x, y, &l);
! 117: expected = i + j + k + 1;
! 118: if (l != expected)
! 119: {
! 120: fprintf (stderr, "Error in mpfr_cmp2:\nx=");
! 121: mpfr_out_str (stderr, 2, 0, x, GMP_RNDN);
! 122: fprintf (stderr, "\ny=");
! 123: mpfr_out_str (stderr, 2, 0, y, GMP_RNDN);
! 124: fprintf (stderr, "\ngot %lu instead of %u\n", l, expected);
! 125: exit(1);
! 126: }
! 127:
! 128: set_bit (x, i + j + k + 2, 0);
! 129: set_bit (x, i + j + k + 3, 0);
! 130: set_bit (y, i + j + k + 3, 1);
! 131: l = 0; mpfr_cmp2 (x, y, &l);
! 132: expected = i + j + k + 2;
! 133: if (l != expected)
! 134: {
! 135: fprintf (stderr, "Error in mpfr_cmp2:\nx=");
! 136: mpfr_out_str (stderr, 2, 0, x, GMP_RNDN);
! 137: fprintf (stderr, "\ny=");
! 138: mpfr_out_str (stderr, 2, 0, y, GMP_RNDN);
! 139: fprintf (stderr, "\ngot %lu instead of %u\n", l, expected);
! 140: exit(1);
! 141: }
! 142: }
! 143: }
! 144: }
1.1 maekawa 145:
1.1.1.2 ! ohara 146: mpfr_clear (x);
! 147: mpfr_clear (y);
! 148: }
! 149:
! 150: void
! 151: tcmp2 (double x, double y, int i)
1.1 maekawa 152: {
1.1.1.2 ! ohara 153: mpfr_t xx, yy;
! 154: mp_prec_t j;
1.1 maekawa 155:
156: if (i==-1) {
157: if (x==y) i=53;
158: else i = (int) floor(log(x)/log(2.0)) - (int) floor(log(x-y)/log(2.0));
159: }
160: mpfr_init2(xx, 53); mpfr_init2(yy, 53);
1.1.1.2 ! ohara 161: mpfr_set_d (xx, x, GMP_RNDN);
! 162: mpfr_set_d (yy, y, GMP_RNDN);
! 163: j = 0;
! 164: if (mpfr_cmp2 (xx, yy, &j) == 0)
! 165: {
! 166: if (x != y)
! 167: {
! 168: fprintf (stderr, "Error in mpfr_cmp2 for\nx=");
! 169: mpfr_out_str (stderr, 2, 0, xx, GMP_RNDN);
! 170: fprintf (stderr, "\ny=");
! 171: mpfr_out_str (stderr, 2, 0, yy, GMP_RNDN);
! 172: fprintf (stderr, "\ngot sign 0 for x != y\n");
! 173: exit(1);
! 174: }
! 175: }
! 176: else if (j != i) {
! 177: fprintf (stderr, "Error in mpfr_cmp2 for\nx=");
! 178: mpfr_out_str (stderr, 2, 0, xx, GMP_RNDN);
! 179: fprintf (stderr, "\ny=");
! 180: mpfr_out_str (stderr, 2, 0, yy, GMP_RNDN);
! 181: fprintf (stderr, "\ngot %lu instead of %u\n", j, i);
1.1 maekawa 182: exit(1);
183: }
1.1.1.2 ! ohara 184: mpfr_clear(xx); mpfr_clear(yy);
! 185: }
! 186:
! 187: void
! 188: special (void)
! 189: {
! 190: mpfr_t x, y;
! 191: mp_prec_t j;
! 192:
! 193: mpfr_init (x); mpfr_init (y);
! 194:
! 195: /* bug found by Nathalie Revol, 21 March 2001 */
! 196: mpfr_set_prec (x, 65);
! 197: mpfr_set_prec (y, 65);
! 198: mpfr_set_str_raw (x, "0.10000000000000000000000000000000000001110010010110100110011110000E1");
! 199: mpfr_set_str_raw (y, "0.11100100101101001100111011111111110001101001000011101001001010010E-35");
! 200: j = 0;
! 201: if (mpfr_cmp2 (x, y, &j) <= 0 || j != 1) {
! 202: printf ("Error in mpfr_cmp2:\n");
! 203: printf ("x=");
! 204: mpfr_print_binary (x);
! 205: putchar ('\n');
! 206: printf ("y=");
! 207: mpfr_print_binary (y);
! 208: putchar ('\n');
! 209: printf ("got %lu, expected 1\n", j);
! 210: exit (1);
! 211: }
! 212:
! 213: mpfr_set_prec(x, 127); mpfr_set_prec(y, 127);
! 214: mpfr_set_str_raw(x, "0.1011010000110111111000000101011110110001000101101011011110010010011110010000101101000010011001100110010000000010110000101000101E6");
! 215: mpfr_set_str_raw(y, "0.1011010000110111111000000101011011111100011101000011001111000010100010100110110100110010011001100110010000110010010110000010110E6");
! 216: j = 0;
! 217: if (mpfr_cmp2(x, y, &j) <= 0 || j != 32) {
1.1 maekawa 218: printf("Error in mpfr_cmp2:\n");
1.1.1.2 ! ohara 219: printf("x="); mpfr_print_binary(x); putchar('\n');
! 220: printf("y="); mpfr_print_binary(y); putchar('\n');
! 221: printf("got %lu, expected 32\n", j);
1.1 maekawa 222: exit(1);
223: }
1.1.1.2 ! ohara 224:
! 225: mpfr_set_prec (x, 128); mpfr_set_prec (y, 239);
! 226: mpfr_set_str_raw (x, "0.10001000110110000111011000101011111100110010010011001101000011111010010110001000000010100110100111111011011010101100100000000000E167");
! 227: mpfr_set_str_raw (y, "0.10001000110110000111011000101011111100110010010011001101000011111010010110001000000010100110100111111011011010101100011111111111111111111111111111111111111111111111011111100101011100011001101000100111000010000000000101100110000111111000101E167");
! 228: j = 0;
! 229: if (mpfr_cmp2(x, y, &j) <= 0 || j != 164) {
! 230: printf("Error in mpfr_cmp2:\n");
! 231: printf("x="); mpfr_print_binary(x); putchar('\n');
! 232: printf("y="); mpfr_print_binary(y); putchar('\n');
! 233: printf("got %lu, expected 164\n", j);
! 234: exit(1);
! 235: }
! 236:
! 237: /* bug found by Nathalie Revol, 29 March 2001 */
! 238: mpfr_set_prec (x, 130); mpfr_set_prec (y, 130);
! 239: mpfr_set_str_raw (x, "0.1100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000E2");
! 240: mpfr_set_str_raw (y, "0.1011111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111111100E2");
! 241: j = 0;
! 242: if (mpfr_cmp2(x, y, &j) <= 0 || j != 127) {
! 243: printf("Error in mpfr_cmp2:\n");
! 244: printf("x="); mpfr_print_binary(x); putchar('\n');
! 245: printf("y="); mpfr_print_binary(y); putchar('\n');
! 246: printf("got %lu, expected 127\n", j);
! 247: exit(1);
! 248: }
! 249:
! 250: /* bug found by Nathalie Revol, 2 Apr 2001 */
! 251: mpfr_set_prec (x, 65); mpfr_set_prec (y, 65);
! 252: mpfr_set_ui (x, 5, GMP_RNDN);
! 253: mpfr_set_str_raw (y, "0.10011111111111111111111111111111111111111111111111111111111111101E3");
! 254: j = 0;
! 255: if (mpfr_cmp2(x, y, &j) <= 0 || j != 63) {
! 256: printf("Error in mpfr_cmp2:\n");
! 257: printf("x="); mpfr_print_binary(x); putchar('\n');
! 258: printf("y="); mpfr_print_binary(y); putchar('\n');
! 259: printf("got %lu, expected 63\n", j);
! 260: exit(1);
! 261: }
! 262:
! 263: /* bug found by Nathalie Revol, 2 Apr 2001 */
! 264: mpfr_set_prec (x, 65); mpfr_set_prec (y, 65);
! 265: mpfr_set_str_raw (x, "0.10011011111000101001110000000000000000000000000000000000000000000E-69");
! 266: mpfr_set_str_raw (y, "0.10011011111000101001101111111111111111111111111111111111111111101E-69");
! 267: j = 0;
! 268: if (mpfr_cmp2(x, y, &j) <= 0 || j != 63) {
! 269: printf("Error in mpfr_cmp2:\n");
! 270: printf("x="); mpfr_print_binary(x); putchar('\n');
! 271: printf("y="); mpfr_print_binary(y); putchar('\n');
! 272: printf("got %lu, expected 63\n", j);
! 273: exit(1);
! 274: }
! 275:
! 276: mpfr_clear(x); mpfr_clear(y);
1.1 maekawa 277: }
278:
1.1.1.2 ! ohara 279: int
! 280: main (void)
1.1 maekawa 281: {
282: int i,j; double x=1.0, y, z;
283:
1.1.1.2 ! ohara 284: mpfr_test_init ();
! 285:
! 286: worst_cases ();
! 287: special ();
! 288: tcmp2(5.43885304644369510000e+185, -1.87427265794105340000e-57, 1);
1.1 maekawa 289: tcmp2(1.06022698059744327881e+71, 1.05824655795525779205e+71, -1);
290: tcmp2(1.0, 1.0, 53);
291: for (i=0;i<54;i++) {
1.1.1.2 ! ohara 292: tcmp2 (1.0, 1.0-x, i);
1.1 maekawa 293: x /= 2.0;
294: }
1.1.1.2 ! ohara 295: for (x=0.5, i=1; i<100; i++) {
! 296: tcmp2 (1.0, x, 1);
! 297: x /= 2.0;
! 298: }
! 299: for (j=0; j<100000; j++) {
! 300: x = DBL_RAND ();
! 301: y = DBL_RAND ();
1.1 maekawa 302: if (x<y) { z=x; x=y; y=z; }
303: if (y != 0.0 && y != -0.0) tcmp2(x, y, -1);
304: }
305:
1.1.1.2 ! ohara 306: return 0;
! 307: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to tcmp2.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / mpfr / tests