tgeneric.c revision 1.1.1.1.2.1
1 /* Generic test file for functions with one or two arguments (the second being 2 either mpfr_t or double). 3 4 Copyright 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013 Free Software Foundation, Inc. 5 Contributed by the AriC and Caramel projects, INRIA. 6 7 This file is part of the GNU MPFR Library. 8 9 The GNU MPFR Library is free software; you can redistribute it and/or modify 10 it under the terms of the GNU Lesser General Public License as published by 11 the Free Software Foundation; either version 3 of the License, or (at your 12 option) any later version. 13 14 The GNU MPFR Library is distributed in the hope that it will be useful, but 15 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 16 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 17 License for more details. 18 19 You should have received a copy of the GNU Lesser General Public License 20 along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 21 http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 22 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 23 24 /* define TWO_ARGS for two-argument functions like mpfr_pow 25 define DOUBLE_ARG1 or DOUBLE_ARG2 for function with a double operand in 26 first or second place like sub_d or d_sub */ 27 28 #ifndef TEST_RANDOM_POS 29 /* For the random function: one number on two is negative. */ 30 #define TEST_RANDOM_POS 256 31 #endif 32 33 #ifndef TEST_RANDOM_POS2 34 /* For the random function: one number on two is negative. */ 35 #define TEST_RANDOM_POS2 256 36 #endif 37 38 #ifndef TEST_RANDOM_EMIN 39 #define TEST_RANDOM_EMIN -256 40 #endif 41 42 #ifndef TEST_RANDOM_EMAX 43 #define TEST_RANDOM_EMAX 255 44 #endif 45 46 /* If the MPFR_SUSPICIOUS_OVERFLOW test fails but this is not a bug, 47 then define TGENERIC_SO_TEST with an adequate test (possibly 0) to 48 omit this particular case. */ 49 #ifndef TGENERIC_SO_TEST 50 #define TGENERIC_SO_TEST 1 51 #endif 52 53 #define STR(F) #F 54 #define MAKE_STR(S) STR(S) 55 56 /* The (void *) below is needed to avoid a warning with gcc 4.2+ and functions 57 * with 2 arguments. See <http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36299>. 58 */ 59 #define TGENERIC_FAIL(S, X, U) \ 60 do \ 61 { \ 62 printf ("%s\nx = ", (S)); \ 63 mpfr_out_str (stdout, 2, 0, (X), MPFR_RNDN); \ 64 printf ("\n"); \ 65 if ((void *) U != 0) \ 66 { \ 67 printf ("u = "); \ 68 mpfr_out_str (stdout, 2, 0, (U), MPFR_RNDN); \ 69 printf ("\n"); \ 70 } \ 71 printf ("yprec = %u, rnd_mode = %s, inexact = %d, flags = %u\n", \ 72 (unsigned int) yprec, mpfr_print_rnd_mode (rnd), compare, \ 73 (unsigned int) __gmpfr_flags); \ 74 exit (1); \ 75 } \ 76 while (0) 77 78 #define TGENERIC_CHECK_AUX(S, EXPR, U) \ 79 do \ 80 if (!(EXPR)) \ 81 TGENERIC_FAIL (S " in " MAKE_STR(TEST_FUNCTION), x, U); \ 82 while (0) 83 84 #undef TGENERIC_CHECK 85 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 86 #define TGENERIC_CHECK(S, EXPR) TGENERIC_CHECK_AUX(S, EXPR, u) 87 #else 88 #define TGENERIC_CHECK(S, EXPR) TGENERIC_CHECK_AUX(S, EXPR, 0) 89 #endif 90 91 #ifdef DEBUG_TGENERIC 92 #define TGENERIC_IAUX(F,P,X,U) \ 93 do \ 94 { \ 95 printf ("tgeneric: testing function " STR(F) \ 96 ", %s, target prec = %lu\nx = ", \ 97 mpfr_print_rnd_mode (rnd), (unsigned long) (P)); \ 98 mpfr_out_str (stdout, 2, 0, (X), MPFR_RNDN); \ 99 printf ("\n"); \ 100 if (U) \ 101 { \ 102 printf ("u = "); \ 103 mpfr_out_str (stdout, 2, 0, (U), MPFR_RNDN); \ 104 printf ("\n"); \ 105 } \ 106 } \ 107 while (0) 108 #undef TGENERIC_INFO 109 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 110 #define TGENERIC_INFO(F,P) TGENERIC_IAUX(F,P,x,u) 111 #else 112 #define TGENERIC_INFO(F,P) TGENERIC_IAUX(F,P,x,0) 113 #endif 114 #endif 115 116 /* For some functions (for example cos), the argument reduction is too 117 expensive when using mpfr_get_emax(). Then simply define REDUCE_EMAX 118 to some reasonable value before including tgeneric.c. */ 119 #ifndef REDUCE_EMAX 120 #define REDUCE_EMAX mpfr_get_emax () 121 #endif 122 123 static void 124 test_generic (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int nmax) 125 { 126 mpfr_prec_t prec, xprec, yprec; 127 mpfr_t x, y, z, t, w; 128 #ifdef TWO_ARGS 129 mpfr_t u; 130 #elif defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 131 mpfr_t u; 132 double d; 133 #endif 134 mpfr_rnd_t rnd; 135 int inexact, compare, compare2; 136 unsigned int n; 137 unsigned long ctrt = 0, ctrn = 0; 138 mpfr_exp_t old_emin, old_emax; 139 140 old_emin = mpfr_get_emin (); 141 old_emax = mpfr_get_emax (); 142 143 mpfr_inits2 (MPFR_PREC_MIN, x, y, z, t, w, (mpfr_ptr) 0); 144 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 145 mpfr_init2 (u, MPFR_PREC_MIN); 146 #endif 147 148 /* generic test */ 149 for (prec = p0; prec <= p1; prec++) 150 { 151 mpfr_set_prec (z, prec); 152 mpfr_set_prec (t, prec); 153 yprec = prec + 10; 154 mpfr_set_prec (y, yprec); 155 mpfr_set_prec (w, yprec); 156 157 /* Note: in precision p1, we test 4 special cases. */ 158 for (n = 0; n < (prec == p1 ? nmax + 4 : nmax); n++) 159 { 160 int infinite_input = 0; 161 162 xprec = prec; 163 if (randlimb () & 1) 164 { 165 xprec *= (double) randlimb () / MP_LIMB_T_MAX; 166 if (xprec < MPFR_PREC_MIN) 167 xprec = MPFR_PREC_MIN; 168 } 169 mpfr_set_prec (x, xprec); 170 #ifdef TWO_ARGS 171 mpfr_set_prec (u, xprec); 172 #elif defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 173 mpfr_set_prec (u, IEEE_DBL_MANT_DIG); 174 #endif 175 176 if (n > 3 || prec < p1) 177 { 178 #if defined(RAND_FUNCTION) 179 RAND_FUNCTION (x); 180 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 181 RAND_FUNCTION (u); 182 #endif 183 #else 184 tests_default_random (x, TEST_RANDOM_POS, 185 TEST_RANDOM_EMIN, TEST_RANDOM_EMAX); 186 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 187 tests_default_random (u, TEST_RANDOM_POS2, 188 TEST_RANDOM_EMIN, TEST_RANDOM_EMAX); 189 #endif 190 #endif 191 } 192 else 193 { 194 /* Special cases tested in precision p1 if n <= 3. They are 195 useful really in the extended exponent range. */ 196 #if (defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2)) && defined(MPFR_ERRDIVZERO) 197 goto next_n; 198 #endif 199 set_emin (MPFR_EMIN_MIN); 200 set_emax (MPFR_EMAX_MAX); 201 if (n <= 1) 202 { 203 mpfr_set_si (x, n == 0 ? 1 : -1, MPFR_RNDN); 204 mpfr_set_exp (x, mpfr_get_emin ()); 205 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 206 mpfr_set_si (u, randlimb () % 2 == 0 ? 1 : -1, MPFR_RNDN); 207 mpfr_set_exp (u, mpfr_get_emin ()); 208 #endif 209 } 210 else /* 2 <= n <= 3 */ 211 { 212 if (getenv ("MPFR_CHECK_MAX") == NULL) 213 goto next_n; 214 mpfr_set_si (x, n == 0 ? 1 : -1, MPFR_RNDN); 215 mpfr_setmax (x, REDUCE_EMAX); 216 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 217 mpfr_set_si (u, randlimb () % 2 == 0 ? 1 : -1, MPFR_RNDN); 218 mpfr_setmax (u, mpfr_get_emax ()); 219 #endif 220 } 221 } 222 223 rnd = RND_RAND (); 224 mpfr_clear_flags (); 225 #ifdef DEBUG_TGENERIC 226 TGENERIC_INFO (TEST_FUNCTION, MPFR_PREC (y)); 227 #endif 228 #if defined(TWO_ARGS) 229 compare = TEST_FUNCTION (y, x, u, rnd); 230 #elif defined(DOUBLE_ARG1) 231 d = mpfr_get_d (u, rnd); 232 compare = TEST_FUNCTION (y, d, x, rnd); 233 /* d can be infinite due to overflow in mpfr_get_d */ 234 infinite_input |= DOUBLE_ISINF (d); 235 #elif defined(DOUBLE_ARG2) 236 d = mpfr_get_d (u, rnd); 237 compare = TEST_FUNCTION (y, x, d, rnd); 238 /* d can be infinite due to overflow in mpfr_get_d */ 239 infinite_input |= DOUBLE_ISINF (d); 240 #else 241 compare = TEST_FUNCTION (y, x, rnd); 242 #endif 243 TGENERIC_CHECK ("Bad inexact flag", 244 (compare != 0) ^ (mpfr_inexflag_p () == 0)); 245 ctrt++; 246 /* Consistency test in a reduced exponent range. Doing it 247 for the first 10 samples and for prec == p1 (which has 248 some special cases) should be sufficient. */ 249 if (ctrt <= 10 || prec == p1) 250 { 251 unsigned int flags, oldflags = __gmpfr_flags; 252 mpfr_exp_t e, emin, emax, oemin, oemax; 253 254 /* Determine the smallest exponent range containing the 255 exponents of the mpfr_t inputs (x, and u if TWO_ARGS) 256 and output (y). */ 257 emin = MPFR_EMAX_MAX; 258 emax = MPFR_EMIN_MIN; 259 if (MPFR_IS_PURE_FP (x)) 260 { 261 e = MPFR_GET_EXP (x); 262 if (e < emin) 263 emin = e; 264 if (e > emax) 265 emax = e; 266 } 267 if (MPFR_IS_PURE_FP (y)) 268 { 269 e = MPFR_GET_EXP (y); 270 if (e < emin) 271 emin = e; 272 if (e > emax) 273 emax = e; 274 } 275 #if defined(TWO_ARGS) 276 if (MPFR_IS_PURE_FP (u)) 277 { 278 e = MPFR_GET_EXP (u); 279 if (e < emin) 280 emin = e; 281 if (e > emax) 282 emax = e; 283 } 284 #endif 285 if (emin > emax) 286 emin = emax; /* case where all values are singular */ 287 oemin = mpfr_get_emin (); 288 oemax = mpfr_get_emax (); 289 mpfr_set_emin (emin); 290 mpfr_set_emax (emax); 291 #ifdef DEBUG_TGENERIC 292 /* Useful information in case of assertion failure. */ 293 printf ("tgeneric: reduced exponent range [%" 294 MPFR_EXP_FSPEC "d,%" MPFR_EXP_FSPEC "d]\n", 295 (mpfr_eexp_t) emin, (mpfr_eexp_t) emax); 296 #endif 297 mpfr_clear_flags (); 298 #if defined(TWO_ARGS) 299 inexact = TEST_FUNCTION (w, x, u, rnd); 300 #elif defined(DOUBLE_ARG1) 301 inexact = TEST_FUNCTION (w, d, x, rnd); 302 #elif defined(DOUBLE_ARG2) 303 inexact = TEST_FUNCTION (w, x, d, rnd); 304 #else 305 inexact = TEST_FUNCTION (w, x, rnd); 306 #endif 307 flags = __gmpfr_flags; 308 mpfr_set_emin (oemin); 309 mpfr_set_emax (oemax); 310 if (! (SAME_VAL (w, y) && 311 SAME_SIGN (inexact, compare) && 312 flags == oldflags)) 313 { 314 printf ("Error in " MAKE_STR(TEST_FUNCTION) 315 ", reduced exponent range [%" 316 MPFR_EXP_FSPEC "d,%" MPFR_EXP_FSPEC "d] on:\n", 317 (mpfr_eexp_t) emin, (mpfr_eexp_t) emax); 318 printf ("x = "); 319 mpfr_dump (x); 320 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 321 printf ("u = "); 322 mpfr_dump (u); 323 #endif 324 printf ("yprec = %u, rnd_mode = %s\n", 325 (unsigned int) yprec, mpfr_print_rnd_mode (rnd)); 326 printf ("Expected:\n y = "); 327 mpfr_dump (y); 328 printf (" inex = %d, flags = %u\n", 329 SIGN (compare), oldflags); 330 printf ("Got:\n w = "); 331 mpfr_dump (w); 332 printf (" inex = %d, flags = %u\n", 333 SIGN (inexact), flags); 334 exit (1); 335 } 336 } 337 if (MPFR_IS_SINGULAR (y)) 338 { 339 if (MPFR_IS_NAN (y) || mpfr_nanflag_p ()) 340 TGENERIC_CHECK ("Bad NaN flag", 341 MPFR_IS_NAN (y) && mpfr_nanflag_p ()); 342 else if (MPFR_IS_INF (y)) 343 { 344 TGENERIC_CHECK ("Bad overflow flag", 345 (compare != 0) ^ (mpfr_overflow_p () == 0)); 346 TGENERIC_CHECK ("Bad divide-by-zero flag", 347 (compare == 0 && !infinite_input) ^ 348 (mpfr_divby0_p () == 0)); 349 } 350 else if (MPFR_IS_ZERO (y)) 351 TGENERIC_CHECK ("Bad underflow flag", 352 (compare != 0) ^ (mpfr_underflow_p () == 0)); 353 } 354 else if (mpfr_divby0_p ()) 355 { 356 TGENERIC_CHECK ("Both overflow and divide-by-zero", 357 ! mpfr_overflow_p ()); 358 TGENERIC_CHECK ("Both underflow and divide-by-zero", 359 ! mpfr_underflow_p ()); 360 TGENERIC_CHECK ("Bad compare value (divide-by-zero)", 361 compare == 0); 362 } 363 else if (mpfr_overflow_p ()) 364 { 365 TGENERIC_CHECK ("Both underflow and overflow", 366 ! mpfr_underflow_p ()); 367 TGENERIC_CHECK ("Bad compare value (overflow)", compare != 0); 368 mpfr_nexttoinf (y); 369 TGENERIC_CHECK ("Should have been max MPFR number", 370 MPFR_IS_INF (y)); 371 } 372 else if (mpfr_underflow_p ()) 373 { 374 TGENERIC_CHECK ("Bad compare value (underflow)", compare != 0); 375 mpfr_nexttozero (y); 376 TGENERIC_CHECK ("Should have been min MPFR number", 377 MPFR_IS_ZERO (y)); 378 } 379 else if (mpfr_can_round (y, yprec, rnd, rnd, prec)) 380 { 381 ctrn++; 382 mpfr_set (t, y, rnd); 383 /* Risk of failures are known when some flags are already set 384 before the function call. Do not set the erange flag, as 385 it will remain set after the function call and no checks 386 are performed in such a case (see the mpfr_erangeflag_p 387 test below). */ 388 if (randlimb () & 1) 389 __gmpfr_flags = MPFR_FLAGS_ALL ^ MPFR_FLAGS_ERANGE; 390 #ifdef DEBUG_TGENERIC 391 TGENERIC_INFO (TEST_FUNCTION, MPFR_PREC (z)); 392 #endif 393 /* Let's increase the precision of the inputs in a random way. 394 In most cases, this doesn't make any difference, but for 395 the mpfr_fmod bug fixed in r6230, this triggers the bug. */ 396 mpfr_prec_round (x, mpfr_get_prec (x) + (randlimb () & 15), 397 MPFR_RNDN); 398 #if defined(TWO_ARGS) 399 mpfr_prec_round (u, mpfr_get_prec (u) + (randlimb () & 15), 400 MPFR_RNDN); 401 inexact = TEST_FUNCTION (z, x, u, rnd); 402 #elif defined(DOUBLE_ARG1) 403 inexact = TEST_FUNCTION (z, d, x, rnd); 404 #elif defined(DOUBLE_ARG2) 405 inexact = TEST_FUNCTION (z, x, d, rnd); 406 #else 407 inexact = TEST_FUNCTION (z, x, rnd); 408 #endif 409 if (mpfr_erangeflag_p ()) 410 goto next_n; 411 if (mpfr_nan_p (z) || mpfr_cmp (t, z) != 0) 412 { 413 printf ("results differ for x="); 414 mpfr_out_str (stdout, 2, xprec, x, MPFR_RNDN); 415 #ifdef TWO_ARGS 416 printf ("\nu="); 417 mpfr_out_str (stdout, 2, xprec, u, MPFR_RNDN); 418 #elif defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 419 printf ("\nu="); 420 mpfr_out_str (stdout, 2, IEEE_DBL_MANT_DIG, u, MPFR_RNDN); 421 #endif 422 printf (" prec=%u rnd_mode=%s\n", (unsigned) prec, 423 mpfr_print_rnd_mode (rnd)); 424 printf ("got "); 425 mpfr_out_str (stdout, 2, prec, z, MPFR_RNDN); 426 puts (""); 427 printf ("expected "); 428 mpfr_out_str (stdout, 2, prec, t, MPFR_RNDN); 429 puts (""); 430 printf ("approx "); 431 mpfr_print_binary (y); 432 puts (""); 433 exit (1); 434 } 435 compare2 = mpfr_cmp (t, y); 436 /* if rounding to nearest, cannot know the sign of t - f(x) 437 because of composed rounding: y = o(f(x)) and t = o(y) */ 438 if (compare * compare2 >= 0) 439 compare = compare + compare2; 440 else 441 compare = inexact; /* cannot determine sign(t-f(x)) */ 442 if (((inexact == 0) && (compare != 0)) || 443 ((inexact > 0) && (compare <= 0)) || 444 ((inexact < 0) && (compare >= 0))) 445 { 446 printf ("Wrong inexact flag for rnd=%s: expected %d, got %d" 447 "\n", mpfr_print_rnd_mode (rnd), compare, inexact); 448 printf ("x="); mpfr_print_binary (x); puts (""); 449 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 450 printf ("u="); mpfr_print_binary (u); puts (""); 451 #endif 452 printf ("y="); mpfr_print_binary (y); puts (""); 453 printf ("t="); mpfr_print_binary (t); puts (""); 454 exit (1); 455 } 456 } 457 else if (getenv ("MPFR_SUSPICIOUS_OVERFLOW") != NULL) 458 { 459 /* For developers only! */ 460 MPFR_ASSERTN (MPFR_IS_PURE_FP (y)); 461 mpfr_nexttoinf (y); 462 if (MPFR_IS_INF (y) && MPFR_IS_LIKE_RNDZ (rnd, MPFR_IS_NEG (y)) 463 && !mpfr_overflow_p () && TGENERIC_SO_TEST) 464 { 465 printf ("Possible bug! |y| is the maximum finite number " 466 "and has been obtained when\nrounding toward zero" 467 " (%s). Thus there is a very probable overflow,\n" 468 "but the overflow flag is not set!\n", 469 mpfr_print_rnd_mode (rnd)); 470 printf ("x="); mpfr_print_binary (x); puts (""); 471 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 472 printf ("u="); mpfr_print_binary (u); puts (""); 473 #endif 474 exit (1); 475 } 476 } 477 478 next_n: 479 /* In case the exponent range has been changed by 480 tests_default_random() or for special values... */ 481 mpfr_set_emin (old_emin); 482 mpfr_set_emax (old_emax); 483 } 484 } 485 486 #ifndef TGENERIC_NOWARNING 487 if (3 * ctrn < 2 * ctrt) 488 printf ("Warning! Too few normal cases in generic tests (%lu / %lu)\n", 489 ctrn, ctrt); 490 #endif 491 492 mpfr_clears (x, y, z, t, w, (mpfr_ptr) 0); 493 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 494 mpfr_clear (u); 495 #endif 496 } 497 498 #undef TEST_RANDOM_POS 499 #undef TEST_RANDOM_POS2 500 #undef TEST_RANDOM_EMIN 501 #undef TEST_RANDOM_EMAX 502 #undef RAND_FUNCTION 503 #undef TWO_ARGS 504 #undef TWO_ARGS_UI 505 #undef TEST_FUNCTION 506 #undef test_generic 507
Indexes created Sun Jul 12 00:26:10 UTC 2026