tgeneric.c revision 1.1.1.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 Free Software Foundation, Inc. 5 Contributed by the Arenaire and Cacao 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 /* The (void *) below is needed to avoid a warning with gcc 4.2+ and functions 47 * with 2 arguments. See <http://gcc.gnu.org/bugzilla/show_bug.cgi?id=36299>. 48 */ 49 #define TGENERIC_FAIL(S, X, U) \ 50 do \ 51 { \ 52 printf ("%s\nx = ", (S)); \ 53 mpfr_out_str (stdout, 2, 0, (X), MPFR_RNDN); \ 54 printf ("\n"); \ 55 if ((void *) U != 0) \ 56 { \ 57 printf ("u = "); \ 58 mpfr_out_str (stdout, 2, 0, (U), MPFR_RNDN); \ 59 printf ("\n"); \ 60 } \ 61 printf ("yprec = %u, rnd_mode = %s, inexact = %d, flags = %u\n", \ 62 (unsigned int) yprec, mpfr_print_rnd_mode (rnd), compare, \ 63 (unsigned int) __gmpfr_flags); \ 64 exit (1); \ 65 } \ 66 while (0) 67 68 #undef TGENERIC_CHECK 69 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 70 #define TGENERIC_CHECK(S, EXPR) \ 71 do if (!(EXPR)) TGENERIC_FAIL (S, x, u); while (0) 72 #else 73 #define TGENERIC_CHECK(S, EXPR) \ 74 do if (!(EXPR)) TGENERIC_FAIL (S, x, 0); while (0) 75 #endif 76 77 #ifdef DEBUG_TGENERIC 78 #define STR(F) #F 79 #define TGENERIC_IAUX(F,P,X,U) \ 80 do \ 81 { \ 82 printf ("tgeneric: testing function " STR(F) \ 83 ", %s, target prec = %lu\nx = ", \ 84 mpfr_print_rnd_mode (rnd), (unsigned long) (P)); \ 85 mpfr_out_str (stdout, 2, 0, (X), MPFR_RNDN); \ 86 printf ("\n"); \ 87 if (U) \ 88 { \ 89 printf ("u = "); \ 90 mpfr_out_str (stdout, 2, 0, (U), MPFR_RNDN); \ 91 printf ("\n"); \ 92 } \ 93 } \ 94 while (0) 95 #undef TGENERIC_INFO 96 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 97 #define TGENERIC_INFO(F,P) TGENERIC_IAUX(F,P,x,u) 98 #else 99 #define TGENERIC_INFO(F,P) TGENERIC_IAUX(F,P,x,0) 100 #endif 101 #endif 102 103 /* For some functions (for example cos), the argument reduction is too 104 expensive when using mpfr_get_emax(). Then simply define REDUCE_EMAX 105 to some reasonable value before including tgeneric.c. */ 106 #ifndef REDUCE_EMAX 107 #define REDUCE_EMAX mpfr_get_emax () 108 #endif 109 110 static void 111 test_generic (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N) 112 { 113 mpfr_prec_t prec, xprec, yprec; 114 mpfr_t x, y, z, t, w; 115 #ifdef TWO_ARGS 116 mpfr_t u; 117 #elif defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 118 mpfr_t u; 119 double d; 120 #endif 121 mpfr_rnd_t rnd; 122 int inexact, compare, compare2; 123 unsigned int n; 124 unsigned long ctrt = 0, ctrn = 0; 125 mpfr_exp_t old_emin, old_emax; 126 127 old_emin = mpfr_get_emin (); 128 old_emax = mpfr_get_emax (); 129 130 mpfr_inits2 (MPFR_PREC_MIN, x, y, z, t, w, (mpfr_ptr) 0); 131 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 132 mpfr_init2 (u, MPFR_PREC_MIN); 133 #endif 134 135 /* generic test */ 136 for (prec = p0; prec <= p1; prec++) 137 { 138 mpfr_set_prec (z, prec); 139 mpfr_set_prec (t, prec); 140 yprec = prec + 10; 141 mpfr_set_prec (y, yprec); 142 mpfr_set_prec (w, yprec); 143 144 /* Note: in precision p1, we test 4 special cases. */ 145 for (n = 0; n < (prec == p1 ? N + 4 : N); n++) 146 { 147 xprec = prec; 148 if (randlimb () & 1) 149 { 150 xprec *= (double) randlimb () / MP_LIMB_T_MAX; 151 if (xprec < MPFR_PREC_MIN) 152 xprec = MPFR_PREC_MIN; 153 } 154 mpfr_set_prec (x, xprec); 155 #ifdef TWO_ARGS 156 mpfr_set_prec (u, xprec); 157 #elif defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 158 mpfr_set_prec (u, IEEE_DBL_MANT_DIG); 159 #endif 160 161 if (n > 3 || prec < p1) 162 { 163 #if defined(RAND_FUNCTION) 164 RAND_FUNCTION (x); 165 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 166 RAND_FUNCTION (u); 167 #endif 168 #else 169 tests_default_random (x, TEST_RANDOM_POS, 170 TEST_RANDOM_EMIN, TEST_RANDOM_EMAX); 171 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 172 tests_default_random (u, TEST_RANDOM_POS2, 173 TEST_RANDOM_EMIN, TEST_RANDOM_EMAX); 174 #endif 175 #endif 176 } 177 else 178 { 179 /* Special cases tested in precision p1 if n <= 3. They are 180 useful really in the extended exponent range. */ 181 set_emin (MPFR_EMIN_MIN); 182 set_emax (MPFR_EMAX_MAX); 183 if (n <= 1) 184 { 185 mpfr_set_si (x, n == 0 ? 1 : -1, MPFR_RNDN); 186 mpfr_set_exp (x, mpfr_get_emin ()); 187 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 188 mpfr_set_si (u, randlimb () % 2 == 0 ? 1 : -1, MPFR_RNDN); 189 mpfr_set_exp (u, mpfr_get_emin ()); 190 #endif 191 } 192 else /* 2 <= n <= 3 */ 193 { 194 if (getenv ("MPFR_CHECK_MAX") == NULL) 195 goto next_n; 196 mpfr_set_si (x, n == 0 ? 1 : -1, MPFR_RNDN); 197 mpfr_setmax (x, REDUCE_EMAX); 198 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 199 mpfr_set_si (u, randlimb () % 2 == 0 ? 1 : -1, MPFR_RNDN); 200 mpfr_setmax (u, mpfr_get_emax ()); 201 #endif 202 } 203 } 204 205 rnd = RND_RAND (); 206 mpfr_clear_flags (); 207 #ifdef DEBUG_TGENERIC 208 TGENERIC_INFO (TEST_FUNCTION, MPFR_PREC (y)); 209 #endif 210 #if defined(TWO_ARGS) 211 compare = TEST_FUNCTION (y, x, u, rnd); 212 #elif defined(DOUBLE_ARG1) 213 d = mpfr_get_d (u, rnd); 214 compare = TEST_FUNCTION (y, d, x, rnd); 215 #elif defined(DOUBLE_ARG2) 216 d = mpfr_get_d (u, rnd); 217 compare = TEST_FUNCTION (y, x, d, rnd); 218 #else 219 compare = TEST_FUNCTION (y, x, rnd); 220 #endif 221 TGENERIC_CHECK ("Bad inexact flag", 222 (compare != 0) ^ (mpfr_inexflag_p () == 0)); 223 ctrt++; 224 /* Consistency test in a reduced exponent range. Doing it 225 for the first 10 samples and for prec == p1 (which has 226 some special cases) should be sufficient. */ 227 if (ctrt <= 10 || prec == p1) 228 { 229 unsigned int flags, oldflags = __gmpfr_flags; 230 mpfr_exp_t e, emin, emax, oemin, oemax; 231 232 /* Determine the smallest exponent range containing the 233 exponents of the mpfr_t inputs (x, and u if TWO_ARGS) 234 and output (y). */ 235 emin = MPFR_EMAX_MAX; 236 emax = MPFR_EMIN_MIN; 237 if (MPFR_IS_PURE_FP (x)) 238 { 239 e = MPFR_GET_EXP (x); 240 if (e < emin) 241 emin = e; 242 if (e > emax) 243 emax = e; 244 } 245 if (MPFR_IS_PURE_FP (y)) 246 { 247 e = MPFR_GET_EXP (y); 248 if (e < emin) 249 emin = e; 250 if (e > emax) 251 emax = e; 252 } 253 #if defined(TWO_ARGS) 254 if (MPFR_IS_PURE_FP (u)) 255 { 256 e = MPFR_GET_EXP (u); 257 if (e < emin) 258 emin = e; 259 if (e > emax) 260 emax = e; 261 } 262 #endif 263 if (emin > emax) 264 emin = emax; /* case where all values are singular */ 265 oemin = mpfr_get_emin (); 266 oemax = mpfr_get_emax (); 267 mpfr_set_emin (emin); 268 mpfr_set_emax (emax); 269 mpfr_clear_flags (); 270 #if defined(TWO_ARGS) 271 inexact = TEST_FUNCTION (w, x, u, rnd); 272 #elif defined(DOUBLE_ARG1) 273 inexact = TEST_FUNCTION (w, d, x, rnd); 274 #elif defined(DOUBLE_ARG2) 275 inexact = TEST_FUNCTION (w, x, d, rnd); 276 #else 277 inexact = TEST_FUNCTION (w, x, rnd); 278 #endif 279 flags = __gmpfr_flags; 280 mpfr_set_emin (oemin); 281 mpfr_set_emax (oemax); 282 if (! (SAME_VAL (w, y) && 283 SAME_SIGN (inexact, compare) && 284 flags == oldflags)) 285 { 286 printf ("Error in reduced exponent range [%" 287 MPFR_EXP_FSPEC "d,%" MPFR_EXP_FSPEC "d] on:\n", 288 (mpfr_eexp_t) emin, (mpfr_eexp_t) emax); 289 printf ("x = "); 290 mpfr_dump (x); 291 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 292 printf ("u = "); 293 mpfr_dump (u); 294 #endif 295 printf ("yprec = %u, rnd_mode = %s\n", 296 (unsigned int) yprec, mpfr_print_rnd_mode (rnd)); 297 printf ("Expected:\n y = "); 298 mpfr_dump (y); 299 printf (" inex = %d, flags = %u\n", 300 SIGN (compare), oldflags); 301 printf ("Got:\n w = "); 302 mpfr_dump (w); 303 printf (" inex = %d, flags = %u\n", 304 SIGN (inexact), flags); 305 exit (1); 306 } 307 } 308 if (MPFR_IS_SINGULAR (y)) 309 { 310 if (MPFR_IS_NAN (y) || mpfr_nanflag_p ()) 311 TGENERIC_CHECK ("Bad NaN flag", 312 MPFR_IS_NAN (y) && mpfr_nanflag_p ()); 313 else if (MPFR_IS_INF (y)) 314 TGENERIC_CHECK ("Bad overflow flag", 315 (compare != 0) ^ (mpfr_overflow_p () == 0)); 316 else if (MPFR_IS_ZERO (y)) 317 TGENERIC_CHECK ("Bad underflow flag", 318 (compare != 0) ^ (mpfr_underflow_p () == 0)); 319 } 320 else if (mpfr_overflow_p ()) 321 { 322 TGENERIC_CHECK ("Bad compare value (overflow)", compare != 0); 323 mpfr_nexttoinf (y); 324 TGENERIC_CHECK ("Should have been max MPFR number", 325 MPFR_IS_INF (y)); 326 } 327 else if (mpfr_underflow_p ()) 328 { 329 TGENERIC_CHECK ("Bad compare value (underflow)", compare != 0); 330 mpfr_nexttozero (y); 331 TGENERIC_CHECK ("Should have been min MPFR number", 332 MPFR_IS_ZERO (y)); 333 } 334 else if (mpfr_can_round (y, yprec, rnd, rnd, prec)) 335 { 336 ctrn++; 337 mpfr_set (t, y, rnd); 338 /* Risk of failures are known when some flags are already set 339 before the function call. Do not set the erange flag, as 340 it will remain set after the function call and no checks 341 are performed in such a case (see the mpfr_erangeflag_p 342 test below). */ 343 if (randlimb () & 1) 344 __gmpfr_flags = MPFR_FLAGS_ALL ^ MPFR_FLAGS_ERANGE; 345 #ifdef DEBUG_TGENERIC 346 TGENERIC_INFO (TEST_FUNCTION, MPFR_PREC (z)); 347 #endif 348 /* Let's increase the precision of the inputs in a random way. 349 In most cases, this doesn't make any difference, but for 350 the mpfr_fmod bug fixed in r6230, this triggers the bug. */ 351 mpfr_prec_round (x, mpfr_get_prec (x) + (randlimb () & 15), 352 MPFR_RNDN); 353 #if defined(TWO_ARGS) 354 mpfr_prec_round (u, mpfr_get_prec (u) + (randlimb () & 15), 355 MPFR_RNDN); 356 inexact = TEST_FUNCTION (z, x, u, rnd); 357 #elif defined(DOUBLE_ARG1) 358 inexact = TEST_FUNCTION (z, d, x, rnd); 359 #elif defined(DOUBLE_ARG2) 360 inexact = TEST_FUNCTION (z, x, d, rnd); 361 #else 362 inexact = TEST_FUNCTION (z, x, rnd); 363 #endif 364 if (mpfr_erangeflag_p ()) 365 goto next_n; 366 if (mpfr_nan_p (z) || mpfr_cmp (t, z) != 0) 367 { 368 printf ("results differ for x="); 369 mpfr_out_str (stdout, 2, xprec, x, MPFR_RNDN); 370 #ifdef TWO_ARGS 371 printf ("\nu="); 372 mpfr_out_str (stdout, 2, xprec, u, MPFR_RNDN); 373 #elif defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 374 printf ("\nu="); 375 mpfr_out_str (stdout, 2, IEEE_DBL_MANT_DIG, u, MPFR_RNDN); 376 #endif 377 printf (" prec=%u rnd_mode=%s\n", (unsigned) prec, 378 mpfr_print_rnd_mode (rnd)); 379 printf ("got "); 380 mpfr_out_str (stdout, 2, prec, z, MPFR_RNDN); 381 puts (""); 382 printf ("expected "); 383 mpfr_out_str (stdout, 2, prec, t, MPFR_RNDN); 384 puts (""); 385 printf ("approx "); 386 mpfr_print_binary (y); 387 puts (""); 388 exit (1); 389 } 390 compare2 = mpfr_cmp (t, y); 391 /* if rounding to nearest, cannot know the sign of t - f(x) 392 because of composed rounding: y = o(f(x)) and t = o(y) */ 393 if (compare * compare2 >= 0) 394 compare = compare + compare2; 395 else 396 compare = inexact; /* cannot determine sign(t-f(x)) */ 397 if (((inexact == 0) && (compare != 0)) || 398 ((inexact > 0) && (compare <= 0)) || 399 ((inexact < 0) && (compare >= 0))) 400 { 401 printf ("Wrong inexact flag for rnd=%s: expected %d, got %d" 402 "\n", mpfr_print_rnd_mode (rnd), compare, inexact); 403 printf ("x="); mpfr_print_binary (x); puts (""); 404 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 405 printf ("u="); mpfr_print_binary (u); puts (""); 406 #endif 407 printf ("y="); mpfr_print_binary (y); puts (""); 408 printf ("t="); mpfr_print_binary (t); puts (""); 409 exit (1); 410 } 411 } 412 else if (getenv ("MPFR_SUSPICIOUS_OVERFLOW") != NULL) 413 { 414 /* For developers only! */ 415 MPFR_ASSERTN (MPFR_IS_PURE_FP (y)); 416 mpfr_nexttoinf (y); 417 if (MPFR_IS_INF (y) && MPFR_IS_LIKE_RNDZ (rnd, MPFR_IS_NEG (y)) 418 && !mpfr_overflow_p ()) 419 { 420 printf ("Possible bug! |y| is the maximum finite number " 421 "and has been obtained when\nrounding toward zero" 422 " (%s). Thus there is a very probable overflow,\n" 423 "but the overflow flag is not set!\n", 424 mpfr_print_rnd_mode (rnd)); 425 printf ("x="); mpfr_print_binary (x); puts (""); 426 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 427 printf ("u="); mpfr_print_binary (u); puts (""); 428 #endif 429 exit (1); 430 } 431 } 432 433 next_n: 434 /* In case the exponent range has been changed by 435 tests_default_random() or for special values... */ 436 mpfr_set_emin (old_emin); 437 mpfr_set_emax (old_emax); 438 } 439 } 440 441 #ifndef TGENERIC_NOWARNING 442 if (3 * ctrn < 2 * ctrt) 443 printf ("Warning! Too few normal cases in generic tests (%lu / %lu)\n", 444 ctrn, ctrt); 445 #endif 446 447 mpfr_clears (x, y, z, t, w, (mpfr_ptr) 0); 448 #if defined(TWO_ARGS) || defined(DOUBLE_ARG1) || defined(DOUBLE_ARG2) 449 mpfr_clear (u); 450 #endif 451 } 452 453 #undef TEST_RANDOM_POS 454 #undef TEST_RANDOM_POS2 455 #undef TEST_RANDOM_EMIN 456 #undef TEST_RANDOM_EMAX 457 #undef RAND_FUNCTION 458 #undef TWO_ARGS 459 #undef TWO_ARGS_UI 460 #undef TEST_FUNCTION 461 #undef test_generic 462
Indexes created Sun Jul 12 00:26:10 UTC 2026