math/test-misc.c

   1 /* Miscellaneous tests which don't fit anywhere else.
   2    Copyright (C) 2000-2023 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    <https://www.gnu.org/licenses/>.  */
  18
  19 #include <fenv.h>
  20 #include <float.h>
  21 #include <ieee754.h>
  22 #include <math.h>
  23 #include <stdio.h>
  24 #include <string.h>
  25 #include <math-tests.h>
  26
  27
  28 static int
  29 do_test (void)
  30 {
  31   int result = 0;
  32
  33 #if LDBL_MANT_DIG > DBL_MANT_DIG
  34   {
  35     long double x = 0x100000001ll + (long double) 0.5;
  36     long double q;
  37     long double r;
  38
  39     r = modfl (x, &q);
  40     if (q != (long double) 0x100000001ll || r != 0.5)
  41       {
  42         printf ("modfl (%Lg, ...) failed\n", x);
  43         result = 1;
  44       }
  45   }
  46
  47   {
  48     long double x;
  49     long double m;
  50     long double r;
  51     int e;
  52     int i;
  53
  54 # if LDBL_MANT_DIG == 64
  55     m = 0xf.fffffffffffffffp-4L;
  56 # elif LDBL_MANT_DIG == 106
  57     /* This has to match the mantissa of LDBL_MAX which actually does have a
  58        missing bit in the middle.  */
  59     m = 0x1.fffffffffffff7ffffffffffff8p-1L;
  60 # elif LDBL_MANT_DIG == 113
  61     m = 0x1.ffffffffffffffffffffffffffffp-1L;
  62 # else
  63 #  error "Unsupported LDBL_MANT_DIG, please adjust"
  64 # endif
  65
  66     for (i = LDBL_MAX_EXP, x = LDBL_MAX; i >= LDBL_MIN_EXP; --i, x /= 2.0L)
  67       {
  68         printf ("2^%d: ", i);
  69
  70         r = frexpl (x, &e);
  71         if (r != m)
  72           {
  73             printf ("mantissa incorrect: %.20La\n", r);
  74             result = 1;
  75             continue;
  76           }
  77         if (e != i)
  78           {
  79             printf ("exponent wrong %d (%.20Lg)\n", e, x);
  80             result = 1;
  81             continue;
  82           }
  83         puts ("ok");
  84       }
  85
  86     for (i = LDBL_MIN_EXP, x = LDBL_MIN; i >= LDBL_MIN_EXP - LDBL_MANT_DIG + 1;
  87          --i, x /= 2.0L)
  88       {
  89         printf ("2^%d: ", i);
  90
  91         r = frexpl (x, &e);
  92         if (r != 0.5L)
  93           {
  94             printf ("mantissa incorrect: %.20La\n", r);
  95             result = 1;
  96             continue;
  97           }
  98         if (e != i)
  99           {
 100             printf ("exponent wrong %d (%.20Lg)\n", e, x);
 101             result = 1;
 102             continue;
 103           }
 104         puts ("ok");
 105       }
 106
 107   }
 108
 109 # if 0
 110   {
 111     int e;
 112     long double r = frexpl (LDBL_MIN * LDBL_EPSILON, &e);
 113
 114     if (r != 0.5)
 115       {
 116         printf ("frexpl (LDBL_MIN * LDBL_EPSILON, ...): mantissa wrong: %Lg\n",
 117                 r);
 118         result = 1;
 119       }
 120     else if (e != -16444)
 121       {
 122         printf ("frexpl (LDBL_MIN * LDBL_EPSILON, ...): exponent wrong: %d\n",
 123                 e);
 124         result = 1;
 125       }
 126   }
 127 # endif
 128 #endif
 129
 130   {
 131     double x = 0x100000001ll + (double) 0.5;
 132     double q;
 133     double r;
 134
 135     r = modf (x, &q);
 136     if (q != (double) 0x100000001ll || r != 0.5)
 137       {
 138         printf ("modf (%g, ...) failed\n", x);
 139         result = 1;
 140       }
 141   }
 142
 143   {
 144     union ieee754_float v1;
 145     union ieee754_float v2;
 146     float f;
 147
 148     v1.f = f = FLT_MIN;
 149     if (fpclassify (f) != FP_NORMAL)
 150       {
 151         printf ("fpclassify (FLT_MIN) failed: %d\n", fpclassify (f));
 152         result = 1;
 153       }
 154     f = nextafterf (f, FLT_MIN / 2.0f);
 155     if (fpclassify (f) != FP_SUBNORMAL)
 156       {
 157         printf ("fpclassify (FLT_MIN-epsilon) failed: %d\n", fpclassify (f));
 158         result = 1;
 159       }
 160     v2.f = f = nextafterf (f, FLT_MIN);
 161     if (fpclassify (f) != FP_NORMAL)
 162       {
 163         printf ("fpclassify (FLT_MIN-epsilon+epsilon) failed: %d\n",
 164                 fpclassify (f));
 165         result = 1;
 166       }
 167
 168     if (v1.ieee.mantissa != v2.ieee.mantissa)
 169       {
 170         printf ("FLT_MIN: mantissa differs: %8x vs %8x\n",
 171                 v1.ieee.mantissa, v2.ieee.mantissa);
 172         result = 1;
 173       }
 174     if (v1.ieee.exponent != v2.ieee.exponent)
 175       {
 176         printf ("FLT_MIN: exponent differs: %4x vs %4x\n",
 177                 v1.ieee.exponent, v2.ieee.exponent);
 178         result = 1;
 179       }
 180     if (v1.ieee.negative != v2.ieee.negative)
 181       {
 182         printf ("FLT_MIN: negative differs: %d vs %d\n",
 183                 v1.ieee.negative, v2.ieee.negative);
 184         result = 1;
 185       }
 186
 187     v1.f = f = -FLT_MIN;
 188     if (fpclassify (f) != FP_NORMAL)
 189       {
 190         printf ("fpclassify (-FLT_MIN) failed: %d\n", fpclassify (f));
 191         result = 1;
 192       }
 193     f = nextafterf (f, -FLT_MIN / 2.0f);
 194     if (fpclassify (f) != FP_SUBNORMAL)
 195       {
 196         printf ("fpclassify (-FLT_MIN-epsilon) failed: %d\n", fpclassify (f));
 197         result = 1;
 198       }
 199     v2.f = f = nextafterf (f, -FLT_MIN);
 200     if (fpclassify (f) != FP_NORMAL)
 201       {
 202         printf ("fpclassify (-FLT_MIN-epsilon+epsilon) failed: %d\n",
 203                 fpclassify (f));
 204         result = 1;
 205       }
 206
 207     if (v1.ieee.mantissa != v2.ieee.mantissa)
 208       {
 209         printf ("-FLT_MIN: mantissa differs: %8x vs %8x\n",
 210                 v1.ieee.mantissa, v2.ieee.mantissa);
 211         result = 1;
 212       }
 213     if (v1.ieee.exponent != v2.ieee.exponent)
 214       {
 215         printf ("-FLT_MIN: exponent differs: %4x vs %4x\n",
 216                 v1.ieee.exponent, v2.ieee.exponent);
 217         result = 1;
 218       }
 219     if (v1.ieee.negative != v2.ieee.negative)
 220       {
 221         printf ("-FLT_MIN: negative differs: %d vs %d\n",
 222                 v1.ieee.negative, v2.ieee.negative);
 223         result = 1;
 224       }
 225
 226     f = FLT_MAX;
 227     if (fpclassify (f) != FP_NORMAL)
 228       {
 229         printf ("fpclassify (FLT_MAX) failed: %d\n", fpclassify (f));
 230         result = 1;
 231       }
 232     f = nextafterf (f, INFINITY);
 233     if (fpclassify (f) != FP_INFINITE)
 234       {
 235         printf ("fpclassify (FLT_MAX+epsilon) failed: %d\n", fpclassify (f));
 236         result = 1;
 237       }
 238
 239     f = -FLT_MAX;
 240     if (fpclassify (f) != FP_NORMAL)
 241       {
 242         printf ("fpclassify (-FLT_MAX) failed: %d\n", fpclassify (f));
 243         result = 1;
 244       }
 245     f = nextafterf (f, -INFINITY);
 246     if (fpclassify (f) != FP_INFINITE)
 247       {
 248         printf ("fpclassify (-FLT_MAX-epsilon) failed: %d\n", fpclassify (f));
 249         result = 1;
 250       }
 251
 252     v1.f = f = 0.0625;
 253     f = nextafterf (f, 0.0);
 254     v2.f = f = nextafterf (f, 1.0);
 255
 256     if (v1.ieee.mantissa != v2.ieee.mantissa)
 257       {
 258         printf ("0.0625f down: mantissa differs: %8x vs %8x\n",
 259                 v1.ieee.mantissa, v2.ieee.mantissa);
 260         result = 1;
 261       }
 262     if (v1.ieee.exponent != v2.ieee.exponent)
 263       {
 264         printf ("0.0625f down: exponent differs: %4x vs %4x\n",
 265                 v1.ieee.exponent, v2.ieee.exponent);
 266         result = 1;
 267       }
 268     if (v1.ieee.negative != v2.ieee.negative)
 269       {
 270         printf ("0.0625f down: negative differs: %d vs %d\n",
 271                 v1.ieee.negative, v2.ieee.negative);
 272         result = 1;
 273       }
 274
 275     v1.f = f = 0.0625;
 276     f = nextafterf (f, 1.0);
 277     v2.f = f = nextafterf (f, 0.0);
 278
 279     if (v1.ieee.mantissa != v2.ieee.mantissa)
 280       {
 281         printf ("0.0625f up: mantissa differs: %8x vs %8x\n",
 282                 v1.ieee.mantissa, v2.ieee.mantissa);
 283         result = 1;
 284       }
 285     if (v1.ieee.exponent != v2.ieee.exponent)
 286       {
 287         printf ("0.0625f up: exponent differs: %4x vs %4x\n",
 288                 v1.ieee.exponent, v2.ieee.exponent);
 289         result = 1;
 290       }
 291     if (v1.ieee.negative != v2.ieee.negative)
 292       {
 293         printf ("0.0625f up: negative differs: %d vs %d\n",
 294                 v1.ieee.negative, v2.ieee.negative);
 295         result = 1;
 296       }
 297
 298     v1.f = f = -0.0625;
 299     f = nextafterf (f, 0.0);
 300     v2.f = f = nextafterf (f, -1.0);
 301
 302     if (v1.ieee.mantissa != v2.ieee.mantissa)
 303       {
 304         printf ("-0.0625f up: mantissa differs: %8x vs %8x\n",
 305                 v1.ieee.mantissa, v2.ieee.mantissa);
 306         result = 1;
 307       }
 308     if (v1.ieee.exponent != v2.ieee.exponent)
 309       {
 310         printf ("-0.0625f up: exponent differs: %4x vs %4x\n",
 311                 v1.ieee.exponent, v2.ieee.exponent);
 312         result = 1;
 313       }
 314     if (v1.ieee.negative != v2.ieee.negative)
 315       {
 316         printf ("-0.0625f up: negative differs: %d vs %d\n",
 317                 v1.ieee.negative, v2.ieee.negative);
 318         result = 1;
 319       }
 320
 321     v1.f = f = -0.0625;
 322     f = nextafterf (f, -1.0);
 323     v2.f = f = nextafterf (f, 0.0);
 324
 325     if (v1.ieee.mantissa != v2.ieee.mantissa)
 326       {
 327         printf ("-0.0625f down: mantissa differs: %8x vs %8x\n",
 328                 v1.ieee.mantissa, v2.ieee.mantissa);
 329         result = 1;
 330       }
 331     if (v1.ieee.exponent != v2.ieee.exponent)
 332       {
 333         printf ("-0.0625f down: exponent differs: %4x vs %4x\n",
 334                 v1.ieee.exponent, v2.ieee.exponent);
 335         result = 1;
 336       }
 337     if (v1.ieee.negative != v2.ieee.negative)
 338       {
 339         printf ("-0.0625f down: negative differs: %d vs %d\n",
 340                 v1.ieee.negative, v2.ieee.negative);
 341         result = 1;
 342       }
 343
 344     v1.f = f = 0.0f;
 345     f = nextafterf (f, 1.0);
 346     v2.f = nextafterf (f, -1.0);
 347
 348     if (v1.ieee.mantissa != v2.ieee.mantissa)
 349       {
 350         printf ("0.0f up: mantissa differs: %8x vs %8x\n",
 351                 v1.ieee.mantissa, v2.ieee.mantissa);
 352         result = 1;
 353       }
 354     if (v1.ieee.exponent != v2.ieee.exponent)
 355       {
 356         printf ("0.0f up: exponent differs: %4x vs %4x\n",
 357                 v1.ieee.exponent, v2.ieee.exponent);
 358         result = 1;
 359       }
 360     if (0 != v2.ieee.negative)
 361       {
 362         printf ("0.0f up: negative differs: 0 vs %d\n",
 363                 v2.ieee.negative);
 364         result = 1;
 365       }
 366
 367     v1.f = f = 0.0f;
 368     f = nextafterf (f, -1.0);
 369     v2.f = nextafterf (f, 1.0);
 370
 371     if (v1.ieee.mantissa != v2.ieee.mantissa)
 372       {
 373         printf ("0.0f down: mantissa differs: %8x vs %8x\n",
 374                 v1.ieee.mantissa, v2.ieee.mantissa);
 375         result = 1;
 376       }
 377     if (v1.ieee.exponent != v2.ieee.exponent)
 378       {
 379         printf ("0.0f down: exponent differs: %4x vs %4x\n",
 380                 v1.ieee.exponent, v2.ieee.exponent);
 381         result = 1;
 382       }
 383     if (1 != v2.ieee.negative)
 384       {
 385         printf ("0.0f down: negative differs: 1 vs %d\n",
 386                 v2.ieee.negative);
 387         result = 1;
 388       }
 389
 390     if (nextafterf (0.0f, INFINITY) != nextafterf (0.0f, 1.0f)
 391         || nextafterf (-0.0f, INFINITY) != nextafterf (-0.0f, 1.0f)
 392         || nextafterf (0.0f, -INFINITY) != nextafterf (0.0f, -1.0f)
 393         || nextafterf (-0.0f, -INFINITY) != nextafterf (-0.0f, -1.0f))
 394       {
 395         printf ("nextafterf (+-0, +-Inf) != nextafterf (+-0, +-1)\n");
 396         result = 1;
 397       }
 398
 399     if (nexttowardf (0.0f, INFINITY) != nexttowardf (0.0f, 1.0f)
 400         || nexttowardf (-0.0f, INFINITY) != nexttowardf (-0.0f, 1.0f)
 401         || nexttowardf (0.0f, -INFINITY) != nexttowardf (0.0f, -1.0f)
 402         || nexttowardf (-0.0f, -INFINITY) != nexttowardf (-0.0f, -1.0f))
 403       {
 404         printf ("nexttowardf (+-0, +-Inf) != nexttowardf (+-0, +-1)\n");
 405         result = 1;
 406       }
 407   }
 408
 409   {
 410     union ieee754_double v1;
 411     union ieee754_double v2;
 412     double d;
 413
 414     v1.d = d = DBL_MIN;
 415     if (fpclassify (d) != FP_NORMAL)
 416       {
 417         printf ("fpclassify (DBL_MIN) failed: %d\n", fpclassify (d));
 418         result = 1;
 419       }
 420     d = nextafter (d, DBL_MIN / 2.0);
 421     if (fpclassify (d) != FP_SUBNORMAL)
 422       {
 423         printf ("fpclassify (DBL_MIN-epsilon) failed: %d\n", fpclassify (d));
 424         result = 1;
 425       }
 426     v2.d = d = nextafter (d, DBL_MIN);
 427     if (fpclassify (d) != FP_NORMAL)
 428       {
 429         printf ("fpclassify (DBL_MIN-epsilon+epsilon) failed: %d\n",
 430                 fpclassify (d));
 431         result = 1;
 432       }
 433
 434     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 435       {
 436         printf ("DBL_MIN: mantissa0 differs: %8x vs %8x\n",
 437                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 438         result = 1;
 439       }
 440     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 441       {
 442         printf ("DBL_MIN: mantissa1 differs: %8x vs %8x\n",
 443                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 444         result = 1;
 445       }
 446     if (v1.ieee.exponent != v2.ieee.exponent)
 447       {
 448         printf ("DBL_MIN: exponent differs: %4x vs %4x\n",
 449                 v1.ieee.exponent, v2.ieee.exponent);
 450         result = 1;
 451       }
 452     if (v1.ieee.negative != v2.ieee.negative)
 453       {
 454         printf ("DBL_MIN: negative differs: %d vs %d\n",
 455                 v1.ieee.negative, v2.ieee.negative);
 456         result = 1;
 457       }
 458
 459     v1.d = d = -DBL_MIN;
 460     if (fpclassify (d) != FP_NORMAL)
 461       {
 462         printf ("fpclassify (-DBL_MIN) failed: %d\n", fpclassify (d));
 463         result = 1;
 464       }
 465     d = nextafter (d, -DBL_MIN / 2.0);
 466     if (fpclassify (d) != FP_SUBNORMAL)
 467       {
 468         printf ("fpclassify (-DBL_MIN-epsilon) failed: %d\n", fpclassify (d));
 469         result = 1;
 470       }
 471     v2.d = d = nextafter (d, -DBL_MIN);
 472     if (fpclassify (d) != FP_NORMAL)
 473       {
 474         printf ("fpclassify (-DBL_MIN-epsilon+epsilon) failed: %d\n",
 475                 fpclassify (d));
 476         result = 1;
 477       }
 478
 479     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 480       {
 481         printf ("-DBL_MIN: mantissa0 differs: %8x vs %8x\n",
 482                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 483         result = 1;
 484       }
 485     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 486       {
 487         printf ("-DBL_MIN: mantissa1 differs: %8x vs %8x\n",
 488                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 489         result = 1;
 490       }
 491     if (v1.ieee.exponent != v2.ieee.exponent)
 492       {
 493         printf ("-DBL_MIN: exponent differs: %4x vs %4x\n",
 494                 v1.ieee.exponent, v2.ieee.exponent);
 495         result = 1;
 496       }
 497     if (v1.ieee.negative != v2.ieee.negative)
 498       {
 499         printf ("-DBL_MIN: negative differs: %d vs %d\n",
 500                 v1.ieee.negative, v2.ieee.negative);
 501         result = 1;
 502       }
 503
 504     d = DBL_MAX;
 505     if (fpclassify (d) != FP_NORMAL)
 506       {
 507         printf ("fpclassify (DBL_MAX) failed: %d\n", fpclassify (d));
 508         result = 1;
 509       }
 510     d = nextafter (d, INFINITY);
 511     if (fpclassify (d) != FP_INFINITE)
 512       {
 513         printf ("fpclassify (DBL_MAX+epsilon) failed: %d\n", fpclassify (d));
 514         result = 1;
 515       }
 516
 517     d = -DBL_MAX;
 518     if (fpclassify (d) != FP_NORMAL)
 519       {
 520         printf ("fpclassify (-DBL_MAX) failed: %d\n", fpclassify (d));
 521         result = 1;
 522       }
 523     d = nextafter (d, -INFINITY);
 524     if (fpclassify (d) != FP_INFINITE)
 525       {
 526         printf ("fpclassify (-DBL_MAX-epsilon) failed: %d\n", fpclassify (d));
 527         result = 1;
 528       }
 529
 530     v1.d = d = 0.0625;
 531     d = nextafter (d, 0.0);
 532     v2.d = d = nextafter (d, 1.0);
 533
 534     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 535       {
 536         printf ("0.0625 down: mantissa0 differs: %8x vs %8x\n",
 537                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 538         result = 1;
 539       }
 540     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 541       {
 542         printf ("0.0625 down: mantissa1 differs: %8x vs %8x\n",
 543                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 544         result = 1;
 545       }
 546     if (v1.ieee.exponent != v2.ieee.exponent)
 547       {
 548         printf ("0.0625 down: exponent differs: %4x vs %4x\n",
 549                 v1.ieee.exponent, v2.ieee.exponent);
 550         result = 1;
 551       }
 552     if (v1.ieee.negative != v2.ieee.negative)
 553       {
 554         printf ("0.0625 down: negative differs: %d vs %d\n",
 555                 v1.ieee.negative, v2.ieee.negative);
 556         result = 1;
 557       }
 558
 559     v1.d = d = 0.0625;
 560     d = nextafter (d, 1.0);
 561     v2.d = d = nextafter (d, 0.0);
 562
 563     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 564       {
 565         printf ("0.0625 up: mantissa0 differs: %8x vs %8x\n",
 566                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 567         result = 1;
 568       }
 569     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 570       {
 571         printf ("0.0625 up: mantissa1 differs: %8x vs %8x\n",
 572                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 573         result = 1;
 574       }
 575     if (v1.ieee.exponent != v2.ieee.exponent)
 576       {
 577         printf ("0.0625 up: exponent differs: %4x vs %4x\n",
 578                 v1.ieee.exponent, v2.ieee.exponent);
 579         result = 1;
 580       }
 581     if (v1.ieee.negative != v2.ieee.negative)
 582       {
 583         printf ("0.0625 up: negative differs: %d vs %d\n",
 584                 v1.ieee.negative, v2.ieee.negative);
 585         result = 1;
 586       }
 587
 588     v1.d = d = -0.0625;
 589     d = nextafter (d, 0.0);
 590     v2.d = d = nextafter (d, -1.0);
 591
 592     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 593       {
 594         printf ("-0.0625 up: mantissa0 differs: %8x vs %8x\n",
 595                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 596         result = 1;
 597       }
 598     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 599       {
 600         printf ("-0.0625 up: mantissa1 differs: %8x vs %8x\n",
 601                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 602         result = 1;
 603       }
 604     if (v1.ieee.exponent != v2.ieee.exponent)
 605       {
 606         printf ("-0.0625 up: exponent differs: %4x vs %4x\n",
 607                 v1.ieee.exponent, v2.ieee.exponent);
 608         result = 1;
 609       }
 610     if (v1.ieee.negative != v2.ieee.negative)
 611       {
 612         printf ("-0.0625 up: negative differs: %d vs %d\n",
 613                 v1.ieee.negative, v2.ieee.negative);
 614         result = 1;
 615       }
 616
 617     v1.d = d = -0.0625;
 618     d = nextafter (d, -1.0);
 619     v2.d = d = nextafter (d, 0.0);
 620
 621     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 622       {
 623         printf ("-0.0625 down: mantissa0 differs: %8x vs %8x\n",
 624                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 625         result = 1;
 626       }
 627     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 628       {
 629         printf ("-0.0625 down: mantissa1 differs: %8x vs %8x\n",
 630                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 631         result = 1;
 632       }
 633     if (v1.ieee.exponent != v2.ieee.exponent)
 634       {
 635         printf ("-0.0625 down: exponent differs: %4x vs %4x\n",
 636                 v1.ieee.exponent, v2.ieee.exponent);
 637         result = 1;
 638       }
 639     if (v1.ieee.negative != v2.ieee.negative)
 640       {
 641         printf ("-0.0625 down: negative differs: %d vs %d\n",
 642                 v1.ieee.negative, v2.ieee.negative);
 643         result = 1;
 644       }
 645
 646     v1.d = d = 0.0;
 647     d = nextafter (d, 1.0);
 648     v2.d = nextafter (d, -1.0);
 649
 650     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 651       {
 652         printf ("0.0 up: mantissa0 differs: %8x vs %8x\n",
 653                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 654         result = 1;
 655       }
 656     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 657       {
 658         printf ("0.0 up: mantissa1 differs: %8x vs %8x\n",
 659                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 660         result = 1;
 661       }
 662     if (v1.ieee.exponent != v2.ieee.exponent)
 663       {
 664         printf ("0.0 up: exponent differs: %4x vs %4x\n",
 665                 v1.ieee.exponent, v2.ieee.exponent);
 666         result = 1;
 667       }
 668     if (0 != v2.ieee.negative)
 669       {
 670         printf ("0.0 up: negative differs: 0 vs %d\n",
 671                 v2.ieee.negative);
 672         result = 1;
 673       }
 674
 675     v1.d = d = 0.0;
 676     d = nextafter (d, -1.0);
 677     v2.d = nextafter (d, 1.0);
 678
 679     if (v1.ieee.mantissa0 != v2.ieee.mantissa0)
 680       {
 681         printf ("0.0 down: mantissa0 differs: %8x vs %8x\n",
 682                 v1.ieee.mantissa0, v2.ieee.mantissa0);
 683         result = 1;
 684       }
 685     if (v1.ieee.mantissa1 != v2.ieee.mantissa1)
 686       {
 687         printf ("0.0 down: mantissa1 differs: %8x vs %8x\n",
 688                 v1.ieee.mantissa1, v2.ieee.mantissa1);
 689         result = 1;
 690       }
 691     if (v1.ieee.exponent != v2.ieee.exponent)
 692       {
 693         printf ("0.0 down: exponent differs: %4x vs %4x\n",
 694                 v1.ieee.exponent, v2.ieee.exponent);
 695         result = 1;
 696       }
 697     if (1 != v2.ieee.negative)
 698       {
 699         printf ("0.0 down: negative differs: 1 vs %d\n",
 700                 v2.ieee.negative);
 701         result = 1;
 702       }
 703
 704     if (nextafter (0.0, INFINITY) != nextafter (0.0, 1.0)
 705         || nextafter (-0.0, INFINITY) != nextafter (-0.0, 1.0)
 706         || nextafter (0.0, -INFINITY) != nextafter (0.0, -1.0)
 707         || nextafter (-0.0, -INFINITY) != nextafter (-0.0, -1.0))
 708       {
 709         printf ("nextafter (+-0, +-Inf) != nextafter (+-0, +-1)\n");
 710         result = 1;
 711       }
 712
 713     if (nexttoward (0.0, INFINITY) != nexttoward (0.0, 1.0)
 714         || nexttoward (-0.0, INFINITY) != nexttoward (-0.0, 1.0)
 715         || nexttoward (0.0, -INFINITY) != nexttoward (0.0, -1.0)
 716         || nexttoward (-0.0, -INFINITY) != nexttoward (-0.0, -1.0))
 717       {
 718         printf ("nexttoward (+-0, +-Inf) != nexttoward (+-0, +-1)\n");
 719         result = 1;
 720       }
 721   }
 722
 723 #if LDBL_MANT_DIG > DBL_MANT_DIG
 724   {
 725     long double v1, v2;
 726
 727     v1 = LDBL_MIN;
 728     if (fpclassify (v1) != FP_NORMAL)
 729       {
 730         printf ("fpclassify (LDBL_MIN) failed: %d (%La)\n",
 731                 fpclassify (v1), v1);
 732         result = 1;
 733       }
 734     v2 = nextafterl (v1, LDBL_MIN / 2.0);
 735     if (fpclassify (v2) != FP_SUBNORMAL)
 736       {
 737         printf ("fpclassify (LDBL_MIN-epsilon) failed: %d (%La)\n",
 738                 fpclassify (v2), v2);
 739         result = 1;
 740       }
 741     v2 = nextafterl (v2, LDBL_MIN);
 742     if (fpclassify (v2) != FP_NORMAL)
 743       {
 744         printf ("fpclassify (LDBL_MIN-epsilon+epsilon) failed: %d (%La)\n",
 745                 fpclassify (v2), v2);
 746         result = 1;
 747       }
 748
 749     if (v1 != v2)
 750       {
 751         printf ("LDBL_MIN-epsilon+epsilon != LDBL_MIN: %La vs %La\n", v2, v1);
 752         result = 1;
 753       }
 754
 755     v1 = -LDBL_MIN;
 756     if (fpclassify (v1) != FP_NORMAL)
 757       {
 758         printf ("fpclassify (-LDBL_MIN) failed: %d (%La)\n",
 759                 fpclassify (v1), v1);
 760         result = 1;
 761       }
 762     v2 = nextafterl (v1, -LDBL_MIN / 2.0);
 763     if (fpclassify (v2) != FP_SUBNORMAL)
 764       {
 765         printf ("fpclassify (-LDBL_MIN-epsilon) failed: %d (%La)\n",
 766                 fpclassify (v2), v2);
 767         result = 1;
 768       }
 769     v2 = nextafterl (v2, -LDBL_MIN);
 770     if (fpclassify (v2) != FP_NORMAL)
 771       {
 772         printf ("fpclassify (-LDBL_MIN-epsilon+epsilon) failed: %d (%La)\n",
 773                 fpclassify (v2), v2);
 774         result = 1;
 775       }
 776
 777     if (v1 != v2)
 778       {
 779         printf ("-LDBL_MIN-epsilon+epsilon != -LDBL_MIN: %La vs %La\n", v2, v1);
 780         result = 1;
 781       }
 782
 783     v1 = LDBL_MAX;
 784     if (fpclassify (v1) != FP_NORMAL)
 785       {
 786         printf ("fpclassify (LDBL_MAX) failed: %d (%La)\n",
 787                 fpclassify (v1), v1);
 788         result = 1;
 789       }
 790     v2 = nextafterl (v1, INFINITY);
 791     if (fpclassify (v2) != FP_INFINITE)
 792       {
 793         printf ("fpclassify (LDBL_MAX+epsilon) failed: %d (%La)\n",
 794                 fpclassify (v2), v2);
 795         result = 1;
 796       }
 797
 798     v1 = -LDBL_MAX;
 799     if (fpclassify (v1) != FP_NORMAL)
 800       {
 801         printf ("fpclassify (-LDBL_MAX) failed: %d (%La)\n",
 802                 fpclassify (v1), v1);
 803         result = 1;
 804       }
 805     v2 = nextafterl (v1, -INFINITY);
 806     if (fpclassify (v2) != FP_INFINITE)
 807       {
 808         printf ("fpclassify (-LDBL_MAX-epsilon) failed: %d (%La)\n",
 809                 fpclassify (v2), v2);
 810         result = 1;
 811       }
 812
 813     v1 = 0.0625;
 814     v2 = nextafterl (v1, 0.0);
 815     v2 = nextafterl (v2, 1.0);
 816
 817     if (v1 != v2)
 818       {
 819         printf ("0.0625L-epsilon+epsilon != 0.0625L: %La vs %La\n", v2, v1);
 820         result = 1;
 821       }
 822
 823     v1 = 0.0625;
 824     v2 = nextafterl (v1, 1.0);
 825     v2 = nextafterl (v2, 0.0);
 826
 827     if (v1 != v2)
 828       {
 829         printf ("0.0625L+epsilon-epsilon != 0.0625L: %La vs %La\n", v2, v1);
 830         result = 1;
 831       }
 832
 833     v1 = -0.0625;
 834     v2 = nextafterl (v1, 0.0);
 835     v2 = nextafterl (v2, -1.0);
 836
 837     if (v1 != v2)
 838       {
 839         printf ("-0.0625L+epsilon-epsilon != -0.0625L: %La vs %La\n", v2, v1);
 840         result = 1;
 841       }
 842
 843     v1 = -0.0625;
 844     v2 = nextafterl (v1, -1.0);
 845     v2 = nextafterl (v2, 0.0);
 846
 847     if (v1 != v2)
 848       {
 849         printf ("-0.0625L-epsilon+epsilon != -0.0625L: %La vs %La\n", v2, v1);
 850         result = 1;
 851       }
 852
 853     v1 = 0.0;
 854     v2 = nextafterl (v1, 1.0);
 855     v2 = nextafterl (v2, -1.0);
 856
 857     if (v1 != v2)
 858       {
 859         printf ("0.0+epsilon-epsilon != 0.0L: %La vs %La\n", v2, v1);
 860         result = 1;
 861       }
 862     if (signbit (v2))
 863       {
 864         printf ("0.0+epsilon-epsilon is negative\n");
 865         result = 1;
 866       }
 867
 868     v1 = 0.0;
 869     v2 = nextafterl (v1, -1.0);
 870     v2 = nextafterl (v2, 1.0);
 871
 872     if (v1 != v2)
 873       {
 874         printf ("0.0-epsilon+epsilon != 0.0L: %La vs %La\n", v2, v1);
 875         result = 1;
 876       }
 877     if (!signbit (v2))
 878       {
 879         printf ("0.0-epsilon+epsilon is positive\n");
 880         result = 1;
 881       }
 882
 883     if (nextafterl (0.0, INFINITY) != nextafterl (0.0, 1.0)
 884         || nextafterl (-0.0, INFINITY) != nextafterl (-0.0, 1.0)
 885         || nextafterl (0.0, -INFINITY) != nextafterl (0.0, -1.0)
 886         || nextafterl (-0.0, -INFINITY) != nextafterl (-0.0, -1.0))
 887       {
 888         printf ("nextafterl (+-0, +-Inf) != nextafterl (+-0, +-1)\n");
 889         result = 1;
 890       }
 891
 892     if (nexttowardl (0.0L, INFINITY) != nexttowardl (0.0L, 1.0L)
 893         || nexttowardl (-0.0L, INFINITY) != nexttowardl (-0.0L, 1.0L)
 894         || nexttowardl (0.0L, -INFINITY) != nexttowardl (0.0L, -1.0L)
 895         || nexttowardl (-0.0L, -INFINITY) != nexttowardl (-0.0L, -1.0L))
 896       {
 897         printf ("nexttowardl (+-0, +-Inf) != nexttowardl (+-0, +-1)\n");
 898         result = 1;
 899       }
 900   }
 901 #endif
 902
 903   if (! isnormal (FLT_MIN))
 904     {
 905       puts ("isnormal (FLT_MIN) failed");
 906       result = 1;
 907     }
 908   if (! isnormal (DBL_MIN))
 909     {
 910       puts ("isnormal (DBL_MIN) failed");
 911       result = 1;
 912     }
 913 #if LDBL_MANT_DIG > DBL_MANT_DIG
 914   if (! isnormal (LDBL_MIN))
 915     {
 916       puts ("isnormal (LDBL_MIN) failed");
 917       result = 1;
 918     }
 919 #endif
 920
 921 #if defined (__i386__) || defined (__x86_64__)
 922   /* This is a test for the strange long doubles in x86 FPUs.  */
 923   {
 924     union
 925     {
 926       char b[10];
 927       long double d;
 928     } u =
 929       { .b = { 0, 0, 0, 0, 0, 0, 0, 0x80, 0, 0 } };
 930
 931     if (fpclassify (u.d) != FP_NORMAL)
 932       {
 933         printf ("fpclassify (0x00008000000000000000) failed: %d (%Lg)\n",
 934                 fpclassify (u.d), u.d);
 935         result = 1;
 936       }
 937   }
 938
 939   /* Special qNaNs in x86 long double.  Test for scalbl.  */
 940   {
 941     union
 942     {
 943       char b[10];
 944       long double d;
 945     } u =
 946       { .b = { 0, 1, 0, 0, 0, 0, 0, 0xc0, 0xff, 0x7f } };
 947     long double r;
 948
 949     r = scalbl (u.d, 0.0);
 950     if (!isnan (r))
 951       {
 952         puts ("scalbl (qNaN, 0) does not return NaN");
 953         result = 1;
 954       }
 955     else if (memcmp (&r, &u.d, sizeof (double)) != 0)
 956       {
 957         puts ("scalbl (qNaN, 0) does not return the same NaN");
 958         result = 1;
 959       }
 960   }
 961 #endif
 962
 963 #if LDBL_MANT_DIG > DBL_MANT_DIG
 964   {
 965     long double r;
 966
 967     feclearexcept (FE_ALL_EXCEPT);
 968     r = scalbl (LDBL_MIN, 2147483647);
 969     if (! isinf (r))
 970       {
 971         puts ("scalbl (LDBL_MIN, 2147483647) does not return Inf");
 972         result = 1;
 973       }
 974     else if (signbit (r) != 0)
 975       {
 976         puts ("scalbl (LDBL_MIN, 2147483647) returns -Inf");
 977         result = 1;
 978       }
 979 # ifdef FE_UNDERFLOW
 980     else if (fetestexcept (FE_UNDERFLOW))
 981       {
 982         puts ("scalbl (LDBL_MIN, 2147483647) raises underflow exception");
 983         result = 1;
 984       }
 985 # endif
 986
 987     feclearexcept (FE_ALL_EXCEPT);
 988     r = scalbl (LDBL_MAX, -2147483647);
 989     if (r != 0.0)
 990       {
 991         puts ("scalbl (LDBL_MAX, -2147483647) does not return 0");
 992         result = 1;
 993       }
 994     else if (signbit (r) != 0)
 995       {
 996         puts ("scalbl (LDBL_MAX, -2147483647) returns -Inf");
 997         result = 1;
 998       }
 999 # ifdef FE_OVERFLOW
1000     else if (fetestexcept (FE_OVERFLOW))
1001       {
1002         puts ("scalbl (LDBL_MAX, -2147483647) raises overflow exception");
1003         result = 1;
1004       }
1005 # endif
1006   }
1007 #endif
1008
1009   /* The tests here are very similar to tests earlier in this file,
1010      the important difference is just that there are no intervening
1011      union variables that cause some GCC versions to hide possible
1012      bugs in nextafter* implementation.  */
1013   if (nextafterf (nextafterf (FLT_MIN, FLT_MIN / 2.0), FLT_MIN) != FLT_MIN)
1014     {
1015       puts ("nextafterf FLT_MIN test failed");
1016       result = 1;
1017     }
1018   if (nextafterf (nextafterf (-FLT_MIN, -FLT_MIN / 2.0), -FLT_MIN)
1019       != -FLT_MIN)
1020     {
1021       puts ("nextafterf -FLT_MIN test failed");
1022       result = 1;
1023     }
1024   if (nextafter (nextafter (DBL_MIN, DBL_MIN / 2.0), DBL_MIN) != DBL_MIN)
1025     {
1026       puts ("nextafter DBL_MIN test failed");
1027       result = 1;
1028     }
1029   if (nextafter (nextafter (-DBL_MIN, -DBL_MIN / 2.0), -DBL_MIN) != -DBL_MIN)
1030     {
1031       puts ("nextafter -DBL_MIN test failed");
1032       result = 1;
1033     }
1034 #if LDBL_MANT_DIG > DBL_MANT_DIG
1035   if (nextafterl (nextafterl (LDBL_MIN, LDBL_MIN / 2.0), LDBL_MIN)
1036       != LDBL_MIN)
1037     {
1038       puts ("nextafterl LDBL_MIN test failed");
1039       result = 1;
1040     }
1041   if (nextafterl (nextafterl (-LDBL_MIN, -LDBL_MIN / 2.0), -LDBL_MIN)
1042       != -LDBL_MIN)
1043     {
1044       puts ("nextafterl -LDBL_MIN test failed");
1045       result = 1;
1046     }
1047 #endif
1048
1049   volatile float f1 = FLT_MAX;
1050   volatile float f2 = FLT_MAX / 2;
1051   (void) &f1;
1052   (void) &f2;
1053   feclearexcept (FE_ALL_EXCEPT);
1054   f2 += f1;
1055 #if defined(FE_OVERFLOW) && defined(FE_INEXACT)
1056   int fe = fetestexcept (FE_ALL_EXCEPT);
1057   if (EXCEPTION_TESTS (float) && fe != (FE_OVERFLOW | FE_INEXACT))
1058     {
1059       printf ("float overflow test failed: %x\n", fe);
1060       result = 1;
1061     }
1062 #endif
1063
1064   volatile double d1 = DBL_MAX;
1065   volatile double d2 = DBL_MAX / 2;
1066   (void) &d1;
1067   (void) &d2;
1068   feclearexcept (FE_ALL_EXCEPT);
1069   d2 += d1;
1070 #if defined(FE_OVERFLOW) && defined(FE_INEXACT)
1071   fe = fetestexcept (FE_ALL_EXCEPT);
1072   if (EXCEPTION_TESTS (double) && fe != (FE_OVERFLOW | FE_INEXACT))
1073     {
1074       printf ("double overflow test failed: %x\n", fe);
1075       result = 1;
1076     }
1077 #endif
1078
1079 #if LDBL_MANT_DIG > DBL_MANT_DIG
1080   volatile long double ld1 = LDBL_MAX;
1081   volatile long double ld2 = LDBL_MAX / 2;
1082   (void) &ld1;
1083   (void) &ld2;
1084   feclearexcept (FE_ALL_EXCEPT);
1085   ld2 += ld1;
1086 # if defined(FE_OVERFLOW) && defined(FE_INEXACT)
1087   fe = fetestexcept (FE_ALL_EXCEPT);
1088   if (EXCEPTION_TESTS (long double) && fe != (FE_OVERFLOW | FE_INEXACT))
1089     {
1090       printf ("long double overflow test failed: %x\n", fe);
1091       result = 1;
1092     }
1093 # endif
1094
1095 # if LDBL_MANT_DIG == 113
1096   volatile long double ld3 = 0x1.0000000000010000000100000001p+1;
1097   volatile long double ld4 = 0x1.0000000000000000000000000001p+1;
1098   (void) &ld3;
1099   (void) &ld4;
1100   ld3 -= ld4;
1101   if (ld3 != 0x1.0p-47)
1102     {
1103       printf ("long double subtraction test failed %.28La\n", ld3);
1104       result = 1;
1105     }
1106 # endif
1107
1108 /* Skip testing IBM long double format, for 2 reasons:
1109    1) it only supports FE_TONEAREST
1110    2) nextafter (0.0, 1.0) == nextafterl (0.0L, 1.0L), so
1111       nextafter (0.0, 1.0) / 16.0L will be 0.0L.  */
1112 # if LDBL_MANT_DIG >= DBL_MANT_DIG + 4 && LDBL_MANT_DIG != 106
1113   int oldmode = fegetround ();
1114   int j;
1115   for (j = 0; j < 4; j++)
1116     {
1117       int mode;
1118       int i;
1119       int k = 0;
1120       const char *mstr;
1121       switch (j)
1122         {
1123 #ifdef FE_TONEAREST
1124         case 0:
1125           mode = FE_TONEAREST;
1126           mstr = "nearest";
1127           k = 8;
1128           break;
1129 #endif
1130 #ifdef FE_DOWNWARD
1131         case 1:
1132           mode = FE_DOWNWARD;
1133           mstr = "-inf";
1134           break;
1135 #endif
1136 #ifdef FE_UPWARD
1137         case 2:
1138           mode = FE_UPWARD;
1139           mstr = "+inf";
1140           k = 15;
1141           break;
1142 #endif
1143 #ifdef FE_TOWARDZERO
1144         case 3:
1145           mode = FE_TOWARDZERO;
1146           mstr = "0";
1147           break;
1148 #endif
1149         default:
1150           continue;
1151         }
1152
1153       volatile long double ld5 = nextafter (0.0, 1.0) / 16.0L;
1154       volatile double d5;
1155       (void) &ld5;
1156       for (i = 0; i <= 32; i++)
1157         {
1158           if (fesetround (mode))
1159             {
1160               printf ("failed to set rounding mode to %s\n", mstr);
1161               if (ROUNDING_TESTS (long double, mode)
1162                   && ROUNDING_TESTS (double, mode))
1163                 result = 1;
1164               else
1165                 puts ("ignoring this failure");
1166               break;
1167             }
1168           d5 = ld5 * i;
1169           (void) &d5;
1170           fesetround (oldmode);
1171           if (d5 != ((j == 0 && i == 8) ? 0 : (i + k) / 16)
1172                     * nextafter (0.0, 1.0))
1173             {
1174               printf ("%La incorrectly rounded to %s as %a\n",
1175                       ld5 * i, mstr, d5);
1176               if (ROUNDING_TESTS (long double, mode)
1177                   && ROUNDING_TESTS (double, mode))
1178                 result = 1;
1179               else
1180                 puts ("ignoring this failure");
1181             }
1182         }
1183     }
1184
1185 #  ifdef FE_UPWARD
1186   volatile long double ld7 = nextafterl (0.0L, 1.0L);
1187   volatile double d7;
1188   (void) &ld7;
1189   fesetround (FE_UPWARD);
1190   d7 = ld7;
1191   (void) &d7;
1192   fesetround (oldmode);
1193
1194   if (d7 != nextafter (0.0, 1.0))
1195     {
1196       printf ("%La incorrectly rounded upward to %a\n", ld7, d7);
1197       if (ROUNDING_TESTS (long double, FE_UPWARD)
1198           && ROUNDING_TESTS (double, FE_UPWARD))
1199         result = 1;
1200       else
1201         puts ("ignoring this failure");
1202     }
1203 #  endif
1204 # endif
1205 #endif
1206
1207   return result;
1208 }
1209
1210 #define TEST_FUNCTION do_test ()
1211 #include "../test-skeleton.c"