math/libm-test-support.c

   1 /* Support code for testing libm functions (compiled once per type).
   2    Copyright (C) 1997-2018 Free Software Foundation, Inc.
   3    This file is part of the GNU C Library.
   4
   5    The GNU C Library is free software; you can redistribute it and/or
   6    modify it under the terms of the GNU Lesser General Public
   7    License as published by the Free Software Foundation; either
   8    version 2.1 of the License, or (at your option) any later version.
   9
  10    The GNU C Library is distributed in the hope that it will be useful,
  11    but WITHOUT ANY WARRANTY; without even the implied warranty of
  12    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13    Lesser General Public License for more details.
  14
  15    You should have received a copy of the GNU Lesser General Public
  16    License along with the GNU C Library; if not, see
  17    <http://www.gnu.org/licenses/>.  */
  18
  19 /* Part of testsuite for libm.
  20
  21    libm-test-support.c contains functions shared by tests of different
  22    libm functions and types; it is compiled once per type.
  23    libm-test-driver.c defines the main function, and various variables
  24    that are used to configure the code in libm-test-support.c for
  25    different types and for variants such as testing inline functions.
  26
  27    The tests of individual functions are in .inc files processed by
  28    gen-libm-test.pl, with the resulting files included together with
  29    libm-test-driver.c.
  30
  31    The per-type headers included both before libm-test-support.c and
  32    for the tests of individual functions must define the following
  33    macros:
  34
  35    FUNC(function): Convert general function name (like cos) to name
  36    with correct suffix (e.g. cosl or cosf).
  37
  38    FLOAT: Floating-point type to test.
  39
  40    BUILD_COMPLEX(real, imag): Create a complex number by calling a
  41    macro such as CMPLX.
  42
  43    PREFIX: The prefix for <float.h> macros for the type (e.g. LDBL,
  44    DBL, or FLT).
  45
  46    TYPE_STR: The name of the type as used in ulps files, as a string.
  47
  48    ULP_IDX, ULP_I_IDX: The array indexes for ulps values for this
  49    function.
  50
  51    LIT: Append the correct suffix to a literal.
  52
  53    LITM: Append the correct suffix to an M_* macro name.
  54
  55    FTOSTR: A function similar in type to strfromf which converts a
  56    FLOAT to a string.
  57
  58    snan_value_MACRO: The macro such as SNAN for a signaling NaN for
  59    the type.
  60
  61 */
  62
  63 /* Parameter handling is primitive in the moment:
  64    --verbose=[0..3] for different levels of output:
  65    0: only error count
  66    1: basic report on failed tests (default)
  67    2: full report on all tests
  68    -v for full output (equals --verbose=3)
  69    -u for generation of an ULPs file
  70  */
  71
  72 /* "Philosophy":
  73
  74    This suite tests some aspects of the correct implementation of
  75    mathematical functions in libm.  Some simple, specific parameters
  76    are tested for correctness but there's no exhaustive
  77    testing.  Handling of specific inputs (e.g. infinity, not-a-number)
  78    is also tested.  Correct handling of exceptions is checked
  79    against.  These implemented tests should check all cases that are
  80    specified in ISO C99.
  81
  82    NaN values: The payload of NaNs is set in inputs for functions
  83    where it is significant, and is examined in the outputs of some
  84    functions.
  85
  86    Inline functions: Inlining functions should give an improvement in
  87    speed - but not in precission.  The inlined functions return
  88    reasonable values for a reasonable range of input values.  The
  89    result is not necessarily correct for all values and exceptions are
  90    not correctly raised in all cases.  Problematic input and return
  91    values are infinity, not-a-number and minus zero.  This suite
  92    therefore does not check these specific inputs and the exception
  93    handling for inlined mathematical functions - just the "reasonable"
  94    values are checked.
  95
  96    Beware: The tests might fail for any of the following reasons:
  97    - Tests are wrong
  98    - Functions are wrong
  99    - Floating Point Unit not working properly
 100    - Compiler has errors
 101
 102    With e.g. gcc 2.7.2.2 the test for cexp fails because of a compiler error.
 103
 104
 105    To Do: All parameter should be numbers that can be represented as
 106    exact floating point values.  Currently some values cannot be
 107    represented exactly and therefore the result is not the expected
 108    result.  For this we will use 36 digits so that numbers can be
 109    represented exactly.  */
 110
 111 #include "libm-test-support.h"
 112
 113 #include <argp.h>
 114 #include <errno.h>
 115 #include <string.h>
 116
 117 /* This header defines func_ulps, func_real_ulps and func_imag_ulps
 118    arrays.  */
 119 #include "libm-test-ulps.h"
 120
 121 /* Maximum character buffer to store a stringitized FLOAT value.  */
 122 #define FSTR_MAX (128)
 123
 124 #define ulps_file_name "ULPs"   /* Name of the ULPs file.  */
 125 static FILE *ulps_file;         /* File to document difference.  */
 126 static int output_ulps;         /* Should ulps printed?  */
 127 static char *output_dir;        /* Directory where generated files will be written.  */
 128
 129 static int noErrors;    /* number of errors */
 130 static int noTests;     /* number of tests (without testing exceptions) */
 131 static int noExcTests;  /* number of tests for exception flags */
 132 static int noErrnoTests;/* number of tests for errno values */
 133
 134 static int verbose;
 135 static int output_max_error;    /* Should the maximal errors printed?  */
 136 static int output_points;       /* Should the single function results printed?  */
 137 static int ignore_max_ulp;      /* Should we ignore max_ulp?  */
 138 static int test_ibm128;         /* Is argument or result IBM long double?  */
 139
 140 static FLOAT max_error, real_max_error, imag_max_error;
 141
 142 static FLOAT prev_max_error, prev_real_max_error, prev_imag_max_error;
 143
 144 static FLOAT max_valid_error;
 145
 146 /* Sufficient numbers of digits to represent any floating-point value
 147    unambiguously (for any choice of the number of bits in the first
 148    hex digit, in the case of TYPE_HEX_DIG).  When used with printf
 149    formats where the precision counts only digits after the point, 1
 150    is subtracted from these values. */
 151 #define TYPE_DECIMAL_DIG __CONCATX (PREFIX, _DECIMAL_DIG)
 152 #define TYPE_HEX_DIG ((MANT_DIG + 6) / 4)
 153
 154 /* Converts VALUE (a floating-point number) to string and writes it to DEST.
 155    PRECISION specifies the number of fractional digits that should be printed.
 156    CONVERSION is the conversion specifier, such as in printf, e.g. 'f' or 'a'.
 157    The output is prepended with an empty space if VALUE is non-negative.  */
 158 static void
 159 fmt_ftostr (char *dest, size_t size, int precision, const char *conversion,
 160             FLOAT value)
 161 {
 162   char format[64];
 163   char *ptr_format;
 164   int ret;
 165
 166   /* Generate the format string.  */
 167   ptr_format = stpcpy (format, "%.");
 168   ret = sprintf (ptr_format, "%d", precision);
 169   ptr_format += ret;
 170   ptr_format = stpcpy (ptr_format, conversion);
 171
 172   /* Add a space to the beginning of the output string, if the floating-point
 173      number is non-negative.  This mimics the behavior of the space (' ') flag
 174      in snprintf, which is not available on strfrom.  */
 175   if (! signbit (value))
 176     {
 177       *dest = ' ';
 178       dest++;
 179       size--;
 180     }
 181
 182   /* Call the float to string conversion function, e.g.: strfromd.  */
 183   FTOSTR (dest, size, format, value);
 184 }
 185
 186 /* Compare KEY (a string, with the name of a function) with ULP (a
 187    pointer to a struct ulp_data structure), returning a value less
 188    than, equal to or greater than zero for use in bsearch.  */
 189
 190 static int
 191 compare_ulp_data (const void *key, const void *ulp)
 192 {
 193   const char *keystr = key;
 194   const struct ulp_data *ulpdat = ulp;
 195   return strcmp (keystr, ulpdat->name);
 196 }
 197
 198 static const int ulp_i_idx = ULP_I_IDX;
 199 static const int ulp_idx = ULP_IDX;
 200
 201 /* Return the ulps for NAME in array DATA with NMEMB elements, or 0 if
 202    no ulps listed.  */
 203
 204 static FLOAT
 205 find_ulps (const char *name, const struct ulp_data *data, size_t nmemb)
 206 {
 207   const struct ulp_data *entry = bsearch (name, data, nmemb, sizeof (*data),
 208                                           compare_ulp_data);
 209   if (entry == NULL)
 210     return 0;
 211   else
 212     return entry->max_ulp[(flag_test_inline ? ulp_i_idx : ulp_idx)];
 213 }
 214
 215 void
 216 init_max_error (const char *name, int exact, int testing_ibm128)
 217 {
 218   max_error = 0;
 219   real_max_error = 0;
 220   imag_max_error = 0;
 221   test_ibm128 = testing_ibm128;
 222   prev_max_error = find_ulps (name, func_ulps,
 223                               sizeof (func_ulps) / sizeof (func_ulps[0]));
 224   prev_real_max_error = find_ulps (name, func_real_ulps,
 225                                    (sizeof (func_real_ulps)
 226                                     / sizeof (func_real_ulps[0])));
 227   prev_imag_max_error = find_ulps (name, func_imag_ulps,
 228                                    (sizeof (func_imag_ulps)
 229                                     / sizeof (func_imag_ulps[0])));
 230   if (testing_ibm128)
 231     /* The documented accuracy of IBM long double division is 3ulp
 232        (see libgcc/config/rs6000/ibm-ldouble-format), so do not
 233        require better accuracy for libm functions that are exactly
 234        defined for other formats.  */
 235     max_valid_error = exact ? 3 : 16;
 236   else
 237     max_valid_error = exact ? 0 : 9;
 238   prev_max_error = (prev_max_error <= max_valid_error
 239                     ? prev_max_error
 240                     : max_valid_error);
 241   prev_real_max_error = (prev_real_max_error <= max_valid_error
 242                          ? prev_real_max_error
 243                          : max_valid_error);
 244   prev_imag_max_error = (prev_imag_max_error <= max_valid_error
 245                          ? prev_imag_max_error
 246                          : max_valid_error);
 247   feclearexcept (FE_ALL_EXCEPT);
 248   errno = 0;
 249 }
 250
 251 static void
 252 set_max_error (FLOAT current, FLOAT *curr_max_error)
 253 {
 254   if (current > *curr_max_error && current <= max_valid_error)
 255     *curr_max_error = current;
 256 }
 257
 258
 259 /* Print a FLOAT.  */
 260 static void
 261 print_float (FLOAT f)
 262 {
 263   /* As printf doesn't differ between a sNaN and a qNaN, do this manually.  */
 264   if (issignaling (f))
 265     printf ("sNaN\n");
 266   else if (isnan (f))
 267     printf ("qNaN\n");
 268   else
 269     {
 270       char fstrn[FSTR_MAX], fstrx[FSTR_MAX];
 271       fmt_ftostr (fstrn, FSTR_MAX, TYPE_DECIMAL_DIG - 1, "e", f);
 272       fmt_ftostr (fstrx, FSTR_MAX, TYPE_HEX_DIG - 1, "a", f);
 273       printf ("%s  %s\n", fstrn, fstrx);
 274     }
 275 }
 276
 277 /* Should the message print to screen?  This depends on the verbose flag,
 278    and the test status.  */
 279 static int
 280 print_screen (int ok)
 281 {
 282   if (output_points
 283       && (verbose > 1
 284           || (verbose == 1 && ok == 0)))
 285     return 1;
 286   return 0;
 287 }
 288
 289
 290 /* Should the message print to screen?  This depends on the verbose flag,
 291    and the test status.  */
 292 static int
 293 print_screen_max_error (int ok)
 294 {
 295   if (output_max_error
 296       && (verbose > 1
 297           || ((verbose == 1) && (ok == 0))))
 298     return 1;
 299   return 0;
 300 }
 301
 302 /* Update statistic counters.  */
 303 static void
 304 update_stats (int ok)
 305 {
 306   ++noTests;
 307   if (!ok)
 308     ++noErrors;
 309 }
 310
 311 static void
 312 print_function_ulps (const char *function_name, FLOAT ulp)
 313 {
 314   if (output_ulps)
 315     {
 316       char ustrn[FSTR_MAX];
 317       FTOSTR (ustrn, FSTR_MAX, "%.0f", FUNC (ceil) (ulp));
 318       fprintf (ulps_file, "Function: \"%s\":\n", function_name);
 319       fprintf (ulps_file, "%s: %s\n", qtype_str, ustrn);
 320     }
 321 }
 322
 323
 324 static void
 325 print_complex_function_ulps (const char *function_name, FLOAT real_ulp,
 326                              FLOAT imag_ulp)
 327 {
 328   if (output_ulps)
 329     {
 330       char fstrn[FSTR_MAX];
 331       if (real_ulp != 0.0)
 332         {
 333           FTOSTR (fstrn, FSTR_MAX, "%.0f", FUNC (ceil) (real_ulp));
 334           fprintf (ulps_file, "Function: Real part of \"%s\":\n", function_name);
 335           fprintf (ulps_file, "%s: %s\n", qtype_str, fstrn);
 336         }
 337       if (imag_ulp != 0.0)
 338         {
 339           FTOSTR (fstrn, FSTR_MAX, "%.0f", FUNC (ceil) (imag_ulp));
 340           fprintf (ulps_file, "Function: Imaginary part of \"%s\":\n", function_name);
 341           fprintf (ulps_file, "%s: %s\n", qtype_str, fstrn);
 342         }
 343
 344
 345     }
 346 }
 347
 348
 349
 350 /* Test if Floating-Point stack hasn't changed */
 351 static void
 352 fpstack_test (const char *test_name)
 353 {
 354 #if defined (__i386__) || defined (__x86_64__)
 355   static int old_stack;
 356   int sw;
 357
 358   asm ("fnstsw" : "=a" (sw));
 359   sw >>= 11;
 360   sw &= 7;
 361
 362   if (sw != old_stack)
 363     {
 364       printf ("FP-Stack wrong after test %s (%d, should be %d)\n",
 365               test_name, sw, old_stack);
 366       ++noErrors;
 367       old_stack = sw;
 368     }
 369 #endif
 370 }
 371
 372
 373 void
 374 print_max_error (const char *func_name)
 375 {
 376   int ok = 0;
 377
 378   if (max_error == 0.0 || (max_error <= prev_max_error && !ignore_max_ulp))
 379     {
 380       ok = 1;
 381     }
 382
 383   if (!ok)
 384     print_function_ulps (func_name, max_error);
 385
 386
 387   if (print_screen_max_error (ok))
 388     {
 389       char mestr[FSTR_MAX], pmestr[FSTR_MAX];
 390       FTOSTR (mestr, FSTR_MAX, "%.0f", FUNC (ceil) (max_error));
 391       FTOSTR (pmestr, FSTR_MAX, "%.0f", FUNC (ceil) (prev_max_error));
 392       printf ("Maximal error of `%s'\n", func_name);
 393       printf (" is      : %s ulp\n", mestr);
 394       printf (" accepted: %s ulp\n", pmestr);
 395     }
 396
 397   update_stats (ok);
 398 }
 399
 400
 401 void
 402 print_complex_max_error (const char *func_name)
 403 {
 404   int real_ok = 0, imag_ok = 0, ok;
 405
 406   if (real_max_error == 0
 407       || (real_max_error <= prev_real_max_error && !ignore_max_ulp))
 408     {
 409       real_ok = 1;
 410     }
 411
 412   if (imag_max_error == 0
 413       || (imag_max_error <= prev_imag_max_error && !ignore_max_ulp))
 414     {
 415       imag_ok = 1;
 416     }
 417
 418   ok = real_ok && imag_ok;
 419
 420   if (!ok)
 421     print_complex_function_ulps (func_name,
 422                                  real_ok ? 0 : real_max_error,
 423                                  imag_ok ? 0 : imag_max_error);
 424
 425   if (print_screen_max_error (ok))
 426     {
 427       char rmestr[FSTR_MAX], prmestr[FSTR_MAX];
 428       char imestr[FSTR_MAX], pimestr[FSTR_MAX];
 429       FTOSTR (rmestr, FSTR_MAX, "%.0f", FUNC (ceil) (real_max_error));
 430       FTOSTR (prmestr, FSTR_MAX, "%.0f", FUNC (ceil) (prev_real_max_error));
 431       FTOSTR (imestr, FSTR_MAX, "%.0f", FUNC (ceil) (imag_max_error));
 432       FTOSTR (pimestr, FSTR_MAX, "%.0f", FUNC (ceil) (prev_imag_max_error));
 433       printf ("Maximal error of real part of: %s\n", func_name);
 434       printf (" is      : %s ulp\n", rmestr);
 435       printf (" accepted: %s ulp\n", prmestr);
 436       printf ("Maximal error of imaginary part of: %s\n", func_name);
 437       printf (" is      : %s ulp\n", imestr);
 438       printf (" accepted: %s ulp\n", pimestr);
 439     }
 440
 441   update_stats (ok);
 442 }
 443
 444
 445 #if FE_ALL_EXCEPT
 446 /* Test whether a given exception was raised.  */
 447 static void
 448 test_single_exception (const char *test_name,
 449                        int exception,
 450                        int exc_flag,
 451                        int fe_flag,
 452                        const char *flag_name)
 453 {
 454   int ok = 1;
 455   if (exception & exc_flag)
 456     {
 457       if (fetestexcept (fe_flag))
 458         {
 459           if (print_screen (1))
 460             printf ("Pass: %s: Exception \"%s\" set\n", test_name, flag_name);
 461         }
 462       else
 463         {
 464           ok = 0;
 465           if (print_screen (0))
 466             printf ("Failure: %s: Exception \"%s\" not set\n",
 467                     test_name, flag_name);
 468         }
 469     }
 470   else
 471     {
 472       if (fetestexcept (fe_flag))
 473         {
 474           ok = 0;
 475           if (print_screen (0))
 476             printf ("Failure: %s: Exception \"%s\" set\n",
 477                     test_name, flag_name);
 478         }
 479       else
 480         {
 481           if (print_screen (1))
 482             printf ("%s: Exception \"%s\" not set\n", test_name,
 483                     flag_name);
 484         }
 485     }
 486   if (!ok)
 487     ++noErrors;
 488 }
 489 #endif
 490
 491 /* Test whether exceptions given by EXCEPTION are raised.  Ignore thereby
 492    allowed but not required exceptions.
 493 */
 494 static void
 495 test_exceptions (const char *test_name, int exception)
 496 {
 497   if (flag_test_exceptions && EXCEPTION_TESTS (FLOAT))
 498     {
 499       ++noExcTests;
 500 #ifdef FE_DIVBYZERO
 501       if ((exception & DIVIDE_BY_ZERO_EXCEPTION_OK) == 0)
 502         test_single_exception (test_name, exception,
 503                                DIVIDE_BY_ZERO_EXCEPTION, FE_DIVBYZERO,
 504                                "Divide by zero");
 505 #endif
 506 #ifdef FE_INVALID
 507       if ((exception & INVALID_EXCEPTION_OK) == 0)
 508         test_single_exception (test_name, exception,
 509                                INVALID_EXCEPTION, FE_INVALID,
 510                                "Invalid operation");
 511 #endif
 512 #ifdef FE_OVERFLOW
 513       if ((exception & OVERFLOW_EXCEPTION_OK) == 0)
 514         test_single_exception (test_name, exception, OVERFLOW_EXCEPTION,
 515                                FE_OVERFLOW, "Overflow");
 516 #endif
 517       /* Spurious "underflow" and "inexact" exceptions are always
 518          allowed for IBM long double, in line with the underlying
 519          arithmetic.  */
 520 #ifdef FE_UNDERFLOW
 521       if ((exception & UNDERFLOW_EXCEPTION_OK) == 0
 522           && !(test_ibm128
 523                && (exception & UNDERFLOW_EXCEPTION) == 0))
 524         test_single_exception (test_name, exception, UNDERFLOW_EXCEPTION,
 525                                FE_UNDERFLOW, "Underflow");
 526 #endif
 527 #ifdef FE_INEXACT
 528       if ((exception & (INEXACT_EXCEPTION | NO_INEXACT_EXCEPTION)) != 0
 529           && !(test_ibm128
 530                && (exception & NO_INEXACT_EXCEPTION) != 0))
 531         test_single_exception (test_name, exception, INEXACT_EXCEPTION,
 532                                FE_INEXACT, "Inexact");
 533 #endif
 534     }
 535   feclearexcept (FE_ALL_EXCEPT);
 536 }
 537
 538 /* Test whether errno for TEST_NAME, set to ERRNO_VALUE, has value
 539    EXPECTED_VALUE (description EXPECTED_NAME).  */
 540 static void
 541 test_single_errno (const char *test_name, int errno_value,
 542                    int expected_value, const char *expected_name)
 543 {
 544   if (errno_value == expected_value)
 545     {
 546       if (print_screen (1))
 547         printf ("Pass: %s: errno set to %d (%s)\n", test_name, errno_value,
 548                 expected_name);
 549     }
 550   else
 551     {
 552       ++noErrors;
 553       if (print_screen (0))
 554         printf ("Failure: %s: errno set to %d, expected %d (%s)\n",
 555                 test_name, errno_value, expected_value, expected_name);
 556     }
 557 }
 558
 559 /* Test whether errno (value ERRNO_VALUE) has been for TEST_NAME set
 560    as required by EXCEPTIONS.  */
 561 static void
 562 test_errno (const char *test_name, int errno_value, int exceptions)
 563 {
 564   if (flag_test_errno)
 565     {
 566       ++noErrnoTests;
 567       if (exceptions & ERRNO_UNCHANGED)
 568         test_single_errno (test_name, errno_value, 0, "unchanged");
 569       if (exceptions & ERRNO_EDOM)
 570         test_single_errno (test_name, errno_value, EDOM, "EDOM");
 571       if (exceptions & ERRNO_ERANGE)
 572         test_single_errno (test_name, errno_value, ERANGE, "ERANGE");
 573     }
 574 }
 575
 576 /* Returns the number of ulps that GIVEN is away from EXPECTED.  */
 577 #define ULPDIFF(given, expected) \
 578         (FUNC(fabs) ((given) - (expected)) / ulp (expected))
 579
 580 /* Returns the size of an ulp for VALUE.  */
 581 static FLOAT
 582 ulp (FLOAT value)
 583 {
 584   FLOAT ulp;
 585
 586   switch (fpclassify (value))
 587     {
 588       case FP_ZERO:
 589         /* We compute the distance to the next FP which is the same as the
 590            value of the smallest subnormal number. Previously we used
 591            2^-(MANT_DIG - 1) which is too large a value to be useful. Note that we
 592            can't use ilogb(0), since that isn't a valid thing to do. As a point
 593            of comparison Java's ulp returns the next normal value e.g.
 594            2^(1 - MAX_EXP) for ulp(0), but that is not what we want for
 595            glibc.  */
 596         /* Fall through...  */
 597       case FP_SUBNORMAL:
 598         /* The next closest subnormal value is a constant distance away.  */
 599         ulp = FUNC(ldexp) (1.0, MIN_EXP - MANT_DIG);
 600         break;
 601
 602       case FP_NORMAL:
 603         ulp = FUNC(ldexp) (1.0, FUNC(ilogb) (value) - MANT_DIG + 1);
 604         break;
 605
 606       default:
 607         /* It should never happen. */
 608         abort ();
 609         break;
 610     }
 611   return ulp;
 612 }
 613
 614 static void
 615 check_float_internal (const char *test_name, FLOAT computed, FLOAT expected,
 616                       int exceptions,
 617                       FLOAT *curr_max_error, FLOAT max_ulp)
 618 {
 619   int ok = 0;
 620   int print_diff = 0;
 621   FLOAT diff = 0;
 622   FLOAT ulps = 0;
 623   int errno_value = errno;
 624
 625   test_exceptions (test_name, exceptions);
 626   test_errno (test_name, errno_value, exceptions);
 627   if (exceptions & IGNORE_RESULT)
 628     goto out;
 629   if (issignaling (computed) && issignaling (expected))
 630     {
 631       if ((exceptions & TEST_NAN_SIGN) != 0
 632           && signbit (computed) != signbit (expected))
 633         {
 634           ok = 0;
 635           printf ("signaling NaN has wrong sign.\n");
 636         }
 637       else if ((exceptions & TEST_NAN_PAYLOAD) != 0
 638                && (FUNC (getpayload) (&computed)
 639                    != FUNC (getpayload) (&expected)))
 640         {
 641           ok = 0;
 642           printf ("signaling NaN has wrong payload.\n");
 643         }
 644       else
 645         ok = 1;
 646     }
 647   else if (issignaling (computed) || issignaling (expected))
 648     ok = 0;
 649   else if (isnan (computed) && isnan (expected))
 650     {
 651       if ((exceptions & TEST_NAN_SIGN) != 0
 652           && signbit (computed) != signbit (expected))
 653         {
 654           ok = 0;
 655           printf ("quiet NaN has wrong sign.\n");
 656         }
 657       else if ((exceptions & TEST_NAN_PAYLOAD) != 0
 658                && (FUNC (getpayload) (&computed)
 659                    != FUNC (getpayload) (&expected)))
 660         {
 661           ok = 0;
 662           printf ("quiet NaN has wrong payload.\n");
 663         }
 664       else
 665         ok = 1;
 666     }
 667   else if (isinf (computed) && isinf (expected))
 668     {
 669       /* Test for sign of infinities.  */
 670       if ((exceptions & IGNORE_ZERO_INF_SIGN) == 0
 671           && signbit (computed) != signbit (expected))
 672         {
 673           ok = 0;
 674           printf ("infinity has wrong sign.\n");
 675         }
 676       else
 677         ok = 1;
 678     }
 679   /* Don't calculate ULPs for infinities or any kind of NaNs.  */
 680   else if (isinf (computed) || isnan (computed)
 681            || isinf (expected) || isnan (expected))
 682     ok = 0;
 683   else
 684     {
 685       diff = FUNC(fabs) (computed - expected);
 686       ulps = ULPDIFF (computed, expected);
 687       set_max_error (ulps, curr_max_error);
 688       print_diff = 1;
 689       if ((exceptions & IGNORE_ZERO_INF_SIGN) == 0
 690           && computed == 0.0 && expected == 0.0
 691           && signbit(computed) != signbit (expected))
 692         ok = 0;
 693       else if (ulps <= max_ulp && !ignore_max_ulp)
 694         ok = 1;
 695       else
 696         ok = 0;
 697     }
 698   if (print_screen (ok))
 699     {
 700       if (!ok)
 701         printf ("Failure: ");
 702       printf ("Test: %s\n", test_name);
 703       printf ("Result:\n");
 704       printf (" is:         ");
 705       print_float (computed);
 706       printf (" should be:  ");
 707       print_float (expected);
 708       if (print_diff)
 709         {
 710           char dstrn[FSTR_MAX], dstrx[FSTR_MAX];
 711           char ustrn[FSTR_MAX], mustrn[FSTR_MAX];
 712           fmt_ftostr (dstrn, FSTR_MAX, TYPE_DECIMAL_DIG - 1, "e", diff);
 713           fmt_ftostr (dstrx, FSTR_MAX, TYPE_HEX_DIG - 1, "a", diff);
 714           fmt_ftostr (ustrn, FSTR_MAX, 4, "f", ulps);
 715           fmt_ftostr (mustrn, FSTR_MAX, 4, "f", max_ulp);
 716           printf (" difference: %s  %s\n", dstrn, dstrx);
 717           printf (" ulp       : %s\n", ustrn);
 718           printf (" max.ulp   : %s\n", mustrn);
 719         }
 720     }
 721   update_stats (ok);
 722
 723  out:
 724   fpstack_test (test_name);
 725   feclearexcept (FE_ALL_EXCEPT);
 726   errno = 0;
 727 }
 728
 729
 730 void
 731 check_float (const char *test_name, FLOAT computed, FLOAT expected,
 732              int exceptions)
 733 {
 734   check_float_internal (test_name, computed, expected,
 735                         exceptions, &max_error, prev_max_error);
 736 }
 737
 738
 739 void
 740 check_complex (const char *test_name, CFLOAT computed,
 741                CFLOAT expected,
 742                int exception)
 743 {
 744   FLOAT part_comp, part_exp;
 745   char *str;
 746
 747   if (asprintf (&str, "Real part of: %s", test_name) == -1)
 748     abort ();
 749
 750   part_comp = __real__ computed;
 751   part_exp = __real__ expected;
 752
 753   check_float_internal (str, part_comp, part_exp,
 754                         exception, &real_max_error, prev_real_max_error);
 755   free (str);
 756
 757   if (asprintf (&str, "Imaginary part of: %s", test_name) == -1)
 758     abort ();
 759
 760   part_comp = __imag__ computed;
 761   part_exp = __imag__ expected;
 762
 763   /* Don't check again for exceptions or errno, just pass through the
 764      other relevant flags.  */
 765   check_float_internal (str, part_comp, part_exp,
 766                         exception & (IGNORE_ZERO_INF_SIGN
 767                                      | TEST_NAN_SIGN
 768                                      | IGNORE_RESULT),
 769                         &imag_max_error, prev_imag_max_error);
 770   free (str);
 771 }
 772
 773
 774 /* Check that computed and expected values are equal (int values).  */
 775 void
 776 check_int (const char *test_name, int computed, int expected,
 777            int exceptions)
 778 {
 779   int ok = 0;
 780   int errno_value = errno;
 781
 782   test_exceptions (test_name, exceptions);
 783   test_errno (test_name, errno_value, exceptions);
 784   if (exceptions & IGNORE_RESULT)
 785     goto out;
 786   noTests++;
 787   if (computed == expected)
 788     ok = 1;
 789
 790   if (print_screen (ok))
 791     {
 792       if (!ok)
 793         printf ("Failure: ");
 794       printf ("Test: %s\n", test_name);
 795       printf ("Result:\n");
 796       printf (" is:         %d\n", computed);
 797       printf (" should be:  %d\n", expected);
 798     }
 799
 800   update_stats (ok);
 801  out:
 802   fpstack_test (test_name);
 803   feclearexcept (FE_ALL_EXCEPT);
 804   errno = 0;
 805 }
 806
 807
 808 /* Check that computed and expected values are equal (long int values).  */
 809 void
 810 check_long (const char *test_name, long int computed, long int expected,
 811             int exceptions)
 812 {
 813   int ok = 0;
 814   int errno_value = errno;
 815
 816   test_exceptions (test_name, exceptions);
 817   test_errno (test_name, errno_value, exceptions);
 818   if (exceptions & IGNORE_RESULT)
 819     goto out;
 820   noTests++;
 821   if (computed == expected)
 822     ok = 1;
 823
 824   if (print_screen (ok))
 825     {
 826       if (!ok)
 827         printf ("Failure: ");
 828       printf ("Test: %s\n", test_name);
 829       printf ("Result:\n");
 830       printf (" is:         %ld\n", computed);
 831       printf (" should be:  %ld\n", expected);
 832     }
 833
 834   update_stats (ok);
 835  out:
 836   fpstack_test (test_name);
 837   feclearexcept (FE_ALL_EXCEPT);
 838   errno = 0;
 839 }
 840
 841
 842 /* Check that computed value is true/false.  */
 843 void
 844 check_bool (const char *test_name, int computed, int expected,
 845             int exceptions)
 846 {
 847   int ok = 0;
 848   int errno_value = errno;
 849
 850   test_exceptions (test_name, exceptions);
 851   test_errno (test_name, errno_value, exceptions);
 852   if (exceptions & IGNORE_RESULT)
 853     goto out;
 854   noTests++;
 855   if ((computed == 0) == (expected == 0))
 856     ok = 1;
 857
 858   if (print_screen (ok))
 859     {
 860       if (!ok)
 861         printf ("Failure: ");
 862       printf ("Test: %s\n", test_name);
 863       printf ("Result:\n");
 864       printf (" is:         %d\n", computed);
 865       printf (" should be:  %d\n", expected);
 866     }
 867
 868   update_stats (ok);
 869  out:
 870   fpstack_test (test_name);
 871   feclearexcept (FE_ALL_EXCEPT);
 872   errno = 0;
 873 }
 874
 875
 876 /* check that computed and expected values are equal (long int values) */
 877 void
 878 check_longlong (const char *test_name, long long int computed,
 879                 long long int expected,
 880                 int exceptions)
 881 {
 882   int ok = 0;
 883   int errno_value = errno;
 884
 885   test_exceptions (test_name, exceptions);
 886   test_errno (test_name, errno_value, exceptions);
 887   if (exceptions & IGNORE_RESULT)
 888     goto out;
 889   noTests++;
 890   if (computed == expected)
 891     ok = 1;
 892
 893   if (print_screen (ok))
 894     {
 895       if (!ok)
 896         printf ("Failure:");
 897       printf ("Test: %s\n", test_name);
 898       printf ("Result:\n");
 899       printf (" is:         %lld\n", computed);
 900       printf (" should be:  %lld\n", expected);
 901     }
 902
 903   update_stats (ok);
 904  out:
 905   fpstack_test (test_name);
 906   feclearexcept (FE_ALL_EXCEPT);
 907   errno = 0;
 908 }
 909
 910
 911 /* Check that computed and expected values are equal (intmax_t values).  */
 912 void
 913 check_intmax_t (const char *test_name, intmax_t computed,
 914                 intmax_t expected, int exceptions)
 915 {
 916   int ok = 0;
 917   int errno_value = errno;
 918
 919   test_exceptions (test_name, exceptions);
 920   test_errno (test_name, errno_value, exceptions);
 921   if (exceptions & IGNORE_RESULT)
 922     goto out;
 923   noTests++;
 924   if (computed == expected)
 925     ok = 1;
 926
 927   if (print_screen (ok))
 928     {
 929       if (!ok)
 930         printf ("Failure:");
 931       printf ("Test: %s\n", test_name);
 932       printf ("Result:\n");
 933       printf (" is:         %jd\n", computed);
 934       printf (" should be:  %jd\n", expected);
 935     }
 936
 937   update_stats (ok);
 938  out:
 939   fpstack_test (test_name);
 940   feclearexcept (FE_ALL_EXCEPT);
 941   errno = 0;
 942 }
 943
 944
 945 /* Check that computed and expected values are equal (uintmax_t values).  */
 946 void
 947 check_uintmax_t (const char *test_name, uintmax_t computed,
 948                  uintmax_t expected, int exceptions)
 949 {
 950   int ok = 0;
 951   int errno_value = errno;
 952
 953   test_exceptions (test_name, exceptions);
 954   test_errno (test_name, errno_value, exceptions);
 955   if (exceptions & IGNORE_RESULT)
 956     goto out;
 957   noTests++;
 958   if (computed == expected)
 959     ok = 1;
 960
 961   if (print_screen (ok))
 962     {
 963       if (!ok)
 964         printf ("Failure:");
 965       printf ("Test: %s\n", test_name);
 966       printf ("Result:\n");
 967       printf (" is:         %ju\n", computed);
 968       printf (" should be:  %ju\n", expected);
 969     }
 970
 971   update_stats (ok);
 972  out:
 973   fpstack_test (test_name);
 974   feclearexcept (FE_ALL_EXCEPT);
 975   errno = 0;
 976 }
 977
 978 /* Return whether a test with flags EXCEPTIONS should be run.  */
 979 int
 980 enable_test (int exceptions)
 981 {
 982   if (exceptions & XFAIL_TEST)
 983     return 0;
 984   if (flag_test_inline && (exceptions & NO_TEST_INLINE))
 985     return 0;
 986   if (flag_test_finite && (exceptions & NON_FINITE) != 0)
 987     return 0;
 988   if ((!SNAN_TESTS (FLOAT) || !snan_tests_arg)
 989       && (exceptions & TEST_SNAN) != 0)
 990     return 0;
 991   if (flag_test_mathvec && (exceptions & NO_TEST_MATHVEC) != 0)
 992     return 0;
 993
 994   return 1;
 995 }
 996
 997 static void
 998 initialize (void)
 999 {
1000   fpstack_test ("start *init*");
1001
1002   /* Clear all exceptions.  From now on we must not get random exceptions.  */
1003   feclearexcept (FE_ALL_EXCEPT);
1004   errno = 0;
1005
1006   /* Test to make sure we start correctly.  */
1007   fpstack_test ("end *init*");
1008 }
1009
1010 /* Definitions of arguments for argp functions.  */
1011 static const struct argp_option options[] =
1012 {
1013   { "verbose", 'v', "NUMBER", 0, "Level of verbosity (0..3)"},
1014   { "ulps-file", 'u', NULL, 0, "Output ulps to file ULPs"},
1015   { "no-max-error", 'f', NULL, 0,
1016     "Don't output maximal errors of functions"},
1017   { "no-points", 'p', NULL, 0,
1018     "Don't output results of functions invocations"},
1019   { "ignore-max-ulp", 'i', "yes/no", 0,
1020     "Ignore given maximal errors"},
1021   { "output-dir", 'o', "DIR", 0,
1022     "Directory where generated files will be placed"},
1023   { NULL, 0, NULL, 0, NULL }
1024 };
1025
1026 /* Prototype for option handler.  */
1027 static error_t parse_opt (int key, char *arg, struct argp_state *state);
1028
1029 /* Data structure to communicate with argp functions.  */
1030 static struct argp argp =
1031 {
1032   options, parse_opt, NULL, doc,
1033 };
1034
1035
1036 /* Handle program arguments.  */
1037 static error_t
1038 parse_opt (int key, char *arg, struct argp_state *state)
1039 {
1040   switch (key)
1041     {
1042     case 'f':
1043       output_max_error = 0;
1044       break;
1045     case 'i':
1046       if (strcmp (arg, "yes") == 0)
1047         ignore_max_ulp = 1;
1048       else if (strcmp (arg, "no") == 0)
1049         ignore_max_ulp = 0;
1050       break;
1051     case 'o':
1052       output_dir = (char *) malloc (strlen (arg) + 1);
1053       if (output_dir != NULL)
1054         strcpy (output_dir, arg);
1055       else
1056         return errno;
1057       break;
1058     case 'p':
1059       output_points = 0;
1060       break;
1061     case 'u':
1062       output_ulps = 1;
1063       break;
1064     case 'v':
1065       if (optarg)
1066         verbose = (unsigned int) strtoul (optarg, NULL, 0);
1067       else
1068         verbose = 3;
1069       break;
1070     default:
1071       return ARGP_ERR_UNKNOWN;
1072     }
1073   return 0;
1074 }
1075
1076 /* Verify that our ulp () implementation is behaving as expected
1077    or abort.  */
1078 static void
1079 check_ulp (void)
1080 {
1081    FLOAT ulps, ulpx, value;
1082    int i;
1083    /* Check ulp of zero is a subnormal value...  */
1084    ulps = ulp (0x0.0p0);
1085    if (fpclassify (ulps) != FP_SUBNORMAL)
1086      {
1087        fprintf (stderr, "ulp (0x0.0p0) is not FP_SUBNORMAL!\n");
1088        exit (EXIT_FAILURE);
1089      }
1090    /* Check that the ulp of one is a normal value... */
1091    ulps = ulp (LIT(1.0));
1092    if (fpclassify (ulps) != FP_NORMAL)
1093      {
1094        fprintf (stderr, "ulp (1.0L) is not FP_NORMAL\n");
1095        exit (EXIT_FAILURE);
1096      }
1097
1098    /* Compute the next subnormal value using nextafter to validate ulp.
1099       We allow +/- 1 ulp around the represented value.  */
1100    value = FUNC(nextafter) (0, 1);
1101    ulps = ULPDIFF (value, 0);
1102    ulpx = ulp (LIT(1.0));
1103    if (ulps < (LIT(1.0) - ulpx) || ulps > (LIT(1.0) + ulpx))
1104      {
1105        fprintf (stderr, "Value outside of 1 +/- 1ulp.\n");
1106        exit (EXIT_FAILURE);
1107      }
1108    /* Compute the nearest representable number from 10 towards 20.
1109       The result is 10 + 1ulp.  We use this to check the ulp function.
1110       We allow +/- 1 ulp around the represented value.  */
1111    value = FUNC(nextafter) (10, 20);
1112    ulps = ULPDIFF (value, 10);
1113    ulpx = ulp (LIT(1.0));
1114    if (ulps < (LIT(1.0) - ulpx) || ulps > (LIT(1.0) + ulpx))
1115      {
1116        fprintf (stderr, "Value outside of 1 +/- 1ulp.\n");
1117        exit (EXIT_FAILURE);
1118      }
1119    /* This gives one more ulp.  */
1120    value = FUNC(nextafter) (value, 20);
1121    ulps = ULPDIFF (value, 10);
1122    ulpx = ulp (LIT(2.0));
1123    if (ulps < (LIT(2.0) - ulpx) || ulps > (LIT(2.0) + ulpx))
1124      {
1125        fprintf (stderr, "Value outside of 2 +/- 1ulp.\n");
1126        exit (EXIT_FAILURE);
1127      }
1128    /* And now calculate 100 ulp.  */
1129    for (i = 2; i < 100; i++)
1130      value = FUNC(nextafter) (value, 20);
1131    ulps = ULPDIFF (value, 10);
1132    ulpx = ulp (LIT(100.0));
1133    if (ulps < (LIT(100.0) - ulpx) || ulps > (LIT(100.0) + ulpx))
1134      {
1135        fprintf (stderr, "Value outside of 100 +/- 1ulp.\n");
1136        exit (EXIT_FAILURE);
1137      }
1138 }
1139
1140 /* Do all initialization for a test run with arguments given by ARGC
1141    and ARGV.  */
1142 void
1143 libm_test_init (int argc, char **argv)
1144 {
1145   int remaining;
1146   char *ulps_file_path;
1147   size_t dir_len = 0;
1148
1149   verbose = 1;
1150   output_ulps = 0;
1151   output_max_error = 1;
1152   output_points = 1;
1153   output_dir = NULL;
1154   /* XXX set to 0 for releases.  */
1155   ignore_max_ulp = 0;
1156
1157   /* Parse and process arguments.  */
1158   argp_parse (&argp, argc, argv, 0, &remaining, NULL);
1159
1160   if (remaining != argc)
1161     {
1162       fprintf (stderr, "wrong number of arguments");
1163       argp_help (&argp, stdout, ARGP_HELP_SEE, program_invocation_short_name);
1164       exit (EXIT_FAILURE);
1165     }
1166
1167   if (output_ulps)
1168     {
1169       if (output_dir != NULL)
1170         dir_len = strlen (output_dir);
1171       ulps_file_path = (char *) malloc (dir_len + strlen (ulps_file_name) + 1);
1172       if (ulps_file_path == NULL)
1173         {
1174           perror ("can't allocate path for `ULPs' file: ");
1175           exit (1);
1176         }
1177       sprintf (ulps_file_path, "%s%s", output_dir == NULL ? "" : output_dir, ulps_file_name);
1178       ulps_file = fopen (ulps_file_path, "a");
1179       if (ulps_file == NULL)
1180         {
1181           perror ("can't open file `ULPs' for writing: ");
1182           exit (1);
1183         }
1184     }
1185
1186
1187   initialize ();
1188   fputs (test_msg, stdout);
1189
1190   check_ulp ();
1191 }
1192
1193 /* Process the test results, returning the exit status.  */
1194 int
1195 libm_test_finish (void)
1196 {
1197   if (output_ulps)
1198     fclose (ulps_file);
1199
1200   printf ("\nTest suite completed:\n");
1201   printf ("  %d test cases plus %d tests for exception flags and\n"
1202           "    %d tests for errno executed.\n",
1203           noTests, noExcTests, noErrnoTests);
1204   if (noErrors)
1205     {
1206       printf ("  %d errors occurred.\n", noErrors);
1207       return 1;
1208     }
1209   printf ("  All tests passed successfully.\n");
1210
1211   return 0;
1212 }