Annotation of OpenXM_contrib/gmp/mpfr/tests/tpow3.c, Revision 1.1
1.1 ! ohara 1: /* Test file for mpfr_pow.
! 2:
! 3: Copyright 2001 Free Software Foundation.
! 4: Adapted from tarctan.c.
! 5:
! 6: This file is part of the MPFR Library.
! 7:
! 8: The MPFR Library is free software; you can redistribute it and/or modify
! 9: it under the terms of the GNU Lesser General Public License as published by
! 10: the Free Software Foundation; either version 2.1 of the License, or (at your
! 11: option) any later version.
! 12:
! 13: The MPFR Library is distributed in the hope that it will be useful, but
! 14: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
! 15: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
! 16: License for more details.
! 17:
! 18: You should have received a copy of the GNU Lesser General Public License
! 19: along with the MPFR Library; see the file COPYING.LIB. If not, write to
! 20: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
! 21: MA 02111-1307, USA. */
! 22:
! 23: #include <stdio.h>
! 24: #include <limits.h>
! 25: #include <stdlib.h>
! 26: #include "gmp.h"
! 27: #include "gmp-impl.h"
! 28: #include "mpfr.h"
! 29: #include "mpfr-impl.h"
! 30: #include "mpfr-test.h"
! 31:
! 32:
! 33: int
! 34: main (int argc, char *argv[])
! 35: {
! 36: mpfr_t x, y, z, ax;
! 37: long int iy;
! 38: mpfr_init (x);
! 39: mpfr_init (ax);
! 40: mpfr_init2 (y,sizeof(unsigned long int)*CHAR_BIT);
! 41: mpfr_init (z);
! 42:
! 43: MPFR_SET_NAN(x);
! 44: mpfr_random(y);
! 45: mpfr_pow (z, x,y, GMP_RNDN);
! 46: if(!MPFR_IS_NAN(z))
! 47: {
! 48: printf ("evaluation of function in x=NAN does not return NAN");
! 49: exit (1);
! 50: }
! 51:
! 52: MPFR_SET_NAN(y);
! 53: mpfr_random(x);
! 54: mpfr_pow (z, x,y, GMP_RNDN);
! 55: if(!MPFR_IS_NAN(z))
! 56: {
! 57: printf ("evaluation of function in y=NAN does not return NAN");
! 58: exit (1);
! 59: }
! 60:
! 61: MPFR_CLEAR_FLAGS(z);
! 62: MPFR_CLEAR_FLAGS(y);
! 63: MPFR_CLEAR_FLAGS(x);
! 64:
! 65: MPFR_SET_ZERO(y);
! 66: mpfr_random(x);
! 67: mpfr_pow (z, x,y, GMP_RNDN);
! 68: if(mpfr_cmp_ui(z,1)!=0 && !(MPFR_IS_NAN(x)))
! 69: {
! 70: printf ("evaluation of function in y=0 does not return 1\n");
! 71: printf ("x =");
! 72: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 73: printf ("\n y =");
! 74: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 75: printf ("\n result =");
! 76: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 77: exit (1);
! 78: }
! 79:
! 80: MPFR_CLEAR_FLAGS(z);
! 81: MPFR_CLEAR_FLAGS(y);
! 82: MPFR_CLEAR_FLAGS(x);
! 83:
! 84: MPFR_SET_INF(y);
! 85: if (MPFR_SIGN(y) < 0)
! 86: MPFR_CHANGE_SIGN(y);
! 87: mpfr_random(x);
! 88: mpfr_set_prec (ax, MPFR_PREC(x));
! 89: mpfr_abs(ax,x,GMP_RNDN);
! 90: mpfr_pow (z, x,y, GMP_RNDN);
! 91: if( !MPFR_IS_INF(z) && (mpfr_cmp_ui(ax,1) > 0) )
! 92: {
! 93: printf ("evaluation of function in y=INF (|x|>1) does not return INF");
! 94: exit (1);
! 95: }
! 96: if( !MPFR_IS_ZERO(z) && (mpfr_cmp_ui(ax,1) < 0) )
! 97: {
! 98: printf ("\nevaluation of function in y=INF (|x|<1) does not return 0");
! 99: printf ("\nx =");
! 100: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 101: printf ("\n y =");
! 102: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 103: printf ("\n result =");
! 104: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 105: putchar('\n');
! 106: exit (1);
! 107: }
! 108:
! 109:
! 110: MPFR_CLEAR_FLAGS(z);
! 111: MPFR_CLEAR_FLAGS(y);
! 112: MPFR_CLEAR_FLAGS(x);
! 113:
! 114: MPFR_SET_INF(y);
! 115: if (MPFR_SIGN(y) > 0)
! 116: MPFR_CHANGE_SIGN(y);
! 117: mpfr_random(x);
! 118: mpfr_set_prec (ax, MPFR_PREC(x));
! 119: mpfr_abs(ax,x,GMP_RNDN);
! 120: mpfr_pow (z, x,y, GMP_RNDN);
! 121: mpfr_pow (z, x,y, GMP_RNDN);
! 122: if( !MPFR_IS_INF(z) && (mpfr_cmp_ui(ax,1) < 0) )
! 123: {
! 124: printf ("evaluation of function in y=INF (for |x| <0) does not return INF");
! 125: exit (1);
! 126: }
! 127: if( !MPFR_IS_ZERO(z) && (mpfr_cmp_ui(ax,1) > 0) )
! 128: {
! 129: printf ("evaluation of function in y=INF (for |x| >0) does not return 0");
! 130: exit (1);
! 131: }
! 132:
! 133: MPFR_CLEAR_FLAGS(z);
! 134: MPFR_CLEAR_FLAGS(y);
! 135: MPFR_CLEAR_FLAGS(x);
! 136:
! 137: MPFR_SET_INF(x);
! 138: if (MPFR_SIGN(x) < 0)
! 139: MPFR_CHANGE_SIGN(x);
! 140: mpfr_random(y);
! 141: mpfr_pow (z, x,y, GMP_RNDN);
! 142: if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) > 0))
! 143: {
! 144: printf ("evaluation of function in INF does not return INF");
! 145: printf ("\nx =");
! 146: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 147: printf ("\n y =");
! 148: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 149: printf ("\n result =");
! 150: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 151: putchar('\n');
! 152: exit (1);
! 153: }
! 154: if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0))
! 155: {
! 156: printf ("evaluation of function in INF does not return INF");
! 157: printf ("\nx =");
! 158: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 159: printf ("\n y =");
! 160: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 161: printf ("\n result =");
! 162: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 163: putchar('\n');
! 164: exit (1);
! 165: }
! 166:
! 167:
! 168: MPFR_CLEAR_FLAGS(z);
! 169: MPFR_CLEAR_FLAGS(y);
! 170: MPFR_CLEAR_FLAGS(x);
! 171:
! 172: MPFR_SET_INF(x);
! 173: if (MPFR_SIGN(x) > 0)
! 174: MPFR_CHANGE_SIGN(x);
! 175: mpfr_random(y);
! 176: if (random() % 2)
! 177: mpfr_neg (y, y, GMP_RNDN);
! 178: mpfr_pow (z, x,y, GMP_RNDN);
! 179: if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) > 0) && (mpfr_isinteger(y)))
! 180: {
! 181: printf ("evaluation of function in x=-INF does not return INF");
! 182: printf ("\nx =");
! 183: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 184: printf ("\n y =");
! 185: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 186: printf ("\n result =");
! 187: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 188: putchar('\n');
! 189: if(mpfr_isinteger(y))
! 190: printf("y is an integer\n");
! 191: else
! 192: printf("y is not an integer\n");
! 193:
! 194: exit (1);
! 195: }
! 196: if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0) && (mpfr_isinteger(y)))
! 197: {
! 198: printf ("evaluation of function in x=-INF does not return 0");
! 199: printf ("\nx =");
! 200: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 201: printf ("\n y =");
! 202: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 203: printf ("\n result =");
! 204: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 205: putchar('\n');
! 206:
! 207: if(mpfr_isinteger(y))
! 208: printf("y is an integer\n");
! 209: else
! 210: printf("y is not an integer\n");
! 211:
! 212: exit (1);
! 213: }
! 214: MPFR_CLEAR_FLAGS(z);
! 215: MPFR_CLEAR_FLAGS(y);
! 216: MPFR_CLEAR_FLAGS(x);
! 217:
! 218: MPFR_SET_INF(x);
! 219: if (MPFR_SIGN(x) > 0)
! 220: MPFR_CHANGE_SIGN(x);
! 221:
! 222: iy=random();
! 223: mpfr_random(y);
! 224: if (random() % 2)
! 225: iy=-iy;
! 226: mpfr_set_d(y,iy,GMP_RNDN);
! 227: mpfr_pow (z, x,y, GMP_RNDN);
! 228: if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) > 0) && (mpfr_isinteger(y)))
! 229: {
! 230: printf ("evaluation of function in x=-INF does not return INF");
! 231: printf ("\nx =");
! 232: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 233: printf ("\n y =");
! 234: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 235: printf ("\n result =");
! 236: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 237: putchar('\n');
! 238: if(mpfr_isinteger(y))
! 239: printf("y is an integer\n");
! 240: else
! 241: printf("y is not an integer\n");
! 242:
! 243: exit (1);
! 244: }
! 245: if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0) && (mpfr_isinteger(y)))
! 246: {
! 247: printf ("evaluation of function in x=-INF does not return 0");
! 248: printf ("\nx =");
! 249: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 250: printf ("\n y =");
! 251: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 252: printf ("\n result =");
! 253: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 254: putchar('\n');
! 255:
! 256: if(mpfr_isinteger(y))
! 257: printf("y is an integer\n");
! 258: else
! 259: printf("y is not an integer\n");
! 260:
! 261: exit (1);
! 262: }
! 263: MPFR_CLEAR_FLAGS(z);
! 264: MPFR_CLEAR_FLAGS(y);
! 265: MPFR_CLEAR_FLAGS(x);
! 266:
! 267: mpfr_set_ui(x,1,GMP_RNDN);
! 268: MPFR_SET_INF(y);
! 269: mpfr_pow (z, x,y, GMP_RNDN);
! 270: if(!MPFR_IS_NAN(z))
! 271: {
! 272: printf ("evaluation of function in x=1, y=INF does not return NAN");
! 273: printf ("\nx =");
! 274: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 275: printf ("\n y =");
! 276: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 277: printf ("\n result =");
! 278: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 279: putchar('\n');
! 280:
! 281: exit (1);
! 282: }
! 283:
! 284: MPFR_CLEAR_FLAGS(z);
! 285: MPFR_CLEAR_FLAGS(y);
! 286: MPFR_CLEAR_FLAGS(x);
! 287:
! 288: MPFR_SET_ZERO(x);
! 289: mpfr_random(y);
! 290: if (random() % 2)
! 291: mpfr_neg (y, y, GMP_RNDN);
! 292:
! 293: mpfr_pow (z, x,y, GMP_RNDN);
! 294:
! 295: if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) < 0) && !(mpfr_isinteger(y)))
! 296: {
! 297: printf ("evaluation of function in y<0 does not return 0");
! 298: printf ("\nx =");
! 299: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 300: printf ("\n y =");
! 301: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 302: printf ("\n result =");
! 303: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 304: putchar('\n');
! 305:
! 306: exit (1);
! 307: }
! 308: if(!MPFR_IS_INF(z) && (MPFR_SIGN(y) < 0) && (mpfr_isinteger(y)))
! 309: {
! 310: printf ("evaluation of function in y<0 (y integer) does not return INF");
! 311: printf ("\nx =");
! 312: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 313: printf ("\n y =");
! 314: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 315: printf ("\n result =");
! 316: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 317: putchar('\n');
! 318: exit (1);
! 319: }
! 320: if(!MPFR_IS_ZERO(z) && (MPFR_SIGN(y) > 0) && (mpfr_isinteger(y)))
! 321: {
! 322: printf ("evaluation of function in y<0 (y integer) does not return 0");
! 323: printf ("\nx =");
! 324: mpfr_out_str (stdout, 10, MPFR_PREC(x), x, GMP_RNDN);
! 325: printf ("\n y =");
! 326: mpfr_out_str (stdout, 10, MPFR_PREC(y), y, GMP_RNDN);
! 327: printf ("\n result =");
! 328: mpfr_out_str (stdout, 10, MPFR_PREC(z), z, GMP_RNDN);
! 329: putchar('\n');
! 330: exit (1);
! 331: }
! 332:
! 333:
! 334: {
! 335: mp_prec_t prec, yprec;
! 336: mpfr_t t, s;
! 337: mp_rnd_t rnd;
! 338: int inexact, compare, compare2;
! 339: unsigned int n, err;
! 340:
! 341: int p0=2;
! 342: int p1=100;
! 343: int N=100;
! 344:
! 345: mpfr_init (s);
! 346: mpfr_init (t);
! 347:
! 348: /* generic test */
! 349: for (prec = p0; prec <= p1; prec++)
! 350: {
! 351: mpfr_set_prec (x, prec);
! 352: mpfr_set_prec (s, sizeof(unsigned long int)*CHAR_BIT);
! 353: mpfr_set_prec (z, prec);
! 354: mpfr_set_prec (t, prec);
! 355: yprec = prec + 10;
! 356:
! 357: for (n=0; n<N; n++)
! 358: {
! 359:
! 360: mpfr_random (x);
! 361:
! 362: mpfr_random (s);
! 363: if (random() % 2)
! 364: mpfr_neg (s, s, GMP_RNDN);
! 365: rnd = random () % 4;
! 366: mpfr_set_prec (y, yprec);
! 367: compare = mpfr_pow (y, x, s, rnd);
! 368: err = (rnd == GMP_RNDN) ? yprec + 1 : yprec;
! 369: if (mpfr_can_round (y, err, rnd, rnd, prec))
! 370: {
! 371: mpfr_set (t, y, rnd);
! 372: inexact = mpfr_pow (z,x, s, rnd);
! 373: if (mpfr_cmp (t, z))
! 374: {
! 375: printf ("results differ for x=");
! 376: mpfr_out_str (stdout, 2, prec, x, GMP_RNDN);
! 377: printf (" values of the exponential=");
! 378: mpfr_out_str (stdout, 2, prec, s, GMP_RNDN);
! 379: printf (" prec=%u rnd_mode=%s\n", (unsigned) prec,
! 380: mpfr_print_rnd_mode (rnd));
! 381: printf ("got ");
! 382: mpfr_out_str (stdout, 2, prec, z, GMP_RNDN);
! 383: putchar ('\n');
! 384: printf ("expected ");
! 385: mpfr_out_str (stdout, 2, prec, t, GMP_RNDN);
! 386: putchar ('\n');
! 387: printf ("approx ");
! 388: mpfr_print_binary (y);
! 389: putchar ('\n');
! 390: exit (1);
! 391: }
! 392: compare2 = mpfr_cmp (t, y);
! 393: /* if rounding to nearest, cannot know the sign of t - f(x)
! 394: because of composed rounding: y = o(f(x)) and t = o(y) */
! 395: if ((rnd != GMP_RNDN) && (compare * compare2 >= 0))
! 396: compare = compare + compare2;
! 397: else
! 398: compare = inexact; /* cannot determine sign(t-f(x)) */
! 399: if (((inexact == 0) && (compare != 0)) ||
! 400: ((inexact > 0) && (compare <= 0)) ||
! 401: ((inexact < 0) && (compare >= 0)))
! 402: {
! 403: fprintf (stderr, "Wrong inexact flag for rnd=%s: expected %d, got %d\n",
! 404: mpfr_print_rnd_mode (rnd), compare, inexact);
! 405: printf ("x="); mpfr_print_binary (x); putchar ('\n');
! 406: printf ("y="); mpfr_print_binary (y); putchar ('\n');
! 407: printf ("t="); mpfr_print_binary (t); putchar ('\n');
! 408: exit (1);
! 409: }
! 410: }
! 411: }
! 412: }
! 413:
! 414: mpfr_clear (s);
! 415: mpfr_clear (t);
! 416:
! 417: }
! 418: mpfr_clear (x);
! 419: mpfr_clear (y);
! 420: mpfr_clear (z);
! 421: mpfr_clear (ax);
! 422:
! 423: return 0;
! 424: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>