tests.c revision 1.1.1.5
1/* Miscellaneous support for test programs. 2 3Copyright 2001-2020 Free Software Foundation, Inc. 4Contributed by the AriC and Caramba projects, INRIA. 5 6This file is part of the GNU MPFR Library. 7 8The GNU MPFR Library is free software; you can redistribute it and/or modify 9it under the terms of the GNU Lesser General Public License as published by 10the Free Software Foundation; either version 3 of the License, or (at your 11option) any later version. 12 13The GNU MPFR Library is distributed in the hope that it will be useful, but 14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 15or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 16License for more details. 17 18You should have received a copy of the GNU Lesser General Public License 19along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 20https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 2151 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 22 23#ifdef HAVE_CONFIG_H 24# include "config.h" 25#endif 26 27#include <float.h> 28 29#ifdef HAVE_LOCALE_H 30#include <locale.h> 31#endif 32 33#ifdef MPFR_TESTS_FPE_DIV 34# ifdef MPFR_TESTS_FPE_TRAP 35# define _GNU_SOURCE /* for feenableexcept */ 36# endif 37# include <fenv.h> 38#endif 39 40#ifdef TIME_WITH_SYS_TIME 41# include <sys/time.h> /* for struct timeval */ 42# include <time.h> 43#elif defined HAVE_SYS_TIME_H 44# include <sys/time.h> 45#else 46# include <time.h> 47#endif 48 49/* <sys/fpu.h> is needed to have union fpc_csr defined under IRIX64 50 (see below). Let's include it only if need be. */ 51#if defined HAVE_SYS_FPU_H && defined HAVE_FPC_CSR 52# include <sys/fpu.h> 53#endif 54 55#ifdef MPFR_TESTS_TIMEOUT 56#include <sys/resource.h> 57#endif 58 59#if defined(HAVE_SIGNAL) || defined(HAVE_SIGACTION) 60# include <signal.h> 61#endif 62 63#include "mpfr-test.h" 64 65#ifdef MPFR_FPU_PREC 66/* This option allows to test MPFR on x86 processors when the FPU 67 * rounding precision has been changed. As MPFR is a library, this can 68 * occur in practice, either by the calling software or by some other 69 * library or plug-in used by the calling software. This option is 70 * mainly for developers. If it is used, then the <fpu_control.h> 71 * header is assumed to exist and work like under Linux/x86. MPFR does 72 * not need to be recompiled. So, a possible usage is the following: 73 * 74 * cd tests 75 * make clean 76 * make check CFLAGS="-g -O2 -ffloat-store -DMPFR_FPU_PREC=_FPU_SINGLE" 77 * 78 * i.e. just add -DMPFR_FPU_PREC=... to the CFLAGS found in Makefile. 79 * 80 * Notes: 81 * + SSE2 (used to implement double's on x86_64, and possibly on x86 82 * too, depending on the compiler configuration and flags) is not 83 * affected by the dynamic precision. 84 * + When the FPU is set to single precision, the behavior of MPFR 85 * functions that have a native floating-point type (float, double, 86 * long double) as argument or return value is not guaranteed. 87 */ 88 89#include <fpu_control.h> 90 91static void 92set_fpu_prec (void) 93{ 94 fpu_control_t cw; 95 96 _FPU_GETCW(cw); 97 cw &= ~(_FPU_EXTENDED|_FPU_DOUBLE|_FPU_SINGLE); 98 cw |= (MPFR_FPU_PREC); 99 _FPU_SETCW(cw); 100} 101 102#endif 103 104char mpfr_rands_initialized = 0; 105gmp_randstate_t mpfr_rands; 106 107char *locale = NULL; 108 109/* Programs that test GMP's mp_set_memory_functions() need to set 110 tests_memory_disabled = 2 before calling tests_start_mpfr(). */ 111#ifdef MPFR_USE_MINI_GMP 112/* disable since mini-gmp does not keep track of old_size in realloc/free */ 113int tests_memory_disabled = 1; 114#else 115int tests_memory_disabled = 0; 116#endif 117 118static mpfr_exp_t default_emin, default_emax; 119 120static void tests_rand_start (void); 121static void tests_rand_end (void); 122static void tests_limit_start (void); 123 124/* We want to always import the function mpfr_dump inside the test 125 suite, so that we can use it in GDB. But it doesn't work if we build 126 a Windows DLL (initializer element is not a constant) */ 127#if !__GMP_LIBGMP_DLL 128extern void (*dummy_func) (mpfr_srcptr); 129void (*dummy_func)(mpfr_srcptr) = mpfr_dump; 130#endif 131 132/* Various version checks. 133 A mismatch on the GMP version is not regarded as fatal. A mismatch 134 on the MPFR version is regarded as fatal, since this means that we 135 would not check the MPFR library that has just been built (the goal 136 of "make check") but a different library that is already installed, 137 i.e. any test result would be meaningless; in such a case, we exit 138 immediately with an error (exit status = 1). 139 Return value: 0 for no errors, 1 in case of any non-fatal error. 140 Note: If the return value is 0, no data must be sent to stdout. */ 141int 142test_version (void) 143{ 144 const char *version; 145 char buffer[256]; 146 int err = 0; 147 148#ifndef MPFR_USE_MINI_GMP 149 sprintf (buffer, "%d.%d.%d", __GNU_MP_VERSION, __GNU_MP_VERSION_MINOR, 150 __GNU_MP_VERSION_PATCHLEVEL); 151 if (strcmp (buffer, gmp_version) != 0 && 152 (__GNU_MP_VERSION_PATCHLEVEL != 0 || 153 (sprintf (buffer, "%d.%d", __GNU_MP_VERSION, __GNU_MP_VERSION_MINOR), 154 strcmp (buffer, gmp_version) != 0))) 155 err = 1; 156 157 /* In some cases, it may be acceptable to have different versions for 158 the header and the library, in particular when shared libraries are 159 used (e.g., after a bug-fix upgrade of the library, and versioning 160 ensures that this can be done only when the binary interface is 161 compatible). However, when recompiling software like here, this 162 should never happen (except if GMP has been upgraded between two 163 "make check" runs, but there's no reason for that). A difference 164 between the versions of gmp.h and libgmp probably indicates either 165 a bad configuration or some other inconsistency in the development 166 environment, and it is better to fail (in particular for automatic 167 installations). */ 168 if (err) 169 { 170 printf ("ERROR! The versions of gmp.h (%s) and libgmp (%s) do not " 171 "match.\nThe possible causes are:\n", buffer, gmp_version); 172 printf (" * A bad configuration in your include/library search paths.\n" 173 " * An inconsistency in the include/library search paths of\n" 174 " your development environment; an example:\n" 175 " " 176 "https://gcc.gnu.org/legacy-ml/gcc-help/2010-11/msg00359.html\n" 177 " * GMP has been upgraded after the first \"make check\".\n" 178 " In such a case, try again after a \"make clean\".\n" 179 " * A new or non-standard version naming is used in GMP.\n" 180 " In this case, a patch may already be available on the\n" 181 " MPFR web site. Otherwise please report the problem.\n"); 182 printf ("In the first two cases, this may lead to errors, in particular" 183 " with MPFR.\nIf some other tests fail, please solve that" 184 " problem first.\n"); 185 } 186#endif 187 188 /* VL: I get the following error on an OpenSUSE machine, and changing 189 the value of shlibpath_overrides_runpath in the libtool file from 190 'no' to 'yes' fixes the problem. */ 191 version = mpfr_get_version (); 192 if (strcmp (MPFR_VERSION_STRING, version) == 0) 193 { 194 int i; 195 196 sprintf (buffer, "%d.%d.%d", MPFR_VERSION_MAJOR, MPFR_VERSION_MINOR, 197 MPFR_VERSION_PATCHLEVEL); 198 for (i = 0; buffer[i] == version[i]; i++) 199 if (buffer[i] == '\0') 200 return err; 201 if (buffer[i] == '\0' && version[i] == '-') 202 return err; 203 printf ("%sMPFR_VERSION_MAJOR.MPFR_VERSION_MINOR.MPFR_VERSION_PATCHLEVEL" 204 " (%s)\nand MPFR_VERSION_STRING (%s) do not match!\nIt seems " 205 "that the mpfr.h file has been corrupted.\n", err ? "\n" : "", 206 buffer, version); 207 } 208 else 209 printf ( 210 "%sIncorrect MPFR version! (%s header vs %s library)\n" 211 "Nothing else has been tested since for this reason, any other test\n" 212 "may fail. Please fix this problem first, as suggested below. It\n" 213 "probably comes from libtool (included in the MPFR tarball), which\n" 214 "is responsible for setting up the search paths depending on the\n" 215 "platform, or automake.\n" 216 " * On some platforms such as Solaris, $LD_LIBRARY_PATH overrides\n" 217 " the rpath, and if the MPFR library is already installed in a\n" 218 " $LD_LIBRARY_PATH directory, you typically get this error. Do\n" 219 " not use $LD_LIBRARY_PATH permanently on such platforms; it may\n" 220 " also break other things.\n" 221 " * You may have an ld option that specifies a library search path\n" 222 " where MPFR can be found, taking the precedence over the path\n" 223 " added by libtool. Check your environment variables, such as\n" 224 " LD_OPTIONS under Solaris. Moreover, under Solaris, the run path\n" 225 " generated by libtool 2.4.6 may be incorrect: the build directory\n" 226 " may not appear first in the run path; set $LD_LIBRARY_PATH to\n" 227 " /path/to/builddir/src/.libs for the tests as a workaround.\n" 228 " * Then look at https://www.mpfr.org/mpfr-current/ for any update.\n" 229 " * Try again on a completely clean source (some errors might come\n" 230 " from a previous build or previous source changes).\n" 231 " * If the error still occurs, you can try to change the value of\n" 232 " shlibpath_overrides_runpath ('yes' or 'no') in the \"libtool\"\n" 233 " file and rebuild MPFR (make clean && make && make check). You\n" 234 " may want to report the problem to the libtool and/or automake\n" 235 " developers, with the effect of this change.\n", 236 err ? "\n" : "", MPFR_VERSION_STRING, version); 237 /* Note about $LD_LIBRARY_PATH under Solaris: 238 * https://en.wikipedia.org/wiki/Rpath#Solaris_ld.so 239 * This cause has been confirmed by a user who got this error. 240 * And about the libtool 2.4.6 bug also concerning Solaris: 241 * https://debbugs.gnu.org/cgi/bugreport.cgi?bug=30222 242 * https://bugs.debian.org/cgi-bin/bugreport.cgi?bug=888059 243 */ 244 exit (1); 245} 246 247/* The inexact exception occurs very often, and is normal. 248 The underflow exception also might occur, for example in test_generic 249 for mpfr_xxx_d functions. Same for overflow. Thus we only check for 250 the division-by-zero and invalid exceptions, which should not occur 251 inside MPFR. */ 252#define FPE_FLAGS (FE_DIVBYZERO | FE_INVALID) 253 254void 255tests_start_mpfr (void) 256{ 257 /* Don't buffer, so output is not lost if a test causes a segv, etc. 258 For stdout, this is important as it will typically be fully buffered 259 by default with "make check". For stderr, the C standard just says 260 that it is not fully buffered (it may be line buffered by default); 261 disabling buffering completely might be useful in some cases. 262 Warning! No operations must have already been done on stdout/stderr 263 (this is a requirement of ISO C, and this is important on AIX). 264 Thus tests_start_mpfr should be called at the beginning of main(), 265 possibly after some variable settings. */ 266 setbuf (stdout, NULL); 267 setbuf (stderr, NULL); 268 269 test_version (); 270 271#if defined HAVE_LOCALE_H && defined HAVE_SETLOCALE 272 /* Added on 2005-07-09. This allows to test MPFR under various 273 locales. New bugs will probably be found, in particular with 274 LC_ALL="tr_TR.ISO8859-9" because of the i/I character... */ 275 locale = setlocale (LC_ALL, ""); 276#endif 277 278#ifdef MPFR_FPU_PREC 279 set_fpu_prec (); 280#endif 281 282#ifdef MPFR_TESTS_FPE_DIV 283 /* Define to test the use of MPFR_ERRDIVZERO */ 284 feclearexcept (FE_ALL_EXCEPT); 285# ifdef MPFR_TESTS_FPE_TRAP 286 /* to trap the corresponding FP exceptions */ 287 feenableexcept (FPE_FLAGS); 288# endif 289#endif 290 291 if (tests_memory_disabled != 2) 292 { 293 if (tests_memory_disabled == 0) 294 tests_memory_start (); 295 tests_rand_start (); 296 } 297 tests_limit_start (); 298 299 default_emin = mpfr_get_emin (); 300 default_emax = mpfr_get_emax (); 301} 302 303void 304tests_end_mpfr (void) 305{ 306 int err = 0; 307 308 if (mpfr_get_emin () != default_emin) 309 { 310 printf ("Default emin value has not been restored!\n"); 311 err = 1; 312 } 313 314 if (mpfr_get_emax () != default_emax) 315 { 316 printf ("Default emax value has not been restored!\n"); 317 err = 1; 318 } 319 320 mpfr_free_cache (); 321 mpfr_free_cache2 (MPFR_FREE_GLOBAL_CACHE); 322 if (tests_memory_disabled != 2) 323 { 324 tests_rand_end (); 325 if (tests_memory_disabled == 0) 326 tests_memory_end (); 327 } 328 329#ifdef MPFR_TESTS_FPE_DIV 330 /* Define to test the use of MPFR_ERRDIVZERO */ 331 if (fetestexcept (FPE_FLAGS)) 332 { 333 /* With MPFR_ERRDIVZERO, such exceptions should never occur 334 because the purpose of defining MPFR_ERRDIVZERO is to avoid 335 all the FP divisions by 0. */ 336 printf ("Some floating-point exception(s) occurred:"); 337 if (fetestexcept (FE_DIVBYZERO)) 338 printf (" DIVBYZERO"); /* e.g. from 1.0 / 0.0 to generate an inf */ 339 if (fetestexcept (FE_INVALID)) 340 printf (" INVALID"); /* e.g. from 0.0 / 0.0 to generate a NaN */ 341 printf ("\n"); 342 err = 1; 343 } 344#endif 345 346 if (err) 347 exit (err); 348} 349 350static void 351tests_limit_start (void) 352{ 353#ifdef MPFR_TESTS_TIMEOUT 354 struct rlimit rlim[1]; 355 char *timeoutp; 356 int timeout; 357 358 timeoutp = getenv ("MPFR_TESTS_TIMEOUT"); 359 timeout = timeoutp != NULL ? atoi (timeoutp) : MPFR_TESTS_TIMEOUT; 360 if (timeout > 0) 361 { 362 /* We need to call getrlimit first to initialize rlim_max to 363 an acceptable value for setrlimit. When enabled, timeouts 364 are regarded as important: we don't want to take too much 365 CPU time on machines shared with other users. So, if we 366 can't set the timeout, we exit immediately. */ 367 if (getrlimit (RLIMIT_CPU, rlim)) 368 { 369 printf ("Error: getrlimit failed\n"); 370 exit (1); 371 } 372 rlim->rlim_cur = timeout; 373 if (setrlimit (RLIMIT_CPU, rlim)) 374 { 375 printf ("Error: setrlimit failed\n"); 376 exit (1); 377 } 378 } 379#endif 380} 381 382static void 383tests_rand_start (void) 384{ 385 char *perform_seed; 386 unsigned long seed; 387 388 if (mpfr_rands_initialized) 389 { 390 printf ( 391 "Please let tests_start() initialize the global mpfr_rands, i.e.\n" 392 "ensure that function is called before the first use of RANDS.\n"); 393 exit (1); 394 } 395 396 gmp_randinit_default (mpfr_rands); 397 mpfr_rands_initialized = 1; 398 399 perform_seed = getenv ("GMP_CHECK_RANDOMIZE"); 400 if (perform_seed != NULL) 401 { 402 seed = strtoul (perform_seed, NULL, 10); 403 if (! (seed == 0 || seed == 1)) 404 { 405 printf ("Re-seeding with GMP_CHECK_RANDOMIZE=%lu\n", seed); 406 gmp_randseed_ui (mpfr_rands, seed); 407 } 408 else 409 { 410#ifdef HAVE_GETTIMEOFDAY 411 struct timeval tv; 412 gettimeofday (&tv, NULL); 413 /* Note: If time_t is a "floating type" (as allowed by ISO C99), 414 the cast below can yield undefined behavior. But this would 415 be uncommon (gettimeofday() is specified by POSIX only and 416 POSIX requires time_t to be an integer type) and this line 417 is not executed by default. So, this should be OK. Moreover, 418 gettimeofday() is marked obsolescent by POSIX.1-2008. */ 419 seed = 1000000 * (unsigned long) tv.tv_sec + tv.tv_usec; 420#else 421 time_t tv; 422 time (&tv); 423 seed = tv; 424#endif 425 gmp_randseed_ui (mpfr_rands, seed); 426 printf ("Seed GMP_CHECK_RANDOMIZE=%lu " 427 "(include this in bug reports)\n", seed); 428 } 429 } 430 else 431 gmp_randseed_ui (mpfr_rands, 0x2143FEDC); 432} 433 434static void 435tests_rand_end (void) 436{ 437 RANDS_CLEAR (); 438} 439 440/* initialization function for tests using the hardware floats 441 Not very useful now. */ 442void 443mpfr_test_init (void) 444{ 445#ifdef HAVE_FPC_CSR 446 /* to get subnormal numbers on IRIX64 */ 447 union fpc_csr exp; 448 449 exp.fc_word = get_fpc_csr(); 450 exp.fc_struct.flush = 0; 451 set_fpc_csr(exp.fc_word); 452#endif 453 454#ifdef HAVE_SUBNORM_DBL 455 { 456 double d = DBL_MIN; 457 if (2.0 * (d / 2.0) != d) 458 { 459 printf ("Error: HAVE_SUBNORM_DBL defined, but no subnormals.\n"); 460 exit (1); 461 } 462 } 463#endif 464 465 /* generate DBL_EPSILON with a loop to avoid that the compiler 466 optimizes the code below in non-IEEE 754 mode, deciding that 467 c = d is always false. */ 468#if 0 469 for (eps = 1.0; eps != DBL_EPSILON; eps /= 2.0); 470 c = 1.0 + eps; 471 d = eps * (1.0 - eps) / 2.0; 472 d += c; 473 if (c != d) 474 { 475 printf ("Warning: IEEE 754 standard not fully supported\n" 476 " (maybe extended precision not disabled)\n" 477 " Some tests may fail\n"); 478 } 479#endif 480} 481 482 483/* generate a random limb */ 484mp_limb_t 485randlimb (void) 486{ 487 mp_limb_t limb; 488 489 mpfr_rand_raw (&limb, RANDS, GMP_NUMB_BITS); 490 return limb; 491} 492 493/* returns ulp(x) for x a 'normal' double-precision number */ 494double 495Ulp (double x) 496{ 497 double y, eps; 498 499 if (x < 0) x = -x; 500 501 y = x * 2.220446049250313080847263336181640625e-16 ; /* x / 2^52 */ 502 503 /* as ulp(x) <= y = x/2^52 < 2*ulp(x), 504 we have x + ulp(x) <= x + y <= x + 2*ulp(x), 505 therefore o(x + y) = x + ulp(x) or x + 2*ulp(x) */ 506 507 eps = x + y; 508 eps = eps - x; /* ulp(x) or 2*ulp(x) */ 509 510 return (eps > y) ? 0.5 * eps : eps; 511} 512 513/* returns the number of ulp's between a and b, 514 where a and b can be any floating-point number, except NaN 515 */ 516int 517ulp (double a, double b) 518{ 519 double twoa; 520 521 if (a == b) return 0; /* also deals with a=b=inf or -inf */ 522 523 twoa = a + a; 524 if (twoa == a) /* a is +/-0.0 or +/-Inf */ 525 return ((b < a) ? INT_MAX : -INT_MAX); 526 527 return (int) ((a - b) / Ulp (a)); 528} 529 530/* return double m*2^e */ 531double 532dbl (double m, int e) 533{ 534 if (e >=0 ) 535 while (e-- > 0) 536 m *= 2.0; 537 else 538 while (e++ < 0) 539 m /= 2.0; 540 return m; 541} 542 543/* Warning: NaN values cannot be distinguished if MPFR_NANISNAN is defined. */ 544int 545Isnan (double d) 546{ 547 return (d) != (d); 548} 549 550void 551d_trace (const char *name, double d) 552{ 553 union { 554 double d; 555 unsigned char b[sizeof(double)]; 556 } u; 557 int i; 558 559 if (name != NULL && name[0] != '\0') 560 printf ("%s=", name); 561 562 u.d = d; 563 printf ("["); 564 for (i = 0; i < (int) sizeof (u.b); i++) 565 { 566 if (i != 0) 567 printf (" "); 568 printf ("%02X", (int) u.b[i]); 569 } 570 printf ("] %.20g\n", d); 571} 572 573void 574ld_trace (const char *name, long double ld) 575{ 576 union { 577 long double ld; 578 unsigned char b[sizeof(long double)]; 579 } u; 580 int i; 581 582 if (name != NULL && name[0] != '\0') 583 printf ("%s=", name); 584 585 u.ld = ld; 586 printf ("["); 587 for (i = 0; i < (int) sizeof (u.b); i++) 588 { 589 if (i != 0) 590 printf (" "); 591 printf ("%02X", (int) u.b[i]); 592 } 593 printf ("] %.20Lg\n", ld); 594} 595 596void 597n_trace (const char *name, mp_limb_t *p, mp_size_t n) 598{ 599 unsigned char *buf; 600 size_t bufsize; 601 mp_size_t i, m; 602 603 if (name != NULL && name[0] != '\0') 604 printf ("%s=", name); 605 606 /* similar to gmp_printf ("%NX\n",...), which is not available 607 with mini-gmp */ 608 bufsize = 2 + ((mpfr_prec_t) n * GMP_NUMB_BITS - 1) / 4; 609 buf = (unsigned char *) tests_allocate (bufsize); 610 m = mpn_get_str (buf, 16, p, n); 611 i = 0; 612 while (i < m - 1 && buf[i] == 0) 613 i++; /* skip leading zeros (keeping at least one digit) */ 614 while (i < m) 615 putchar ("0123456789ABCDEF"[buf[i++]]); 616 putchar ('\n'); 617 tests_free (buf, bufsize); 618} 619 620/* Open a file in the SRCDIR directory, i.e. the "tests" source directory, 621 which is different from the current directory when objdir is different 622 from srcdir. One should generally use this function instead of fopen 623 directly. */ 624FILE * 625src_fopen (const char *filename, const char *mode) 626{ 627#ifndef SRCDIR 628 return fopen (filename, mode); 629#else 630 const char *srcdir = SRCDIR; 631 char *buffer; 632 size_t buffsize; 633 FILE *f; 634 635 buffsize = strlen (filename) + strlen (srcdir) + 2; 636 buffer = (char *) tests_allocate (buffsize); 637 if (buffer == NULL) 638 { 639 printf ("src_fopen: failed to alloc memory)\n"); 640 exit (1); 641 } 642 sprintf (buffer, "%s/%s", srcdir, filename); 643 f = fopen (buffer, mode); 644 tests_free (buffer, buffsize); 645 return f; 646#endif 647} 648 649void 650set_emin (mpfr_exp_t exponent) 651{ 652 if (mpfr_set_emin (exponent)) 653 { 654 printf ("set_emin: setting emin to %ld failed\n", (long int) exponent); 655 exit (1); 656 } 657} 658 659void 660set_emax (mpfr_exp_t exponent) 661{ 662 if (mpfr_set_emax (exponent)) 663 { 664 printf ("set_emax: setting emax to %ld failed\n", (long int) exponent); 665 exit (1); 666 } 667} 668 669/* pos is 512 times the proportion of negative numbers. 670 If pos=256, half of the numbers are negative. 671 If pos=0, all generated numbers are positive. 672*/ 673void 674tests_default_random (mpfr_ptr x, int pos, mpfr_exp_t emin, mpfr_exp_t emax, 675 int always_scale) 676{ 677 MPFR_ASSERTN (emin <= emax); 678 MPFR_ASSERTN (emin >= MPFR_EMIN_MIN); 679 MPFR_ASSERTN (emax <= MPFR_EMAX_MAX); 680 /* but it isn't required that emin and emax are in the current 681 exponent range (see below), so that underflow/overflow checks 682 can be done on 64-bit machines without a manual change of the 683 exponent range (well, this is a bit ugly...). */ 684 685 mpfr_urandomb (x, RANDS); 686 if (MPFR_IS_PURE_FP (x) && (emin >= 1 || always_scale || (randlimb () & 1))) 687 { 688 mpfr_exp_t e; 689 e = emin + (mpfr_exp_t) (randlimb () % (emax - emin + 1)); 690 /* Note: There should be no overflow here because both terms are 691 between MPFR_EMIN_MIN and MPFR_EMAX_MAX. */ 692 MPFR_ASSERTD (e >= emin && e <= emax); 693 if (mpfr_set_exp (x, e)) 694 { 695 /* The random number doesn't fit in the current exponent range. 696 In this case, test the function in the extended exponent range, 697 which should be restored by the caller. */ 698 mpfr_set_emin (MPFR_EMIN_MIN); 699 mpfr_set_emax (MPFR_EMAX_MAX); 700 mpfr_set_exp (x, e); 701 } 702 } 703 if (randlimb () % 512 < pos) 704 mpfr_neg (x, x, MPFR_RNDN); 705} 706 707/* If sb = -1, then the function is tested in only one rounding mode 708 (the one provided in rnd) and the ternary value is not checked. 709 Otherwise, the function is tested in the 5 rounding modes, rnd must 710 initially be MPFR_RNDZ, y = RNDZ(f(x)), and sb is 0 if f(x) is exact, 711 1 if f(x) is inexact (in which case, y must be a regular number, 712 i.e. not the result of an overflow or an underflow); the successive 713 rounding modes are: 714 * MPFR_RNDZ, MPFR_RNDD, MPFR_RNDA, MPFR_RNDU, MPFR_RNDN for positive y; 715 * MPFR_RNDZ, MPFR_RNDU, MPFR_RNDA, MPFR_RNDD, MPFR_RNDN for negative y; 716 for the last test MPFR_RNDN, the target precision is decreased by 1 in 717 order to be able to deduce the result (anyway, for a hard-to-round case 718 in directed rounding modes, if yprec is chosen to be minimum precision 719 preserving this hard-to-round case, then one has a hard-to-round case 720 in round-to-nearest for precision yprec-1). If the target precision was 721 MPFR_PREC_MIN, then skip the MPFR_RNDN test; thus to test exact special 722 cases, use a target precision larger than MPFR_PREC_MIN. 723 Note: if y is a regular number, sb corresponds to the sticky bit when 724 considering round-to-nearest with precision yprec-1. 725 As examples of use, see the calls to test5rm from the data_check and 726 bad_cases functions. */ 727static void 728test5rm (int (*fct) (FLIST), mpfr_srcptr x, mpfr_ptr y, mpfr_ptr z, 729 mpfr_rnd_t rnd, int sb, const char *name) 730{ 731 mpfr_prec_t yprec = MPFR_PREC (y); 732 mpfr_rnd_t rndnext = MPFR_RND_MAX; /* means uninitialized */ 733 int expected_inex = INT_MIN; 734 735 MPFR_ASSERTN (sb == -1 || rnd == MPFR_RNDZ); 736 mpfr_set_prec (z, yprec); 737 while (1) 738 { 739 int inex; 740 741 MPFR_ASSERTN (rnd != MPFR_RND_MAX); 742 inex = fct (z, x, rnd); 743 if (sb == -1) 744 expected_inex = inex; /* not checked */ 745 else if (rnd != MPFR_RNDN) 746 expected_inex = 747 sb == 0 ? 0 : MPFR_IS_LIKE_RNDD (rnd, MPFR_SIGN (y)) ? -1 : 1; 748 MPFR_ASSERTN (expected_inex != INT_MIN); 749 if (!(SAME_VAL (y, z) || SAME_SIGN (inex, expected_inex))) 750 { 751 printf ("Error for %s with xprec=%lu, yprec=%lu, rnd=%s\nx = ", 752 name, (unsigned long) MPFR_PREC (x), (unsigned long) yprec, 753 mpfr_print_rnd_mode (rnd)); 754 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 755 printf ("\nexpected "); 756 mpfr_out_str (stdout, 16, 0, y, MPFR_RNDN); 757 printf ("\ngot "); 758 mpfr_out_str (stdout, 16, 0, z, MPFR_RNDN); 759 printf ("\n"); 760 if (sb != -1) 761 printf ("Expected inex = %d, got %d\n", expected_inex, inex); 762 exit (1); 763 } 764 765 if (sb == -1 || rnd == MPFR_RNDN) 766 break; 767 else if (rnd == MPFR_RNDZ) 768 { 769 rnd = MPFR_IS_NEG (y) ? MPFR_RNDU : MPFR_RNDD; 770 rndnext = MPFR_RNDA; 771 } 772 else 773 { 774 rnd = rndnext; 775 if (rnd == MPFR_RNDA) 776 { 777 if (sb) 778 mpfr_nexttoinf (y); 779 rndnext = MPFR_IS_NEG (y) ? MPFR_RNDD : MPFR_RNDU; 780 } 781 else if (rndnext != MPFR_RNDN) 782 rndnext = MPFR_RNDN; 783 else 784 { 785 if (yprec == MPFR_PREC_MIN) 786 break; 787 /* If sb = 1, then mpfr_nexttoinf was called on y for the 788 MPFR_RNDA test, i.e. y = RNDA(yprec,f(x)); we use MPFR_RNDZ 789 since one has the property RNDN(p,w) = RNDZ(p,RNDA(p+1,w)) 790 when w is not a midpoint in precision p. If sb = 0, then 791 y = f(x), so that RNDN(yprec-1,f(x)) = RNDN(yprec-1,y). */ 792 inex = mpfr_prec_round (y, --yprec, sb ? MPFR_RNDZ : MPFR_RNDN); 793 expected_inex = sb ? 794 MPFR_SIGN (y) * (inex == 0 ? 1 : -1) : inex; 795 mpfr_set_prec (z, yprec); 796 } 797 } 798 } 799} 800 801/* Check data in file f for function foo, with name 'name'. 802 Each line consists of the file f one: 803 804 xprec yprec rnd x y 805 806 where: 807 808 xprec is the input precision 809 yprec is the output precision 810 rnd is the rounding mode (n, z, u, d, a, Z, *) 811 x is the input (hexadecimal format) 812 y is the expected output (hexadecimal format) for foo(x) with rounding rnd 813 814 If rnd is Z, then y is the expected output in round-toward-zero and 815 it is assumed to be inexact; the four directed rounding modes are 816 tested, and the round-to-nearest mode is tested in precision yprec-1. 817 See details in the description of test5rm above. 818 819 If rnd is *, y must be an exact case (possibly a special case). 820 All the rounding modes are tested and the ternary value is checked 821 (it must be 0). 822 */ 823void 824data_check (const char *f, int (*foo) (FLIST), const char *name) 825{ 826 FILE *fp; 827 long int xprec, yprec; /* not mpfr_prec_t because of the fscanf */ 828 mpfr_t x, y, z; 829 mpfr_rnd_t rnd; 830 char r; 831 int c; 832 833 fp = fopen (f, "r"); 834 if (fp == NULL) 835 fp = src_fopen (f, "r"); 836 if (fp == NULL) 837 { 838 char *v = (char *) MPFR_VERSION_STRING; 839 840 /* In the '-dev' versions, assume that the data file exists and 841 return an error if the file cannot be opened to make sure 842 that such failures are detected. */ 843 while (*v != '\0') 844 v++; 845 if (v[-4] == '-' && v[-3] == 'd' && v[-2] == 'e' && v[-1] == 'v') 846 { 847 printf ("Error: unable to open file '%s'\n", f); 848 exit (1); 849 } 850 else 851 return; 852 } 853 854 mpfr_init (x); 855 mpfr_init (y); 856 mpfr_init (z); 857 858 while (!feof (fp)) 859 { 860 /* skip whitespace, for consistency */ 861 if (fscanf (fp, " ") == EOF) 862 { 863 if (ferror (fp)) 864 { 865 perror ("data_check"); 866 exit (1); 867 } 868 else 869 break; /* end of file */ 870 } 871 872 if ((c = getc (fp)) == EOF) 873 { 874 if (ferror (fp)) 875 { 876 perror ("data_check"); 877 exit (1); 878 } 879 else 880 break; /* end of file */ 881 } 882 883 if (c == '#') /* comment: read entire line */ 884 { 885 do 886 { 887 c = getc (fp); 888 } 889 while (!feof (fp) && c != '\n'); 890 } 891 else 892 { 893 ungetc (c, fp); 894 895 c = fscanf (fp, "%ld %ld %c", &xprec, &yprec, &r); 896 MPFR_ASSERTN (MPFR_PREC_COND (xprec)); 897 MPFR_ASSERTN (MPFR_PREC_COND (yprec)); 898 if (c == EOF) 899 { 900 perror ("data_check"); 901 exit (1); 902 } 903 else if (c != 3) 904 { 905 printf ("Error: corrupted line in file '%s'\n", f); 906 exit (1); 907 } 908 909 switch (r) 910 { 911 case 'n': 912 rnd = MPFR_RNDN; 913 break; 914 case 'z': case 'Z': 915 rnd = MPFR_RNDZ; 916 break; 917 case 'u': 918 rnd = MPFR_RNDU; 919 break; 920 case 'd': 921 rnd = MPFR_RNDD; 922 break; 923 case '*': 924 rnd = MPFR_RND_MAX; /* non-existing rounding mode */ 925 break; 926 default: 927 printf ("Error: unexpected rounding mode" 928 " in file '%s': %c\n", f, (int) r); 929 exit (1); 930 } 931 mpfr_set_prec (x, xprec); 932 mpfr_set_prec (y, yprec); 933 if (mpfr_inp_str (x, fp, 0, MPFR_RNDN) == 0) 934 { 935 printf ("Error: corrupted argument in file '%s'\n", f); 936 exit (1); 937 } 938 if (mpfr_inp_str (y, fp, 0, MPFR_RNDN) == 0) 939 { 940 printf ("Error: corrupted result in file '%s'\n", f); 941 exit (1); 942 } 943 if (getc (fp) != '\n') 944 { 945 printf ("Error: result not followed by \\n in file '%s'\n", f); 946 exit (1); 947 } 948 /* Skip whitespace, in particular at the end of the file. */ 949 if (fscanf (fp, " ") == EOF && ferror (fp)) 950 { 951 perror ("data_check"); 952 exit (1); 953 } 954 if (r == '*') 955 test5rm (foo, x, y, z, MPFR_RNDZ, 0, name); 956 else 957 test5rm (foo, x, y, z, rnd, r == 'Z' ? 1 : -1, name); 958 } 959 } 960 961 mpfr_clear (x); 962 mpfr_clear (y); 963 mpfr_clear (z); 964 965 fclose (fp); 966} 967 968/* Test n random bad cases. A precision py in [pymin,pymax] and 969 * a number y of precision py are chosen randomly. One computes 970 * x = inv(y) in precision px = py + psup (rounded to nearest). 971 * Then (in general), y is a bad case for fct in precision py (in 972 * the directed rounding modes, but also in the rounding-to-nearest 973 * mode for some lower precision: see data_check). 974 * fct, inv, name: data related to the function. 975 * pos, emin, emax: arguments for tests_default_random. 976 * For debugging purpose (e.g. in case of crash or infinite loop), 977 * you can set the MPFR_DEBUG_BADCASES environment variable to 1 in 978 * order to output information about the tested worst cases. You can 979 * also enable logging (when supported), but this may give too much 980 * information. 981 */ 982void 983bad_cases (int (*fct)(FLIST), int (*inv)(FLIST), const char *name, 984 int pos, mpfr_exp_t emin, mpfr_exp_t emax, 985 mpfr_prec_t pymin, mpfr_prec_t pymax, mpfr_prec_t psup, 986 int n) 987{ 988 mpfr_t x, y, z; 989 char *dbgenv; 990 int cnt = 0, i, dbg; 991 mpfr_exp_t old_emin, old_emax; 992 993 old_emin = mpfr_get_emin (); 994 old_emax = mpfr_get_emax (); 995 996 dbgenv = getenv ("MPFR_DEBUG_BADCASES"); 997 dbg = dbgenv != 0 ? atoi (dbgenv) : 0; /* debug level */ 998 mpfr_inits2 (MPFR_PREC_MIN, x, y, z, (mpfr_ptr) 0); 999 for (i = 0; i < n; i++) 1000 { 1001 mpfr_prec_t px, py, pz; 1002 int inex_inv, inex, sb; 1003 1004 if (dbg) 1005 printf ("bad_cases: %s, i = %d\n", name, i); 1006 py = pymin + (randlimb () % (pymax - pymin + 1)); 1007 mpfr_set_prec (y, py); 1008 tests_default_random (y, pos, emin, emax, 0); 1009 if (dbg) 1010 { 1011 printf ("bad_cases: yprec =%4ld, y = ", (long) py); 1012 mpfr_out_str (stdout, 16, 0, y, MPFR_RNDN); 1013 printf ("\n"); 1014 } 1015 px = py + psup; 1016 mpfr_set_prec (x, px); 1017 if (dbg) 1018 printf ("bad_cases: xprec =%4ld\n", (long) px); 1019 mpfr_clear_flags (); 1020 inex_inv = inv (x, y, MPFR_RNDN); 1021 if (mpfr_nanflag_p () || mpfr_overflow_p () || mpfr_underflow_p ()) 1022 { 1023 if (dbg) 1024 printf ("bad_cases: no normal inverse\n"); 1025 goto next_i; 1026 } 1027 if (dbg > 1) 1028 { 1029 printf ("bad_cases: x = "); 1030 mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); 1031 printf ("\n"); 1032 } 1033 pz = px; 1034 do 1035 { 1036 pz += 32; 1037 mpfr_set_prec (z, pz); 1038 sb = fct (z, x, MPFR_RNDN) != 0; 1039 if (!sb) 1040 { 1041 if (dbg) 1042 printf ("bad_cases: exact case\n"); 1043 if (inex_inv) 1044 { 1045 printf ("bad_cases: f exact while f^(-1) inexact,\n" 1046 "due to a poor choice of the parameters.\n"); 1047 exit (1); 1048 /* alternatively, goto next_i */ 1049 } 1050 inex = 0; 1051 break; 1052 } 1053 if (dbg) 1054 { 1055 if (dbg > 1) 1056 { 1057 printf ("bad_cases: %s(x) ~= ", name); 1058 mpfr_out_str (stdout, 16, 0, z, MPFR_RNDN); 1059 } 1060 else 1061 { 1062 printf ("bad_cases: [MPFR_RNDZ] ~= "); 1063 mpfr_out_str (stdout, 16, 40, z, MPFR_RNDZ); 1064 } 1065 printf ("\n"); 1066 } 1067 inex = mpfr_prec_round (z, py, MPFR_RNDN); 1068 if (mpfr_nanflag_p () || mpfr_overflow_p () || mpfr_underflow_p () 1069 || ! mpfr_equal_p (z, y)) 1070 { 1071 if (dbg) 1072 { 1073 printf ("bad_cases: inverse doesn't match for %s\ny = ", 1074 name); 1075 mpfr_out_str (stdout, 16, 0, y, MPFR_RNDN); 1076 printf ("\nz = "); 1077 mpfr_out_str (stdout, 16, 0, z, MPFR_RNDN); 1078 printf ("\n"); 1079 } 1080 goto next_i; 1081 } 1082 } 1083 while (inex == 0); 1084 /* We really have a bad case. */ 1085 do 1086 py--; 1087 while (py >= MPFR_PREC_MIN && mpfr_prec_round (z, py, MPFR_RNDZ) == 0); 1088 py++; 1089 /* py is now the smallest output precision such that we have 1090 a bad case in the directed rounding modes. */ 1091 if (mpfr_prec_round (y, py, MPFR_RNDZ) != 0) 1092 { 1093 printf ("Internal error for i = %d\n", i); 1094 exit (1); 1095 } 1096 if ((inex > 0 && MPFR_IS_POS (z)) || 1097 (inex < 0 && MPFR_IS_NEG (z))) 1098 { 1099 mpfr_nexttozero (y); 1100 if (mpfr_zero_p (y)) 1101 goto next_i; 1102 } 1103 if (dbg) 1104 { 1105 printf ("bad_cases: yprec =%4ld, y = ", (long) py); 1106 mpfr_out_str (stdout, 16, 0, y, MPFR_RNDN); 1107 printf ("\n"); 1108 } 1109 /* Note: y is now the expected result rounded toward zero. */ 1110 test5rm (fct, x, y, z, MPFR_RNDZ, sb, name); 1111 cnt++; 1112 next_i: 1113 /* In case the exponent range has been changed by 1114 tests_default_random()... */ 1115 mpfr_set_emin (old_emin); 1116 mpfr_set_emax (old_emax); 1117 } 1118 mpfr_clears (x, y, z, (mpfr_ptr) 0); 1119 1120 if (dbg) 1121 printf ("bad_cases: %d bad cases over %d generated values\n", cnt, n); 1122 1123 if (getenv ("MPFR_CHECK_BADCASES") && n - cnt > n/10) 1124 { 1125 printf ("bad_cases: too few bad cases (%d over %d generated values)\n", 1126 cnt, n); 1127 exit (1); 1128 } 1129} 1130 1131void 1132flags_out (unsigned int flags) 1133{ 1134 int none = 1; 1135 1136 if (flags & MPFR_FLAGS_UNDERFLOW) 1137 none = 0, printf (" underflow"); 1138 if (flags & MPFR_FLAGS_OVERFLOW) 1139 none = 0, printf (" overflow"); 1140 if (flags & MPFR_FLAGS_NAN) 1141 none = 0, printf (" nan"); 1142 if (flags & MPFR_FLAGS_INEXACT) 1143 none = 0, printf (" inexact"); 1144 if (flags & MPFR_FLAGS_ERANGE) 1145 none = 0, printf (" erange"); 1146 if (none) 1147 printf (" none"); 1148 printf (" (%u)\n", flags); 1149} 1150 1151static void 1152abort_called (int x) 1153{ 1154 /* Ok, abort has been called */ 1155 exit (0); 1156} 1157 1158/* This function has to be called for a test 1159 that will call the abort function */ 1160void 1161tests_expect_abort (void) 1162{ 1163#if defined(HAVE_SIGACTION) 1164 struct sigaction act; 1165 int ret; 1166 1167 memset (&act, 0, sizeof act); 1168 act.sa_handler = abort_called; 1169 ret = sigaction (SIGABRT, &act, NULL); 1170 if (ret != 0) 1171 { 1172 /* Can't register error handler: Skip test */ 1173 exit (77); 1174 } 1175#elif defined(HAVE_SIGNAL) 1176 signal (SIGABRT, abort_called); 1177#else 1178 /* Can't register error handler: Skip test */ 1179 exit (77); 1180#endif 1181} 1182 1183/* Guess whether the test runs within Valgrind. */ 1184int 1185tests_run_within_valgrind (void) 1186{ 1187 char *p; 1188 1189 p = getenv ("LD_PRELOAD"); 1190 if (p == NULL) 1191 return 0; 1192 return (strstr (p, "/valgrind/") != NULL || 1193 strstr (p, "/vgpreload") != NULL); 1194} 1195