1 /*
   2  * CDDL HEADER START
   3  *
   4  * The contents of this file are subject to the terms of the
   5  * Common Development and Distribution License (the "License").
   6  * You may not use this file except in compliance with the License.
   7  *
   8  * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
   9  * or http://www.opensolaris.org/os/licensing.
  10  * See the License for the specific language governing permissions
  11  * and limitations under the License.
  12  *
  13  * When distributing Covered Code, include this CDDL HEADER in each
  14  * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
  15  * If applicable, add the following below this CDDL HEADER, with the
  16  * fields enclosed by brackets "[]" replaced with your own identifying
  17  * information: Portions Copyright [yyyy] [name of copyright owner]
  18  *
  19  * CDDL HEADER END
  20  */
  21 /*
  22  * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
  23  */
  24 /*
  25  * Copyright 2005 Sun Microsystems, Inc.  All rights reserved.
  26  * Use is subject to license terms.
  27  */
  28 
  29 #include "libm.h"
  30 #include "xpg6.h"       /* __xpg6 */
  31 #include <stdio.h>
  32 #include <float.h>                /* DBL_MAX, DBL_MIN */
  33 #include <unistd.h>               /* write */
  34 #if defined(__x86)
  35 #include <ieeefp.h>
  36 #undef  fp_class
  37 #define fp_class        fpclass
  38 #define fp_quiet        FP_QNAN
  39 #endif
  40 #include <errno.h>
  41 #undef fflush
  42 #include <sys/isa_defs.h>
  43 
  44 /* INDENT OFF */
  45 /*
  46  * Report libm exception error according to System V Interface Definition
  47  * (SVID).
  48  * Error mapping:
  49  *      1 -- acos(|x|>1)
  50  *      2 -- asin(|x|>1)
  51  *      3 -- atan2(+-0,+-0)
  52  *      4 -- hypot overflow
  53  *      5 -- cosh overflow
  54  *      6 -- exp overflow
  55  *      7 -- exp underflow
  56  *      8 -- y0(0)
  57  *      9 -- y0(-ve)
  58  *      10-- y1(0)
  59  *      11-- y1(-ve)
  60  *      12-- yn(0)
  61  *      13-- yn(-ve)
  62  *      14-- lgamma(finite) overflow
  63  *      15-- lgamma(-integer)
  64  *      16-- log(0)
  65  *      17-- log(x<0)
  66  *      18-- log10(0)
  67  *      19-- log10(x<0)
  68  *      20-- pow(0.0,0.0)
  69  *      21-- pow(x,y) overflow
  70  *      22-- pow(x,y) underflow
  71  *      23-- pow(0,negative)
  72  *      24-- pow(neg,non-integral)
  73  *      25-- sinh(finite) overflow
  74  *      26-- sqrt(negative)
  75  *      27-- fmod(x,0)
  76  *      28-- remainder(x,0)
  77  *      29-- acosh(x<1)
  78  *      30-- atanh(|x|>1)
  79  *      31-- atanh(|x|=1)
  80  *      32-- scalb overflow
  81  *      33-- scalb underflow
  82  *      34-- j0(|x|>X_TLOSS)
  83  *      35-- y0(x>X_TLOSS)
  84  *      36-- j1(|x|>X_TLOSS)
  85  *      37-- y1(x>X_TLOSS)
  86  *      38-- jn(|x|>X_TLOSS, n)
  87  *      39-- yn(x>X_TLOSS, n)
  88  *      40-- gamma(finite) overflow
  89  *      41-- gamma(-integer)
  90  *      42-- pow(NaN,0.0) return NaN for SVID/XOPEN
  91  *      43-- log1p(-1)
  92  *      44-- log1p(x<-1)
  93  *      45-- logb(0)
  94  *      46-- nextafter overflow
  95  *      47-- scalb(x,inf)
  96  */
  97 /* INDENT ON */
  98 
  99 static double setexception(int, double);
 100 
 101 static const union {
 102         unsigned        x[2];
 103         double          d;
 104 } C[] = {
 105 #ifdef _LITTLE_ENDIAN
 106         { 0xffffffff, 0x7fffffff },
 107         { 0x54442d18, 0x400921fb },
 108 #else
 109         { 0x7fffffff, 0xffffffff },
 110         { 0x400921fb, 0x54442d18 },
 111 #endif
 112 };
 113 
 114 #define NaN     C[0].d
 115 #define PI_RZ   C[1].d
 116 
 117 #define __HI(x) ((unsigned *)&x)[HIWORD]
 118 #define __LO(x) ((unsigned *)&x)[LOWORD]
 119 #undef  Inf
 120 #define Inf     HUGE_VAL
 121 
 122 double
 123 _SVID_libm_err(double x, double y, int type) {
 124         struct exception        exc;
 125         double                  t, w, ieee_retval = 0;
 126         enum version            lib_version = _lib_version;
 127         int                     iy;
 128 
 129         /* force libm_ieee behavior in SUSv3 mode */
 130         if ((__xpg6 & _C99SUSv3_math_errexcept) != 0)
 131                 lib_version = libm_ieee;
 132         if (lib_version == c_issue_4) {
 133                 (void) fflush(stdout);
 134         }
 135         exc.arg1 = x;
 136         exc.arg2 = y;
 137         switch (type) {
 138         case 1:
 139                 /* acos(|x|>1) */
 140                 exc.type = DOMAIN;
 141                 exc.name = "acos";
 142                 ieee_retval = setexception(3, 1.0);
 143                 exc.retval = 0.0;
 144                 if (lib_version == strict_ansi) {
 145                         errno = EDOM;
 146                 } else if (!matherr(&exc)) {
 147                         if (lib_version == c_issue_4) {
 148                                 (void) write(2, "acos: DOMAIN error\n", 19);
 149                         }
 150                         errno = EDOM;
 151                 }
 152                 break;
 153         case 2:
 154                 /* asin(|x|>1) */
 155                 exc.type = DOMAIN;
 156                 exc.name = "asin";
 157                 exc.retval = 0.0;
 158                 ieee_retval = setexception(3, 1.0);
 159                 if (lib_version == strict_ansi) {
 160                         errno = EDOM;
 161                 } else if (!matherr(&exc)) {
 162                         if (lib_version == c_issue_4) {
 163                                 (void) write(2, "asin: DOMAIN error\n", 19);
 164                         }
 165                         errno = EDOM;
 166                 }
 167                 break;
 168         case 3:
 169                 /* atan2(+-0,+-0) */
 170                 exc.arg1 = y;
 171                 exc.arg2 = x;
 172                 exc.type = DOMAIN;
 173                 exc.name = "atan2";
 174                 ieee_retval = copysign(1.0, x) == 1.0 ? y :
 175                         copysign(PI_RZ + DBL_MIN, y);
 176                 exc.retval = 0.0;
 177                 if (lib_version == strict_ansi) {
 178                         errno = EDOM;
 179                 } else if (!matherr(&exc)) {
 180                         if (lib_version == c_issue_4) {
 181                                 (void) write(2, "atan2: DOMAIN error\n", 20);
 182                         }
 183                         errno = EDOM;
 184                 }
 185                 break;
 186         case 4:
 187                 /* hypot(finite,finite) overflow */
 188                 exc.type = OVERFLOW;
 189                 exc.name = "hypot";
 190                 ieee_retval = Inf;
 191                 if (lib_version == c_issue_4)
 192                         exc.retval = HUGE;
 193                 else
 194                         exc.retval = HUGE_VAL;
 195                 if (lib_version == strict_ansi)
 196                         errno = ERANGE;
 197                 else if (!matherr(&exc))
 198                         errno = ERANGE;
 199                 break;
 200         case 5:
 201                 /* cosh(finite) overflow */
 202                 exc.type = OVERFLOW;
 203                 exc.name = "cosh";
 204                 ieee_retval = setexception(2, 1.0);
 205                 if (lib_version == c_issue_4)
 206                         exc.retval = HUGE;
 207                 else
 208                         exc.retval = HUGE_VAL;
 209                 if (lib_version == strict_ansi)
 210                         errno = ERANGE;
 211                 else if (!matherr(&exc))
 212                         errno = ERANGE;
 213                 break;
 214         case 6:
 215                 /* exp(finite) overflow */
 216                 exc.type = OVERFLOW;
 217                 exc.name = "exp";
 218                 ieee_retval = setexception(2, 1.0);
 219                 if (lib_version == c_issue_4)
 220                         exc.retval = HUGE;
 221                 else
 222                         exc.retval = HUGE_VAL;
 223                 if (lib_version == strict_ansi)
 224                         errno = ERANGE;
 225                 else if (!matherr(&exc))
 226                         errno = ERANGE;
 227                 break;
 228         case 7:
 229                 /* exp(finite) underflow */
 230                 exc.type = UNDERFLOW;
 231                 exc.name = "exp";
 232                 ieee_retval = setexception(1, 1.0);
 233                 exc.retval = 0.0;
 234                 if (lib_version == strict_ansi)
 235                         errno = ERANGE;
 236                 else if (!matherr(&exc))
 237                         errno = ERANGE;
 238                 break;
 239         case 8:
 240                 /* y0(0) = -inf */
 241                 exc.type = DOMAIN;      /* should be SING for IEEE */
 242                 exc.name = "y0";
 243                 ieee_retval = setexception(0, -1.0);
 244                 if (lib_version == c_issue_4)
 245                         exc.retval = -HUGE;
 246                 else
 247                         exc.retval = -HUGE_VAL;
 248                 if (lib_version == strict_ansi) {
 249                         errno = EDOM;
 250                 } else if (!matherr(&exc)) {
 251                         if (lib_version == c_issue_4) {
 252                                 (void) write(2, "y0: DOMAIN error\n", 17);
 253                         }
 254                         errno = EDOM;
 255                 }
 256                 break;
 257         case 9:
 258                 /* y0(x<0) = NaN */
 259                 exc.type = DOMAIN;
 260                 exc.name = "y0";
 261                 ieee_retval = setexception(3, 1.0);
 262                 if (lib_version == c_issue_4)
 263                         exc.retval = -HUGE;
 264                 else
 265                         exc.retval = -HUGE_VAL;
 266                 if (lib_version == strict_ansi) {
 267                         errno = EDOM;
 268                 } else if (!matherr(&exc)) {
 269                         if (lib_version == c_issue_4) {
 270                                 (void) write(2, "y0: DOMAIN error\n", 17);
 271                         }
 272                         errno = EDOM;
 273                 }
 274                 break;
 275         case 10:
 276                 /* y1(0) = -inf */
 277                 exc.type = DOMAIN;      /* should be SING for IEEE */
 278                 exc.name = "y1";
 279                 ieee_retval = setexception(0, -1.0);
 280                 if (lib_version == c_issue_4)
 281                         exc.retval = -HUGE;
 282                 else
 283                         exc.retval = -HUGE_VAL;
 284                 if (lib_version == strict_ansi) {
 285                         errno = EDOM;
 286                 } else if (!matherr(&exc)) {
 287                         if (lib_version == c_issue_4) {
 288                                 (void) write(2, "y1: DOMAIN error\n", 17);
 289                         }
 290                         errno = EDOM;
 291                 }
 292                 break;
 293         case 11:
 294                 /* y1(x<0) = NaN */
 295                 exc.type = DOMAIN;
 296                 exc.name = "y1";
 297                 ieee_retval = setexception(3, 1.0);
 298                 if (lib_version == c_issue_4)
 299                         exc.retval = -HUGE;
 300                 else
 301                         exc.retval = -HUGE_VAL;
 302                 if (lib_version == strict_ansi) {
 303                         errno = EDOM;
 304                 } else if (!matherr(&exc)) {
 305                         if (lib_version == c_issue_4) {
 306                                 (void) write(2, "y1: DOMAIN error\n", 17);
 307                         }
 308                         errno = EDOM;
 309                 }
 310                 break;
 311         case 12:
 312                 /* yn(n,0) = -inf */
 313                 exc.type = DOMAIN;      /* should be SING for IEEE */
 314                 exc.name = "yn";
 315                 ieee_retval = setexception(0, -1.0);
 316                 if (lib_version == c_issue_4)
 317                         exc.retval = -HUGE;
 318                 else
 319                         exc.retval = -HUGE_VAL;
 320                 if (lib_version == strict_ansi) {
 321                         errno = EDOM;
 322                 } else if (!matherr(&exc)) {
 323                         if (lib_version == c_issue_4) {
 324                                 (void) write(2, "yn: DOMAIN error\n", 17);
 325                         }
 326                         errno = EDOM;
 327                 }
 328                 break;
 329         case 13:
 330                 /* yn(x<0) = NaN */
 331                 exc.type = DOMAIN;
 332                 exc.name = "yn";
 333                 ieee_retval = setexception(3, 1.0);
 334                 if (lib_version == c_issue_4)
 335                         exc.retval = -HUGE;
 336                 else
 337                         exc.retval = -HUGE_VAL;
 338                 if (lib_version == strict_ansi) {
 339                         errno = EDOM;
 340                 } else if (!matherr(&exc)) {
 341                         if (lib_version == c_issue_4) {
 342                                 (void) write(2, "yn: DOMAIN error\n", 17);
 343                         }
 344                         errno = EDOM;
 345                 }
 346                 break;
 347         case 14:
 348                 /* lgamma(finite) overflow */
 349                 exc.type = OVERFLOW;
 350                 exc.name = "lgamma";
 351                 ieee_retval = setexception(2, 1.0);
 352                 if (lib_version == c_issue_4)
 353                         exc.retval = HUGE;
 354                 else
 355                         exc.retval = HUGE_VAL;
 356                 if (lib_version == strict_ansi)
 357                         errno = ERANGE;
 358                 else if (!matherr(&exc))
 359                         errno = ERANGE;
 360                 break;
 361         case 15:
 362                 /* lgamma(-integer) or lgamma(0) */
 363                 exc.type = SING;
 364                 exc.name = "lgamma";
 365                 ieee_retval = setexception(0, 1.0);
 366                 if (lib_version == c_issue_4)
 367                         exc.retval = HUGE;
 368                 else
 369                         exc.retval = HUGE_VAL;
 370                 if (lib_version == strict_ansi) {
 371                         errno = EDOM;
 372                 } else if (!matherr(&exc)) {
 373                         if (lib_version == c_issue_4) {
 374                                 (void) write(2, "lgamma: SING error\n", 19);
 375                         }
 376                         errno = EDOM;
 377                 }
 378                 break;
 379         case 16:
 380                 /* log(0) */
 381                 exc.type = SING;
 382                 exc.name = "log";
 383                 ieee_retval = setexception(0, -1.0);
 384                 if (lib_version == c_issue_4)
 385                         exc.retval = -HUGE;
 386                 else
 387                         exc.retval = -HUGE_VAL;
 388                 if (lib_version == strict_ansi) {
 389                         errno = ERANGE;
 390                 } else if (!matherr(&exc)) {
 391                         if (lib_version == c_issue_4) {
 392                                 (void) write(2, "log: SING error\n", 16);
 393                                 errno = EDOM;
 394                         } else {
 395                                 errno = ERANGE;
 396                         }
 397                 }
 398                 break;
 399         case 17:
 400                 /* log(x<0) */
 401                 exc.type = DOMAIN;
 402                 exc.name = "log";
 403                 ieee_retval = setexception(3, 1.0);
 404                 if (lib_version == c_issue_4)
 405                         exc.retval = -HUGE;
 406                 else
 407                         exc.retval = -HUGE_VAL;
 408                 if (lib_version == strict_ansi) {
 409                         errno = EDOM;
 410                 } else if (!matherr(&exc)) {
 411                         if (lib_version == c_issue_4) {
 412                                 (void) write(2, "log: DOMAIN error\n", 18);
 413                         }
 414                         errno = EDOM;
 415                 }
 416                 break;
 417         case 18:
 418                 /* log10(0) */
 419                 exc.type = SING;
 420                 exc.name = "log10";
 421                 ieee_retval = setexception(0, -1.0);
 422                 if (lib_version == c_issue_4)
 423                         exc.retval = -HUGE;
 424                 else
 425                         exc.retval = -HUGE_VAL;
 426                 if (lib_version == strict_ansi) {
 427                         errno = ERANGE;
 428                 } else if (!matherr(&exc)) {
 429                         if (lib_version == c_issue_4) {
 430                                 (void) write(2, "log10: SING error\n", 18);
 431                                 errno = EDOM;
 432                         } else {
 433                                 errno = ERANGE;
 434                         }
 435                 }
 436                 break;
 437         case 19:
 438                 /* log10(x<0) */
 439                 exc.type = DOMAIN;
 440                 exc.name = "log10";
 441                 ieee_retval = setexception(3, 1.0);
 442                 if (lib_version == c_issue_4)
 443                         exc.retval = -HUGE;
 444                 else
 445                         exc.retval = -HUGE_VAL;
 446                 if (lib_version == strict_ansi) {
 447                         errno = EDOM;
 448                 } else if (!matherr(&exc)) {
 449                         if (lib_version == c_issue_4) {
 450                                 (void) write(2, "log10: DOMAIN error\n", 20);
 451                         }
 452                         errno = EDOM;
 453                 }
 454                 break;
 455         case 20:
 456                 /* pow(0.0,0.0) */
 457                 /* error only if lib_version == c_issue_4 */
 458                 exc.type = DOMAIN;
 459                 exc.name = "pow";
 460                 exc.retval = 0.0;
 461                 ieee_retval = 1.0;
 462                 if (lib_version != c_issue_4) {
 463                         exc.retval = 1.0;
 464                 } else if (!matherr(&exc)) {
 465                         (void) write(2, "pow(0,0): DOMAIN error\n", 23);
 466                         errno = EDOM;
 467                 }
 468                 break;
 469         case 21:
 470                 /* pow(x,y) overflow */
 471                 exc.type = OVERFLOW;
 472                 exc.name = "pow";
 473                 exc.retval = (lib_version == c_issue_4)? HUGE : HUGE_VAL;
 474                 if (signbit(x)) {
 475                         t = rint(y);
 476                         if (t == y) {
 477                                 w = rint(0.5 * y);
 478                                 if (t != w + w) {       /* y is odd */
 479                                         exc.retval = -exc.retval;
 480                                 }
 481                         }
 482                 }
 483                 ieee_retval = setexception(2, exc.retval);
 484                 if (lib_version == strict_ansi)
 485                         errno = ERANGE;
 486                 else if (!matherr(&exc))
 487                         errno = ERANGE;
 488                 break;
 489         case 22:
 490                 /* pow(x,y) underflow */
 491                 exc.type = UNDERFLOW;
 492                 exc.name = "pow";
 493                 exc.retval = 0.0;
 494                 if (signbit(x)) {
 495                         t = rint(y);
 496                         if (t == y) {
 497                                 w = rint(0.5 * y);
 498                                 if (t != w + w) /* y is odd */
 499                                         exc.retval = -exc.retval;
 500                         }
 501                 }
 502                 ieee_retval = setexception(1, exc.retval);
 503                 if (lib_version == strict_ansi)
 504                         errno = ERANGE;
 505                 else if (!matherr(&exc))
 506                         errno = ERANGE;
 507                 break;
 508         case 23:
 509                 /* (+-0)**neg */
 510                 exc.type = DOMAIN;
 511                 exc.name = "pow";
 512                 ieee_retval = setexception(0, 1.0);
 513                 {
 514                         int ahy, k, j, yisint, ly, hx;
 515                         /* INDENT OFF */
 516                         /*
 517                          * determine if y is an odd int when x = -0
 518                          * yisint = 0       ... y is not an integer
 519                          * yisint = 1       ... y is an odd int
 520                          * yisint = 2       ... y is an even int
 521                          */
 522                         /* INDENT ON */
 523                         hx  = __HI(x);
 524                         ahy = __HI(y)&0x7fffffff;
 525                         ly  = __LO(y);
 526 
 527                         yisint = 0;
 528                         if (ahy >= 0x43400000) {
 529                                 yisint = 2;     /* even integer y */
 530                         } else if (ahy >= 0x3ff00000) {
 531                                 k = (ahy >> 20) - 0x3ff;  /* exponent */
 532                                 if (k > 20) {
 533                                         j = ly >> (52 - k);
 534                                         if ((j << (52 - k)) == ly)
 535                                                 yisint = 2 - (j & 1);
 536                                 } else if (ly == 0) {
 537                                         j = ahy >> (20 - k);
 538                                         if ((j << (20 - k)) == ahy)
 539                                                 yisint = 2 - (j & 1);
 540                                 }
 541                         }
 542                         if (hx < 0 && yisint == 1)
 543                                 ieee_retval = -ieee_retval;
 544                 }
 545                 if (lib_version == c_issue_4)
 546                         exc.retval = 0.0;
 547                 else
 548                         exc.retval = -HUGE_VAL;
 549                 if (lib_version == strict_ansi) {
 550                         errno = EDOM;
 551                 } else if (!matherr(&exc)) {
 552                         if (lib_version == c_issue_4) {
 553                                 (void) write(2, "pow(0,neg): DOMAIN error\n",
 554                                     25);
 555                         }
 556                         errno = EDOM;
 557                 }
 558                 break;
 559         case 24:
 560                 /* neg**non-integral */
 561                 exc.type = DOMAIN;
 562                 exc.name = "pow";
 563                 ieee_retval = setexception(3, 1.0);
 564                 if (lib_version == c_issue_4)
 565                         exc.retval = 0.0;
 566                 else
 567                         exc.retval = ieee_retval;       /* X/Open allow NaN */
 568                 if (lib_version == strict_ansi) {
 569                         errno = EDOM;
 570                 } else if (!matherr(&exc)) {
 571                         if (lib_version == c_issue_4) {
 572                                 (void) write(2,
 573                                     "neg**non-integral: DOMAIN error\n", 32);
 574                         }
 575                         errno = EDOM;
 576                 }
 577                 break;
 578         case 25:
 579                 /* sinh(finite) overflow */
 580                 exc.type = OVERFLOW;
 581                 exc.name = "sinh";
 582                 ieee_retval = copysign(Inf, x);
 583                 if (lib_version == c_issue_4)
 584                         exc.retval = x > 0.0 ? HUGE : -HUGE;
 585                 else
 586                         exc.retval = x > 0.0 ? HUGE_VAL : -HUGE_VAL;
 587                 if (lib_version == strict_ansi)
 588                         errno = ERANGE;
 589                 else if (!matherr(&exc))
 590                         errno = ERANGE;
 591                 break;
 592         case 26:
 593                 /* sqrt(x<0) */
 594                 exc.type = DOMAIN;
 595                 exc.name = "sqrt";
 596                 ieee_retval = setexception(3, 1.0);
 597                 if (lib_version == c_issue_4)
 598                         exc.retval = 0.0;
 599                 else
 600                         exc.retval = ieee_retval;       /* quiet NaN */
 601                 if (lib_version == strict_ansi) {
 602                         errno = EDOM;
 603                 } else if (!matherr(&exc)) {
 604                         if (lib_version == c_issue_4) {
 605                                 (void) write(2, "sqrt: DOMAIN error\n", 19);
 606                         }
 607                         errno = EDOM;
 608                 }
 609                 break;
 610         case 27:
 611                 /* fmod(x,0) */
 612                 exc.type = DOMAIN;
 613                 exc.name = "fmod";
 614                 if (fp_class(x) == fp_quiet)
 615                         ieee_retval = NaN;
 616                 else
 617                         ieee_retval = setexception(3, 1.0);
 618                 if (lib_version == c_issue_4)
 619                         exc.retval = x;
 620                 else
 621                         exc.retval = ieee_retval;
 622                 if (lib_version == strict_ansi) {
 623                         errno = EDOM;
 624                 } else if (!matherr(&exc)) {
 625                         if (lib_version == c_issue_4) {
 626                                 (void) write(2, "fmod:  DOMAIN error\n", 20);
 627                         }
 628                         errno = EDOM;
 629                 }
 630                 break;
 631         case 28:
 632                 /* remainder(x,0) */
 633                 exc.type = DOMAIN;
 634                 exc.name = "remainder";
 635                 if (fp_class(x) == fp_quiet)
 636                         ieee_retval = NaN;
 637                 else
 638                         ieee_retval = setexception(3, 1.0);
 639                 exc.retval = NaN;
 640                 if (lib_version == strict_ansi) {
 641                         errno = EDOM;
 642                 } else if (!matherr(&exc)) {
 643                         if (lib_version == c_issue_4) {
 644                                 (void) write(2, "remainder: DOMAIN error\n",
 645                                     24);
 646                         }
 647                         errno = EDOM;
 648                 }
 649                 break;
 650         case 29:
 651                 /* acosh(x<1) */
 652                 exc.type = DOMAIN;
 653                 exc.name = "acosh";
 654                 ieee_retval = setexception(3, 1.0);
 655                 exc.retval = NaN;
 656                 if (lib_version == strict_ansi) {
 657                         errno = EDOM;
 658                 } else if (!matherr(&exc)) {
 659                         if (lib_version == c_issue_4) {
 660                                 (void) write(2, "acosh: DOMAIN error\n", 20);
 661                         }
 662                         errno = EDOM;
 663                 }
 664                 break;
 665         case 30:
 666                 /* atanh(|x|>1) */
 667                 exc.type = DOMAIN;
 668                 exc.name = "atanh";
 669                 ieee_retval = setexception(3, 1.0);
 670                 exc.retval = NaN;
 671                 if (lib_version == strict_ansi) {
 672                         errno = EDOM;
 673                 } else if (!matherr(&exc)) {
 674                         if (lib_version == c_issue_4) {
 675                                 (void) write(2, "atanh: DOMAIN error\n", 20);
 676                         }
 677                         errno = EDOM;
 678                 }
 679                 break;
 680         case 31:
 681                 /* atanh(|x|=1) */
 682                 exc.type = SING;
 683                 exc.name = "atanh";
 684                 ieee_retval = setexception(0, x);
 685                 exc.retval = ieee_retval;
 686                 if (lib_version == strict_ansi) {
 687                         errno = ERANGE;
 688                 } else if (!matherr(&exc)) {
 689                         if (lib_version == c_issue_4) {
 690                                 (void) write(2, "atanh: SING error\n", 18);
 691                                 errno = EDOM;
 692                         } else {
 693                                 errno = ERANGE;
 694                         }
 695                 }
 696                 break;
 697         case 32:
 698                 /* scalb overflow; SVID also returns +-HUGE_VAL */
 699                 exc.type = OVERFLOW;
 700                 exc.name = "scalb";
 701                 ieee_retval = setexception(2, x);
 702                 exc.retval = x > 0.0 ? HUGE_VAL : -HUGE_VAL;
 703                 if (lib_version == strict_ansi)
 704                         errno = ERANGE;
 705                 else if (!matherr(&exc))
 706                         errno = ERANGE;
 707                 break;
 708         case 33:
 709                 /* scalb underflow */
 710                 exc.type = UNDERFLOW;
 711                 exc.name = "scalb";
 712                 ieee_retval = setexception(1, x);
 713                 exc.retval = ieee_retval;       /* +-0.0 */
 714                 if (lib_version == strict_ansi)
 715                         errno = ERANGE;
 716                 else if (!matherr(&exc))
 717                         errno = ERANGE;
 718                 break;
 719         case 34:
 720                 /* j0(|x|>X_TLOSS) */
 721                 exc.type = TLOSS;
 722                 exc.name = "j0";
 723                 exc.retval = 0.0;
 724                 ieee_retval = y;
 725                 if (lib_version == strict_ansi) {
 726                         errno = ERANGE;
 727                 } else if (!matherr(&exc)) {
 728                         if (lib_version == c_issue_4) {
 729                                 (void) write(2, exc.name, 2);
 730                                 (void) write(2, ": TLOSS error\n", 14);
 731                         }
 732                         errno = ERANGE;
 733                 }
 734                 break;
 735         case 35:
 736                 /* y0(x>X_TLOSS) */
 737                 exc.type = TLOSS;
 738                 exc.name = "y0";
 739                 exc.retval = 0.0;
 740                 ieee_retval = y;
 741                 if (lib_version == strict_ansi) {
 742                         errno = ERANGE;
 743                 } else if (!matherr(&exc)) {
 744                         if (lib_version == c_issue_4) {
 745                                 (void) write(2, exc.name, 2);
 746                                 (void) write(2, ": TLOSS error\n", 14);
 747                         }
 748                         errno = ERANGE;
 749                 }
 750                 break;
 751         case 36:
 752                 /* j1(|x|>X_TLOSS) */
 753                 exc.type = TLOSS;
 754                 exc.name = "j1";
 755                 exc.retval = 0.0;
 756                 ieee_retval = y;
 757                 if (lib_version == strict_ansi) {
 758                         errno = ERANGE;
 759                 } else if (!matherr(&exc)) {
 760                         if (lib_version == c_issue_4) {
 761                                 (void) write(2, exc.name, 2);
 762                                 (void) write(2, ": TLOSS error\n", 14);
 763                         }
 764                         errno = ERANGE;
 765                 }
 766                 break;
 767         case 37:
 768                 /* y1(x>X_TLOSS) */
 769                 exc.type = TLOSS;
 770                 exc.name = "y1";
 771                 exc.retval = 0.0;
 772                 ieee_retval = y;
 773                 if (lib_version == strict_ansi) {
 774                         errno = ERANGE;
 775                 } else if (!matherr(&exc)) {
 776                         if (lib_version == c_issue_4) {
 777                                 (void) write(2, exc.name, 2);
 778                                 (void) write(2, ": TLOSS error\n", 14);
 779                         }
 780                         errno = ERANGE;
 781                 }
 782                 break;
 783         case 38:
 784                 /* jn(|x|>X_TLOSS) */
 785                 /* incorrect ieee value: ieee should never be here */
 786                 exc.type = TLOSS;
 787                 exc.name = "jn";
 788                 exc.retval = 0.0;
 789                 ieee_retval = 0.0;      /* shall not be used */
 790                 if (lib_version == strict_ansi) {
 791                         errno = ERANGE;
 792                 } else if (!matherr(&exc)) {
 793                         if (lib_version == c_issue_4) {
 794                                 (void) write(2, exc.name, 2);
 795                                 (void) write(2, ": TLOSS error\n", 14);
 796                         }
 797                         errno = ERANGE;
 798                 }
 799                 break;
 800         case 39:
 801                 /* yn(x>X_TLOSS) */
 802                 /* incorrect ieee value: ieee should never be here */
 803                 exc.type = TLOSS;
 804                 exc.name = "yn";
 805                 exc.retval = 0.0;
 806                 ieee_retval = 0.0;      /* shall not be used */
 807                 if (lib_version == strict_ansi) {
 808                         errno = ERANGE;
 809                 } else if (!matherr(&exc)) {
 810                         if (lib_version == c_issue_4) {
 811                                 (void) write(2, exc.name, 2);
 812                                 (void) write(2, ": TLOSS error\n", 14);
 813                         }
 814                         errno = ERANGE;
 815                 }
 816                 break;
 817         case 40:
 818                 /* gamma(finite) overflow */
 819                 exc.type = OVERFLOW;
 820                 exc.name = "gamma";
 821                 ieee_retval = setexception(2, 1.0);
 822                 if (lib_version == c_issue_4)
 823                         exc.retval = HUGE;
 824                 else
 825                         exc.retval = HUGE_VAL;
 826                 if (lib_version == strict_ansi)
 827                         errno = ERANGE;
 828                 else if (!matherr(&exc))
 829                         errno = ERANGE;
 830                 break;
 831         case 41:
 832                 /* gamma(-integer) or gamma(0) */
 833                 exc.type = SING;
 834                 exc.name = "gamma";
 835                 ieee_retval = setexception(0, 1.0);
 836                 if (lib_version == c_issue_4)
 837                         exc.retval = HUGE;
 838                 else
 839                         exc.retval = HUGE_VAL;
 840                 if (lib_version == strict_ansi) {
 841                         errno = EDOM;
 842                 } else if (!matherr(&exc)) {
 843                         if (lib_version == c_issue_4) {
 844                                 (void) write(2, "gamma: SING error\n", 18);
 845                         }
 846                         errno = EDOM;
 847                 }
 848                 break;
 849         case 42:
 850                 /* pow(NaN,0.0) */
 851                 /* error if lib_version == c_issue_4 or ansi_1 */
 852                 exc.type = DOMAIN;
 853                 exc.name = "pow";
 854                 exc.retval = x;
 855                 ieee_retval = 1.0;
 856                 if (lib_version == strict_ansi) {
 857                         exc.retval = 1.0;
 858                 } else if (!matherr(&exc)) {
 859                         if ((lib_version == c_issue_4) || (lib_version == ansi_1))
 860                                 errno = EDOM;
 861                 }
 862                 break;
 863         case 43:
 864                 /* log1p(-1) */
 865                 exc.type = SING;
 866                 exc.name = "log1p";
 867                 ieee_retval = setexception(0, -1.0);
 868                 if (lib_version == c_issue_4)
 869                         exc.retval = -HUGE;
 870                 else
 871                         exc.retval = -HUGE_VAL;
 872                 if (lib_version == strict_ansi) {
 873                         errno = ERANGE;
 874                 } else if (!matherr(&exc)) {
 875                         if (lib_version == c_issue_4) {
 876                                 (void) write(2, "log1p: SING error\n", 18);
 877                                 errno = EDOM;
 878                         } else {
 879                                 errno = ERANGE;
 880                         }
 881                 }
 882                 break;
 883         case 44:
 884                 /* log1p(x<-1) */
 885                 exc.type = DOMAIN;
 886                 exc.name = "log1p";
 887                 ieee_retval = setexception(3, 1.0);
 888                 exc.retval = ieee_retval;
 889                 if (lib_version == strict_ansi) {
 890                         errno = EDOM;
 891                 } else if (!matherr(&exc)) {
 892                         if (lib_version == c_issue_4) {
 893                                 (void) write(2, "log1p: DOMAIN error\n", 20);
 894                         }
 895                         errno = EDOM;
 896                 }
 897                 break;
 898         case 45:
 899                 /* logb(0) */
 900                 exc.type = DOMAIN;
 901                 exc.name = "logb";
 902                 ieee_retval = setexception(0, -1.0);
 903                 exc.retval = -HUGE_VAL;
 904                 if (lib_version == strict_ansi)
 905                         errno = EDOM;
 906                 else if (!matherr(&exc))
 907                         errno = EDOM;
 908                 break;
 909         case 46:
 910                 /* nextafter overflow */
 911                 exc.type = OVERFLOW;
 912                 exc.name = "nextafter";
 913                 /*
 914                  * The value as returned by setexception is +/-DBL_MAX in
 915                  * round-to-{zero,-/+Inf} mode respectively, which is not
 916                  * usable.
 917                  */
 918                 (void) setexception(2, x);
 919                 ieee_retval = x > 0 ? Inf : -Inf;
 920                 exc.retval = x > 0 ? HUGE_VAL : -HUGE_VAL;
 921                 if (lib_version == strict_ansi)
 922                         errno = ERANGE;
 923                 else if (!matherr(&exc))
 924                         errno = ERANGE;
 925                 break;
 926         case 47:
 927                 /* scalb(x,inf) */
 928                 iy = ((int *)&y)[HIWORD];
 929                 if (lib_version == c_issue_4)
 930                         /* SVID3: ERANGE in all cases */
 931                         errno = ERANGE;
 932                 else if ((x == 0.0 && iy > 0) || (!finite(x) && iy < 0))
 933                         /* EDOM for scalb(0,+inf) or scalb(inf,-inf) */
 934                         errno = EDOM;
 935                 exc.retval = ieee_retval = ((iy < 0)? x / -y : x * y);
 936                 break;
 937         }
 938         switch (lib_version) {
 939         case c_issue_4:
 940         case ansi_1:
 941         case strict_ansi:
 942                 return (exc.retval);
 943                 /* NOTREACHED */
 944         default:
 945                 return (ieee_retval);
 946         }
 947         /* NOTREACHED */
 948 }
 949 
 950 static double
 951 setexception(int n, double x) {
 952         /*
 953          * n =
 954          * 0    division by zero
 955          * 1    underflow
 956          * 2    overflow
 957          * 3    invalid
 958          */
 959         volatile double one = 1.0, zero = 0.0, retv;
 960 
 961         switch (n) {
 962         case 0:         /* division by zero */
 963                 retv = copysign(one / zero, x);
 964                 break;
 965         case 1:         /* underflow */
 966                 retv = DBL_MIN * copysign(DBL_MIN, x);
 967                 break;
 968         case 2:         /* overflow */
 969                 retv = DBL_MAX * copysign(DBL_MAX, x);
 970                 break;
 971         case 3:         /* invalid */
 972                 retv = zero * Inf;      /* for Cheetah */
 973                 break;
 974         }
 975         return (retv);
 976 }