1 /* Test file for round away. 2 3 Copyright 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013 Free Software Foundation, Inc. 4 Contributed by the AriC and Caramel projects, INRIA. 5 6 This file is part of the GNU MPFR Library. 7 8 The GNU 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 3 of the License, or (at your 11 option) any later version. 12 13 The GNU 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 GNU MPFR Library; see the file COPYING.LESSER. If not, see 20 http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 21 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 22 23 #include <stdio.h> 24 #include <stdlib.h> 25 26 #include "mpfr-test.h" 27 28 #define DISP(s, t) {printf(s); mpfr_out_str(stdout, 2, 0, t, MPFR_RNDN); } 29 #define DISP2(s,t) {DISP(s,t); putchar('\n');} 30 31 #define SPECIAL_MAX 12 32 33 static void 34 set_special (mpfr_ptr x, unsigned int select) 35 { 36 MPFR_ASSERTN (select < SPECIAL_MAX); 37 switch (select) 38 { 39 case 0: 40 MPFR_SET_NAN (x); 41 break; 42 case 1: 43 MPFR_SET_INF (x); 44 MPFR_SET_POS (x); 45 break; 46 case 2: 47 MPFR_SET_INF (x); 48 MPFR_SET_NEG (x); 49 break; 50 case 3: 51 MPFR_SET_ZERO (x); 52 MPFR_SET_POS (x); 53 break; 54 case 4: 55 MPFR_SET_ZERO (x); 56 MPFR_SET_NEG (x); 57 break; 58 case 5: 59 mpfr_set_str_binary (x, "1"); 60 break; 61 case 6: 62 mpfr_set_str_binary (x, "-1"); 63 break; 64 case 7: 65 mpfr_set_str_binary (x, "1e-1"); 66 break; 67 case 8: 68 mpfr_set_str_binary (x, "1e+1"); 69 break; 70 case 9: 71 mpfr_const_pi (x, MPFR_RNDN); 72 break; 73 case 10: 74 mpfr_const_pi (x, MPFR_RNDN); 75 MPFR_SET_EXP (x, MPFR_GET_EXP (x)-1); 76 break; 77 default: 78 mpfr_urandomb (x, RANDS); 79 break; 80 } 81 } 82 /* same than mpfr_cmp, but returns 0 for both NaN's */ 83 static int 84 mpfr_compare (mpfr_srcptr a, mpfr_srcptr b) 85 { 86 return (MPFR_IS_NAN(a)) ? !MPFR_IS_NAN(b) : 87 (MPFR_IS_NAN(b) || mpfr_cmp(a, b)); 88 } 89 90 static void 91 test3 (int (*testfunc)(mpfr_ptr, mpfr_srcptr, mpfr_srcptr, mpfr_rnd_t), 92 const char *foo) 93 { 94 mpfr_t ref1, ref2, ref3; 95 mpfr_t res1; 96 mpfr_prec_t p1, p2, p3; 97 int i, inexa, inexd; 98 mpfr_rnd_t r; 99 100 p1 = (randlimb () % 200) + MPFR_PREC_MIN; 101 p2 = (randlimb () % 200) + MPFR_PREC_MIN; 102 p3 = (randlimb () % 200) + MPFR_PREC_MIN; 103 104 mpfr_init2 (ref1, p1); 105 mpfr_init2 (ref2, p2); 106 mpfr_init2 (ref3, p3); 107 mpfr_init2 (res1, p1); 108 109 /* for each variable, consider each of the following 6 possibilities: 110 NaN, +Infinity, -Infinity, +0, -0 or a random number */ 111 for (i = 0; i < SPECIAL_MAX * SPECIAL_MAX; i++) 112 { 113 set_special (ref2, i%SPECIAL_MAX); 114 set_special (ref3, i/SPECIAL_MAX); 115 116 inexa = testfunc (res1, ref2, ref3, MPFR_RNDA); 117 r = MPFR_SIGN(res1) > 0 ? MPFR_RNDU : MPFR_RNDD; 118 inexd = testfunc (ref1, ref2, ref3, r); 119 120 if (mpfr_compare (res1, ref1) || inexa != inexd) 121 { 122 printf ("Error with RNDA for %s with ", foo); 123 DISP("x=",ref2); DISP2(", y=",ref3); 124 printf ("inexa=%d inexd=%d\n", inexa, inexd); 125 printf ("expected "); mpfr_print_binary (ref1); puts (""); 126 printf ("got "); mpfr_print_binary (res1); puts (""); 127 exit (1); 128 } 129 } 130 131 mpfr_clear (ref1); 132 mpfr_clear (ref2); 133 mpfr_clear (ref3); 134 mpfr_clear (res1); 135 } 136 137 static void 138 test4 (int (*testfunc)(mpfr_ptr, mpfr_srcptr, mpfr_srcptr, mpfr_srcptr, 139 mpfr_rnd_t), const char *foo) 140 { 141 mpfr_t ref, op1, op2, op3; 142 mpfr_prec_t pout, p1, p2, p3; 143 mpfr_t res; 144 int i, j, k, inexa, inexd; 145 mpfr_rnd_t r; 146 147 pout = (randlimb () % 200) + MPFR_PREC_MIN; 148 p1 = (randlimb () % 200) + MPFR_PREC_MIN; 149 p2 = (randlimb () % 200) + MPFR_PREC_MIN; 150 p3 = (randlimb () % 200) + MPFR_PREC_MIN; 151 152 mpfr_init2 (ref, pout); 153 mpfr_init2 (res, pout); 154 mpfr_init2 (op1, p1); 155 mpfr_init2 (op2, p2); 156 mpfr_init2 (op3, p3); 157 158 /* for each variable, consider each of the following 6 possibilities: 159 NaN, +Infinity, -Infinity, +0, -0 or a random number */ 160 161 for (i = 0; i < SPECIAL_MAX; i++) 162 { 163 set_special (op1, i); 164 for (j = 0; j < SPECIAL_MAX; j++) 165 { 166 set_special (op2, j); 167 for (k = 0; k < SPECIAL_MAX; k++) 168 { 169 set_special (op3, k); 170 171 inexa = testfunc (res, op1, op2, op3, MPFR_RNDA); 172 r = MPFR_SIGN(res) > 0 ? MPFR_RNDU : MPFR_RNDD; 173 inexd = testfunc (ref, op1, op2, op3, r); 174 175 if (mpfr_compare (res, ref) || inexa != inexd) 176 { 177 printf ("Error with RNDA for %s with ", foo); 178 DISP("a=", op1); DISP(", b=", op2); DISP2(", c=", op3); 179 printf ("inexa=%d inexd=%d\n", inexa, inexd); 180 DISP("expected ", ref); DISP2(", got ", res); 181 exit (1); 182 } 183 } 184 } 185 } 186 187 mpfr_clear (ref); 188 mpfr_clear (op1); 189 mpfr_clear (op2); 190 mpfr_clear (op3); 191 mpfr_clear (res); 192 } 193 194 static void 195 test2ui (int (*testfunc)(mpfr_ptr, mpfr_srcptr, unsigned long int, mpfr_rnd_t), 196 const char *foo) 197 { 198 mpfr_t ref1, ref2; 199 unsigned int ref3; 200 mpfr_t res1; 201 mpfr_prec_t p1, p2; 202 int i, inexa, inexd; 203 mpfr_rnd_t r; 204 205 p1 = (randlimb () % 200) + MPFR_PREC_MIN; 206 p2 = (randlimb () % 200) + MPFR_PREC_MIN; 207 208 mpfr_init2 (ref1, p1); 209 mpfr_init2 (ref2, p2); 210 mpfr_init2 (res1, p1); 211 212 /* ref2 can be NaN, +Inf, -Inf, +0, -0 or any number 213 ref3 can be 0 or any number */ 214 for (i = 0; i < SPECIAL_MAX * 2; i++) 215 { 216 set_special (ref2, i % SPECIAL_MAX); 217 ref3 = i / SPECIAL_MAX == 0 ? 0 : randlimb (); 218 219 inexa = testfunc (res1, ref2, ref3, MPFR_RNDA); 220 r = MPFR_SIGN(res1) > 0 ? MPFR_RNDU : MPFR_RNDD; 221 inexd = testfunc (ref1, ref2, ref3, r); 222 223 if (mpfr_compare (res1, ref1) || inexa != inexd) 224 { 225 printf ("Error with RNDA for %s for c=%u\n", foo, ref3); 226 DISP2("a=",ref2); 227 printf ("inexa=%d inexd=%d\n", inexa, inexd); 228 printf ("expected "); mpfr_print_binary (ref1); puts (""); 229 printf ("got "); mpfr_print_binary (res1); puts (""); 230 exit (1); 231 } 232 } 233 234 mpfr_clear (ref1); 235 mpfr_clear (ref2); 236 mpfr_clear (res1); 237 } 238 239 static void 240 testui2 (int (*testfunc)(mpfr_ptr, unsigned long int, mpfr_srcptr, mpfr_rnd_t), 241 const char *foo) 242 { 243 mpfr_t ref1, ref3; 244 unsigned int ref2; 245 mpfr_t res1; 246 mpfr_prec_t p1, p3; 247 int i, inexa, inexd; 248 mpfr_rnd_t r; 249 250 p1 = (randlimb () % 200) + MPFR_PREC_MIN; 251 p3 = (randlimb () % 200) + MPFR_PREC_MIN; 252 253 mpfr_init2 (ref1, p1); 254 mpfr_init2 (ref3, p3); 255 mpfr_init2 (res1, p1); 256 257 for (i = 0; i < SPECIAL_MAX * 2; i++) 258 { 259 set_special (ref3, i % SPECIAL_MAX); 260 ref2 = i / SPECIAL_MAX == 0 ? 0 : randlimb (); 261 262 inexa = testfunc (res1, ref2, ref3, MPFR_RNDA); 263 r = MPFR_SIGN(res1) > 0 ? MPFR_RNDU : MPFR_RNDD; 264 inexd = testfunc (ref1, ref2, ref3, r); 265 266 if (mpfr_compare (res1, ref1) || inexa != inexd) 267 { 268 printf ("Error with RNDA for %s for b=%u\n", foo, ref2); 269 DISP2("a=", ref3); 270 printf ("inexa=%d inexd=%d\n", inexa, inexd); 271 DISP("expected ", ref1); DISP2(", got ", res1); 272 exit (1); 273 } 274 } 275 276 mpfr_clear (ref1); 277 mpfr_clear (ref3); 278 mpfr_clear (res1); 279 } 280 281 /* foo(mpfr_ptr, mpfr_srcptr, mp_rndt) */ 282 static void 283 test2 (int (*testfunc)(mpfr_ptr, mpfr_srcptr, mpfr_rnd_t), const char *foo) 284 { 285 mpfr_t ref1, ref2; 286 mpfr_t res1; 287 mpfr_prec_t p1, p2; 288 int i, inexa, inexd; 289 mpfr_rnd_t r; 290 291 p1 = (randlimb () % 200) + MPFR_PREC_MIN; 292 p2 = (randlimb () % 200) + MPFR_PREC_MIN; 293 294 mpfr_init2 (ref1, p1); 295 mpfr_init2 (ref2, p2); 296 mpfr_init2 (res1, p1); 297 298 for (i = 0; i < SPECIAL_MAX; i++) 299 { 300 set_special (ref2, i); 301 302 /* first round to away */ 303 inexa = testfunc (res1, ref2, MPFR_RNDA); 304 305 r = MPFR_SIGN(res1) > 0 ? MPFR_RNDU : MPFR_RNDD; 306 inexd = testfunc (ref1, ref2, r); 307 if (mpfr_compare (res1, ref1) || inexa != inexd) 308 { 309 printf ("Error with RNDA for %s with ", foo); 310 DISP2("x=", ref2); 311 printf ("inexa=%d inexd=%d\n", inexa, inexd); 312 DISP("expected ", ref1); DISP2(", got ", res1); 313 exit (1); 314 } 315 } 316 317 mpfr_clear (ref1); 318 mpfr_clear (ref2); 319 mpfr_clear (res1); 320 } 321 322 /* one operand, two results, like mpfr_sin_cos */ 323 static void 324 test3a (int (*testfunc)(mpfr_ptr, mpfr_ptr, mpfr_srcptr, mpfr_rnd_t), 325 const char *foo) 326 { 327 mpfr_t ref1, ref2, ref3; 328 mpfr_t res1, res2; 329 mpfr_prec_t p1, p2, p3; 330 int i, inexa, inexd; 331 mpfr_rnd_t r; 332 333 p1 = (randlimb () % 200) + MPFR_PREC_MIN; 334 p2 = (randlimb () % 200) + MPFR_PREC_MIN; 335 p3 = (randlimb () % 200) + MPFR_PREC_MIN; 336 337 mpfr_init2 (ref1, p1); 338 mpfr_init2 (ref2, p2); 339 mpfr_init2 (ref3, p3); 340 mpfr_init2 (res1, p1); 341 mpfr_init2 (res2, p2); 342 343 for (i = 0; i < SPECIAL_MAX; i++) 344 { 345 set_special (ref3, i); 346 347 inexa = testfunc (res1, res2, ref3, MPFR_RNDA); 348 349 /* first check wrt the first operand */ 350 r = MPFR_SIGN(res1) > 0 ? MPFR_RNDU : MPFR_RNDD; 351 inexd = testfunc (ref1, ref2, ref3, r); 352 /* the low 2 bits of the inexact flag concern the 1st operand */ 353 if (mpfr_compare (res1, ref1) || (inexa & 3) != (inexd & 3)) 354 { 355 printf ("Error with RNDA for %s (1st operand)\n", foo); 356 DISP2("a=",ref3); 357 DISP("expected ", ref1); printf ("\n"); 358 DISP("got ", res1); printf ("\n"); 359 printf ("inexa=%d inexd=%d\n", inexa & 3, inexd & 3); 360 exit (1); 361 } 362 363 /* now check wrt the second operand */ 364 r = MPFR_SIGN(res2) > 0 ? MPFR_RNDU : MPFR_RNDD; 365 inexd = testfunc (ref1, ref2, ref3, r); 366 /* bits 2..3 of the inexact flag concern the 2nd operand */ 367 if (mpfr_compare (res2, ref2) || (inexa >> 2) != (inexd >> 2)) 368 { 369 printf ("Error with RNDA for %s (2nd operand)\n", foo); 370 DISP2("a=",ref3); 371 DISP("expected ", ref2); printf ("\n"); 372 DISP("got ", res2); printf ("\n"); 373 printf ("inexa=%d inexd=%d\n", inexa >> 2, inexd >> 2); 374 exit (1); 375 } 376 377 } 378 379 mpfr_clear (ref1); 380 mpfr_clear (ref2); 381 mpfr_clear (ref3); 382 mpfr_clear (res1); 383 mpfr_clear (res2); 384 } 385 386 int 387 main (void) 388 { 389 int N = 20; 390 391 tests_start_mpfr (); 392 393 while (N--) 394 { 395 /* no need to test mpfr_round, mpfr_ceil, mpfr_floor, mpfr_trunc since 396 they take no rounding mode */ 397 398 test2ui (mpfr_add_ui, "mpfr_add_ui"); 399 test2ui (mpfr_div_2exp, "mpfr_div_2exp"); 400 test2ui (mpfr_div_ui, "mpfr_div_ui"); 401 test2ui (mpfr_mul_2exp, "mpfr_mul_2exp"); 402 test2ui (mpfr_mul_ui, "mpfr_mul_ui"); 403 test2ui (mpfr_pow_ui, "mpfr_pow_ui"); 404 test2ui (mpfr_sub_ui, "mpfr_sub_ui"); 405 406 testui2 (mpfr_ui_div, "mpfr_ui_div"); 407 testui2 (mpfr_ui_sub, "mpfr_ui_sub"); 408 testui2 (mpfr_ui_pow, "mpfr_ui_pow"); 409 410 test2 (mpfr_sqr, "mpfr_sqr"); 411 test2 (mpfr_sqrt, "mpfr_sqrt"); 412 test2 (mpfr_abs, "mpfr_abs"); 413 test2 (mpfr_neg, "mpfr_neg"); 414 415 test2 (mpfr_log, "mpfr_log"); 416 test2 (mpfr_log2, "mpfr_log2"); 417 test2 (mpfr_log10, "mpfr_log10"); 418 test2 (mpfr_log1p, "mpfr_log1p"); 419 420 test2 (mpfr_exp, "mpfr_exp"); 421 test2 (mpfr_exp2, "mpfr_exp2"); 422 test2 (mpfr_exp10, "mpfr_exp10"); 423 test2 (mpfr_expm1, "mpfr_expm1"); 424 test2 (mpfr_eint, "mpfr_eint"); 425 426 test2 (mpfr_sinh, "mpfr_sinh"); 427 test2 (mpfr_cosh, "mpfr_cosh"); 428 test2 (mpfr_tanh, "mpfr_tanh"); 429 test2 (mpfr_asinh, "mpfr_asinh"); 430 test2 (mpfr_acosh, "mpfr_acosh"); 431 test2 (mpfr_atanh, "mpfr_atanh"); 432 test2 (mpfr_sech, "mpfr_sech"); 433 test2 (mpfr_csch, "mpfr_csch"); 434 test2 (mpfr_coth, "mpfr_coth"); 435 436 test2 (mpfr_asin, "mpfr_asin"); 437 test2 (mpfr_acos, "mpfr_acos"); 438 test2 (mpfr_atan, "mpfr_atan"); 439 test2 (mpfr_cos, "mpfr_cos"); 440 test2 (mpfr_sin, "mpfr_sin"); 441 test2 (mpfr_tan, "mpfr_tan"); 442 test2 (mpfr_sec, "mpfr_sec"); 443 test2 (mpfr_csc, "mpfr_csc"); 444 test2 (mpfr_cot, "mpfr_cot"); 445 446 test2 (mpfr_erf, "mpfr_erf"); 447 test2 (mpfr_erfc, "mpfr_erfc"); 448 test2 (mpfr_j0, "mpfr_j0"); 449 test2 (mpfr_j1, "mpfr_j1"); 450 test2 (mpfr_y0, "mpfr_y0"); 451 test2 (mpfr_y1, "mpfr_y1"); 452 test2 (mpfr_zeta, "mpfr_zeta"); 453 test2 (mpfr_gamma, "mpfr_gamma"); 454 test2 (mpfr_lngamma, "mpfr_lngamma"); 455 456 test2 (mpfr_rint, "mpfr_rint"); 457 test2 (mpfr_rint_ceil, "mpfr_rint_ceil"); 458 test2 (mpfr_rint_floor, "mpfr_rint_floor"); 459 test2 (mpfr_rint_round, "mpfr_rint_round"); 460 test2 (mpfr_rint_trunc, "mpfr_rint_trunc"); 461 test2 (mpfr_frac, "mpfr_frac"); 462 463 test3 (mpfr_add, "mpfr_add"); 464 test3 (mpfr_sub, "mpfr_sub"); 465 test3 (mpfr_mul, "mpfr_mul"); 466 test3 (mpfr_div, "mpfr_div"); 467 468 test3 (mpfr_agm, "mpfr_agm"); 469 test3 (mpfr_min, "mpfr_min"); 470 test3 (mpfr_max, "mpfr_max"); 471 472 /* we don't test reldiff since it does not guarantee correct rounding, 473 thus we can get different results with RNDA and RNDU or RNDD. */ 474 test3 (mpfr_dim, "mpfr_dim"); 475 476 test3 (mpfr_remainder, "mpfr_remainder"); 477 test3 (mpfr_pow, "mpfr_pow"); 478 test3 (mpfr_atan2, "mpfr_atan2"); 479 test3 (mpfr_hypot, "mpfr_hypot"); 480 481 test3a (mpfr_sin_cos, "mpfr_sin_cos"); 482 483 test4 (mpfr_fma, "mpfr_fma"); 484 test4 (mpfr_fms, "mpfr_fms"); 485 486 #if MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0) 487 test2 (mpfr_li2, "mpfr_li2"); 488 test2 (mpfr_rec_sqrt, "mpfr_rec_sqrt"); 489 test3 (mpfr_fmod, "mpfr_fmod"); 490 test3a (mpfr_modf, "mpfr_modf"); 491 test3a (mpfr_sinh_cosh, "mpfr_sinh_cosh"); 492 #endif 493 } 494 495 tests_end_mpfr (); 496 return 0; 497 } 498