Print this page
11210 libm should be cstyle(1ONBLD) clean
   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)


  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 }
   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 /*
  23  * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
  24  */
  25 
  26 /*
  27  * Copyright 2005 Sun Microsystems, Inc.  All rights reserved.
  28  * Use is subject to license terms.
  29  */
  30 
  31 #include "libm.h"
  32 #include "xpg6.h"                       /* __xpg6 */
  33 #include <stdio.h>
  34 #include <float.h>                        /* DBL_MAX, DBL_MIN */
  35 #include <unistd.h>                       /* write */
  36 #if defined(__x86)
  37 #include <ieeefp.h>
  38 #undef  fp_class
  39 #define fp_class        fpclass
  40 #define fp_quiet        FP_QNAN
  41 #endif
  42 #include <errno.h>
  43 #undef fflush
  44 #include <sys/isa_defs.h>
  45 
  46 
  47 /*
  48  * Report libm exception error according to System V Interface Definition
  49  * (SVID).
  50  * Error mapping:
  51  *      1 -- acos(|x|>1)
  52  *      2 -- asin(|x|>1)
  53  *      3 -- atan2(+-0,+-0)
  54  *      4 -- hypot overflow
  55  *      5 -- cosh overflow
  56  *      6 -- exp overflow
  57  *      7 -- exp underflow
  58  *      8 -- y0(0)
  59  *      9 -- y0(-ve)
  60  *      10-- y1(0)
  61  *      11-- y1(-ve)
  62  *      12-- yn(0)
  63  *      13-- yn(-ve)
  64  *      14-- lgamma(finite) overflow
  65  *      15-- lgamma(-integer)
  66  *      16-- log(0)


  79  *      29-- acosh(x<1)
  80  *      30-- atanh(|x|>1)
  81  *      31-- atanh(|x|=1)
  82  *      32-- scalb overflow
  83  *      33-- scalb underflow
  84  *      34-- j0(|x|>X_TLOSS)
  85  *      35-- y0(x>X_TLOSS)
  86  *      36-- j1(|x|>X_TLOSS)
  87  *      37-- y1(x>X_TLOSS)
  88  *      38-- jn(|x|>X_TLOSS, n)
  89  *      39-- yn(x>X_TLOSS, n)
  90  *      40-- gamma(finite) overflow
  91  *      41-- gamma(-integer)
  92  *      42-- pow(NaN,0.0) return NaN for SVID/XOPEN
  93  *      43-- log1p(-1)
  94  *      44-- log1p(x<-1)
  95  *      45-- logb(0)
  96  *      46-- nextafter overflow
  97  *      47-- scalb(x,inf)
  98  */

  99 
 100 static double setexception(int, double);
 101 
 102 static const union {
 103         unsigned x[2];
 104         double d;
 105 } C[] = {
 106 #ifdef _LITTLE_ENDIAN
 107         { 0xffffffff, 0x7fffffff },
 108         { 0x54442d18, 0x400921fb },
 109 #else
 110         { 0x7fffffff, 0xffffffff },
 111         { 0x400921fb, 0x54442d18 },
 112 #endif
 113 };
 114 
 115 #define NaN             C[0].d
 116 #define PI_RZ           C[1].d
 117 
 118 #define __HI(x)         ((unsigned *)&x)[HIWORD]
 119 #define __LO(x)         ((unsigned *)&x)[LOWORD]
 120 #undef  Inf
 121 #define Inf             HUGE_VAL
 122 
 123 double
 124 _SVID_libm_err(double x, double y, int type)
 125 {
 126         struct exception exc;
 127         double t, w, ieee_retval = 0;
 128         enum version lib_version = _lib_version;
 129         int iy;
 130 
 131         /* force libm_ieee behavior in SUSv3 mode */
 132         if ((__xpg6 & _C99SUSv3_math_errexcept) != 0)
 133                 lib_version = libm_ieee;
 134 
 135         if (lib_version == c_issue_4)
 136                 (void) fflush(stdout);
 137 
 138         exc.arg1 = x;
 139         exc.arg2 = y;
 140 
 141         switch (type) {
 142         case 1:
 143                 /* acos(|x|>1) */
 144                 exc.type = DOMAIN;
 145                 exc.name = "acos";
 146                 ieee_retval = setexception(3, 1.0);
 147                 exc.retval = 0.0;
 148 
 149                 if (lib_version == strict_ansi) {
 150                         errno = EDOM;
 151                 } else if (!matherr(&exc)) {
 152                         if (lib_version == c_issue_4)
 153                                 (void) write(2, "acos: DOMAIN error\n", 19);
 154 
 155                         errno = EDOM;
 156                 }
 157 
 158                 break;
 159         case 2:
 160                 /* asin(|x|>1) */
 161                 exc.type = DOMAIN;
 162                 exc.name = "asin";
 163                 exc.retval = 0.0;
 164                 ieee_retval = setexception(3, 1.0);
 165 
 166                 if (lib_version == strict_ansi) {
 167                         errno = EDOM;
 168                 } else if (!matherr(&exc)) {
 169                         if (lib_version == c_issue_4)
 170                                 (void) write(2, "asin: DOMAIN error\n", 19);
 171 
 172                         errno = EDOM;
 173                 }
 174 
 175                 break;
 176         case 3:
 177                 /* atan2(+-0,+-0) */
 178                 exc.arg1 = y;
 179                 exc.arg2 = x;
 180                 exc.type = DOMAIN;
 181                 exc.name = "atan2";
 182                 ieee_retval = copysign(1.0, x) == 1.0 ? y : copysign(PI_RZ +
 183                     DBL_MIN, y);
 184                 exc.retval = 0.0;
 185 
 186                 if (lib_version == strict_ansi) {
 187                         errno = EDOM;
 188                 } else if (!matherr(&exc)) {
 189                         if (lib_version == c_issue_4)
 190                                 (void) write(2, "atan2: DOMAIN error\n", 20);
 191 
 192                         errno = EDOM;
 193                 }
 194 
 195                 break;
 196         case 4:
 197                 /* hypot(finite,finite) overflow */
 198                 exc.type = OVERFLOW;
 199                 exc.name = "hypot";
 200                 ieee_retval = Inf;
 201 
 202                 if (lib_version == c_issue_4)
 203                         exc.retval = HUGE;
 204                 else
 205                         exc.retval = HUGE_VAL;
 206 
 207                 if (lib_version == strict_ansi)
 208                         errno = ERANGE;
 209                 else if (!matherr(&exc))
 210                         errno = ERANGE;
 211 
 212                 break;
 213         case 5:
 214                 /* cosh(finite) overflow */
 215                 exc.type = OVERFLOW;
 216                 exc.name = "cosh";
 217                 ieee_retval = setexception(2, 1.0);
 218 
 219                 if (lib_version == c_issue_4)
 220                         exc.retval = HUGE;
 221                 else
 222                         exc.retval = HUGE_VAL;
 223 
 224                 if (lib_version == strict_ansi)
 225                         errno = ERANGE;
 226                 else if (!matherr(&exc))
 227                         errno = ERANGE;
 228 
 229                 break;
 230         case 6:
 231                 /* exp(finite) overflow */
 232                 exc.type = OVERFLOW;
 233                 exc.name = "exp";
 234                 ieee_retval = setexception(2, 1.0);
 235 
 236                 if (lib_version == c_issue_4)
 237                         exc.retval = HUGE;
 238                 else
 239                         exc.retval = HUGE_VAL;
 240 
 241                 if (lib_version == strict_ansi)
 242                         errno = ERANGE;
 243                 else if (!matherr(&exc))
 244                         errno = ERANGE;
 245 
 246                 break;
 247         case 7:
 248                 /* exp(finite) underflow */
 249                 exc.type = UNDERFLOW;
 250                 exc.name = "exp";
 251                 ieee_retval = setexception(1, 1.0);
 252                 exc.retval = 0.0;
 253 
 254                 if (lib_version == strict_ansi)
 255                         errno = ERANGE;
 256                 else if (!matherr(&exc))
 257                         errno = ERANGE;
 258 
 259                 break;
 260         case 8:
 261                 /* y0(0) = -inf */
 262                 exc.type = DOMAIN;      /* should be SING for IEEE */
 263                 exc.name = "y0";
 264                 ieee_retval = setexception(0, -1.0);
 265 
 266                 if (lib_version == c_issue_4)
 267                         exc.retval = -HUGE;
 268                 else
 269                         exc.retval = -HUGE_VAL;
 270 
 271                 if (lib_version == strict_ansi) {
 272                         errno = EDOM;
 273                 } else if (!matherr(&exc)) {
 274                         if (lib_version == c_issue_4)
 275                                 (void) write(2, "y0: DOMAIN error\n", 17);
 276 
 277                         errno = EDOM;
 278                 }
 279 
 280                 break;
 281         case 9:
 282                 /* y0(x<0) = NaN */
 283                 exc.type = DOMAIN;
 284                 exc.name = "y0";
 285                 ieee_retval = setexception(3, 1.0);
 286 
 287                 if (lib_version == c_issue_4)
 288                         exc.retval = -HUGE;
 289                 else
 290                         exc.retval = -HUGE_VAL;
 291 
 292                 if (lib_version == strict_ansi) {
 293                         errno = EDOM;
 294                 } else if (!matherr(&exc)) {
 295                         if (lib_version == c_issue_4)
 296                                 (void) write(2, "y0: DOMAIN error\n", 17);
 297 
 298                         errno = EDOM;
 299                 }
 300 
 301                 break;
 302         case 10:
 303                 /* y1(0) = -inf */
 304                 exc.type = DOMAIN;      /* should be SING for IEEE */
 305                 exc.name = "y1";
 306                 ieee_retval = setexception(0, -1.0);
 307 
 308                 if (lib_version == c_issue_4)
 309                         exc.retval = -HUGE;
 310                 else
 311                         exc.retval = -HUGE_VAL;
 312 
 313                 if (lib_version == strict_ansi) {
 314                         errno = EDOM;
 315                 } else if (!matherr(&exc)) {
 316                         if (lib_version == c_issue_4)
 317                                 (void) write(2, "y1: DOMAIN error\n", 17);
 318 
 319                         errno = EDOM;
 320                 }
 321 
 322                 break;
 323         case 11:
 324                 /* y1(x<0) = NaN */
 325                 exc.type = DOMAIN;
 326                 exc.name = "y1";
 327                 ieee_retval = setexception(3, 1.0);
 328 
 329                 if (lib_version == c_issue_4)
 330                         exc.retval = -HUGE;
 331                 else
 332                         exc.retval = -HUGE_VAL;
 333 
 334                 if (lib_version == strict_ansi) {
 335                         errno = EDOM;
 336                 } else if (!matherr(&exc)) {
 337                         if (lib_version == c_issue_4)
 338                                 (void) write(2, "y1: DOMAIN error\n", 17);
 339 
 340                         errno = EDOM;
 341                 }
 342 
 343                 break;
 344         case 12:
 345                 /* yn(n,0) = -inf */
 346                 exc.type = DOMAIN;      /* should be SING for IEEE */
 347                 exc.name = "yn";
 348                 ieee_retval = setexception(0, -1.0);
 349 
 350                 if (lib_version == c_issue_4)
 351                         exc.retval = -HUGE;
 352                 else
 353                         exc.retval = -HUGE_VAL;
 354 
 355                 if (lib_version == strict_ansi) {
 356                         errno = EDOM;
 357                 } else if (!matherr(&exc)) {
 358                         if (lib_version == c_issue_4)
 359                                 (void) write(2, "yn: DOMAIN error\n", 17);
 360 
 361                         errno = EDOM;
 362                 }
 363 
 364                 break;
 365         case 13:
 366                 /* yn(x<0) = NaN */
 367                 exc.type = DOMAIN;
 368                 exc.name = "yn";
 369                 ieee_retval = setexception(3, 1.0);
 370 
 371                 if (lib_version == c_issue_4)
 372                         exc.retval = -HUGE;
 373                 else
 374                         exc.retval = -HUGE_VAL;
 375 
 376                 if (lib_version == strict_ansi) {
 377                         errno = EDOM;
 378                 } else if (!matherr(&exc)) {
 379                         if (lib_version == c_issue_4)
 380                                 (void) write(2, "yn: DOMAIN error\n", 17);
 381 
 382                         errno = EDOM;
 383                 }
 384 
 385                 break;
 386         case 14:
 387                 /* lgamma(finite) overflow */
 388                 exc.type = OVERFLOW;
 389                 exc.name = "lgamma";
 390                 ieee_retval = setexception(2, 1.0);
 391 
 392                 if (lib_version == c_issue_4)
 393                         exc.retval = HUGE;
 394                 else
 395                         exc.retval = HUGE_VAL;
 396 
 397                 if (lib_version == strict_ansi)
 398                         errno = ERANGE;
 399                 else if (!matherr(&exc))
 400                         errno = ERANGE;
 401 
 402                 break;
 403         case 15:
 404                 /* lgamma(-integer) or lgamma(0) */
 405                 exc.type = SING;
 406                 exc.name = "lgamma";
 407                 ieee_retval = setexception(0, 1.0);
 408 
 409                 if (lib_version == c_issue_4)
 410                         exc.retval = HUGE;
 411                 else
 412                         exc.retval = HUGE_VAL;
 413 
 414                 if (lib_version == strict_ansi) {
 415                         errno = EDOM;
 416                 } else if (!matherr(&exc)) {
 417                         if (lib_version == c_issue_4)
 418                                 (void) write(2, "lgamma: SING error\n", 19);
 419 
 420                         errno = EDOM;
 421                 }
 422 
 423                 break;
 424         case 16:
 425                 /* log(0) */
 426                 exc.type = SING;
 427                 exc.name = "log";
 428                 ieee_retval = setexception(0, -1.0);
 429 
 430                 if (lib_version == c_issue_4)
 431                         exc.retval = -HUGE;
 432                 else
 433                         exc.retval = -HUGE_VAL;
 434 
 435                 if (lib_version == strict_ansi) {
 436                         errno = ERANGE;
 437                 } else if (!matherr(&exc)) {
 438                         if (lib_version == c_issue_4) {
 439                                 (void) write(2, "log: SING error\n", 16);
 440                                 errno = EDOM;
 441                         } else {
 442                                 errno = ERANGE;
 443                         }
 444                 }
 445 
 446                 break;
 447         case 17:
 448                 /* log(x<0) */
 449                 exc.type = DOMAIN;
 450                 exc.name = "log";
 451                 ieee_retval = setexception(3, 1.0);
 452 
 453                 if (lib_version == c_issue_4)
 454                         exc.retval = -HUGE;
 455                 else
 456                         exc.retval = -HUGE_VAL;
 457 
 458                 if (lib_version == strict_ansi) {
 459                         errno = EDOM;
 460                 } else if (!matherr(&exc)) {
 461                         if (lib_version == c_issue_4)
 462                                 (void) write(2, "log: DOMAIN error\n", 18);
 463 
 464                         errno = EDOM;
 465                 }
 466 
 467                 break;
 468         case 18:
 469                 /* log10(0) */
 470                 exc.type = SING;
 471                 exc.name = "log10";
 472                 ieee_retval = setexception(0, -1.0);
 473 
 474                 if (lib_version == c_issue_4)
 475                         exc.retval = -HUGE;
 476                 else
 477                         exc.retval = -HUGE_VAL;
 478 
 479                 if (lib_version == strict_ansi) {
 480                         errno = ERANGE;
 481                 } else if (!matherr(&exc)) {
 482                         if (lib_version == c_issue_4) {
 483                                 (void) write(2, "log10: SING error\n", 18);
 484                                 errno = EDOM;
 485                         } else {
 486                                 errno = ERANGE;
 487                         }
 488                 }
 489 
 490                 break;
 491         case 19:
 492                 /* log10(x<0) */
 493                 exc.type = DOMAIN;
 494                 exc.name = "log10";
 495                 ieee_retval = setexception(3, 1.0);
 496 
 497                 if (lib_version == c_issue_4)
 498                         exc.retval = -HUGE;
 499                 else
 500                         exc.retval = -HUGE_VAL;
 501 
 502                 if (lib_version == strict_ansi) {
 503                         errno = EDOM;
 504                 } else if (!matherr(&exc)) {
 505                         if (lib_version == c_issue_4)
 506                                 (void) write(2, "log10: DOMAIN error\n", 20);
 507 
 508                         errno = EDOM;
 509                 }
 510 
 511                 break;
 512         case 20:
 513 
 514                 /*
 515                  * pow(0.0,0.0)
 516                  * error only if lib_version == c_issue_4
 517                  */
 518                 exc.type = DOMAIN;
 519                 exc.name = "pow";
 520                 exc.retval = 0.0;
 521                 ieee_retval = 1.0;
 522 
 523                 if (lib_version != c_issue_4) {
 524                         exc.retval = 1.0;
 525                 } else if (!matherr(&exc)) {
 526                         (void) write(2, "pow(0,0): DOMAIN error\n", 23);
 527                         errno = EDOM;
 528                 }
 529 
 530                 break;
 531         case 21:
 532                 /* pow(x,y) overflow */
 533                 exc.type = OVERFLOW;
 534                 exc.name = "pow";
 535                 exc.retval = (lib_version == c_issue_4) ? HUGE : HUGE_VAL;
 536 
 537                 if (signbit(x)) {
 538                         t = rint(y);
 539 
 540                         if (t == y) {
 541                                 w = rint(0.5 * y);
 542 
 543                                 if (t != w + w) /* y is odd */
 544                                         exc.retval = -exc.retval;
 545                         }
 546                 }
 547 
 548                 ieee_retval = setexception(2, exc.retval);
 549 
 550                 if (lib_version == strict_ansi)
 551                         errno = ERANGE;
 552                 else if (!matherr(&exc))
 553                         errno = ERANGE;
 554 
 555                 break;
 556         case 22:
 557                 /* pow(x,y) underflow */
 558                 exc.type = UNDERFLOW;
 559                 exc.name = "pow";
 560                 exc.retval = 0.0;
 561 
 562                 if (signbit(x)) {
 563                         t = rint(y);
 564 
 565                         if (t == y) {
 566                                 w = rint(0.5 * y);
 567 
 568                                 if (t != w + w) /* y is odd */
 569                                         exc.retval = -exc.retval;
 570                         }
 571                 }
 572 
 573                 ieee_retval = setexception(1, exc.retval);
 574 
 575                 if (lib_version == strict_ansi)
 576                         errno = ERANGE;
 577                 else if (!matherr(&exc))
 578                         errno = ERANGE;
 579 
 580                 break;
 581         case 23:
 582                 /* (+-0)**neg */
 583                 exc.type = DOMAIN;
 584                 exc.name = "pow";
 585                 ieee_retval = setexception(0, 1.0);
 586                 {
 587                         int ahy, k, j, yisint, ly, hx;
 588 
 589                         /* BEGIN CSTYLED */
 590                         /*
 591                          * determine if y is an odd int when x = -0
 592                          * yisint = 0       ... y is not an integer
 593                          * yisint = 1       ... y is an odd int
 594                          * yisint = 2       ... y is an even int
 595                          */
 596                         /* END CSTYLED */
 597                         hx = __HI(x);
 598                         ahy = __HI(y) & 0x7fffffff;
 599                         ly = __LO(y);
 600 
 601                         yisint = 0;
 602 
 603                         if (ahy >= 0x43400000) {
 604                                 yisint = 2; /* even integer y */
 605                         } else if (ahy >= 0x3ff00000) {
 606                                 k = (ahy >> 20) - 0x3ff;  /* exponent */
 607 
 608                                 if (k > 20) {
 609                                         j = ly >> (52 - k);
 610 
 611                                         if ((j << (52 - k)) == ly)
 612                                                 yisint = 2 - (j & 1);
 613                                 } else if (ly == 0) {
 614                                         j = ahy >> (20 - k);
 615 
 616                                         if ((j << (20 - k)) == ahy)
 617                                                 yisint = 2 - (j & 1);
 618                                 }
 619                         }
 620 
 621                         if (hx < 0 && yisint == 1)
 622                                 ieee_retval = -ieee_retval;
 623                 }
 624 
 625                 if (lib_version == c_issue_4)
 626                         exc.retval = 0.0;
 627                 else
 628                         exc.retval = -HUGE_VAL;
 629 
 630                 if (lib_version == strict_ansi) {
 631                         errno = EDOM;
 632                 } else if (!matherr(&exc)) {
 633                         if (lib_version == c_issue_4) {
 634                                 (void) write(2, "pow(0,neg): DOMAIN error\n",
 635                                     25);
 636                         }
 637 
 638                         errno = EDOM;
 639                 }
 640 
 641                 break;
 642         case 24:
 643                 /* neg**non-integral */
 644                 exc.type = DOMAIN;
 645                 exc.name = "pow";
 646                 ieee_retval = setexception(3, 1.0);
 647 
 648                 if (lib_version == c_issue_4)
 649                         exc.retval = 0.0;
 650                 else
 651                         exc.retval = ieee_retval;       /* X/Open allow NaN */
 652 
 653                 if (lib_version == strict_ansi) {
 654                         errno = EDOM;
 655                 } else if (!matherr(&exc)) {
 656                         if (lib_version == c_issue_4) {
 657                                 (void) write(2,
 658                                     "neg**non-integral: DOMAIN error\n", 32);
 659                         }
 660 
 661                         errno = EDOM;
 662                 }
 663 
 664                 break;
 665         case 25:
 666                 /* sinh(finite) overflow */
 667                 exc.type = OVERFLOW;
 668                 exc.name = "sinh";
 669                 ieee_retval = copysign(Inf, x);
 670 
 671                 if (lib_version == c_issue_4)
 672                         exc.retval = x > 0.0 ? HUGE : -HUGE;
 673                 else
 674                         exc.retval = x > 0.0 ? HUGE_VAL : -HUGE_VAL;
 675 
 676                 if (lib_version == strict_ansi)
 677                         errno = ERANGE;
 678                 else if (!matherr(&exc))
 679                         errno = ERANGE;
 680 
 681                 break;
 682         case 26:
 683                 /* sqrt(x<0) */
 684                 exc.type = DOMAIN;
 685                 exc.name = "sqrt";
 686                 ieee_retval = setexception(3, 1.0);
 687 
 688                 if (lib_version == c_issue_4)
 689                         exc.retval = 0.0;
 690                 else
 691                         exc.retval = ieee_retval;       /* quiet NaN */
 692 
 693                 if (lib_version == strict_ansi) {
 694                         errno = EDOM;
 695                 } else if (!matherr(&exc)) {
 696                         if (lib_version == c_issue_4)
 697                                 (void) write(2, "sqrt: DOMAIN error\n", 19);
 698 
 699                         errno = EDOM;
 700                 }
 701 
 702                 break;
 703         case 27:
 704                 /* fmod(x,0) */
 705                 exc.type = DOMAIN;
 706                 exc.name = "fmod";
 707 
 708                 if (fp_class(x) == fp_quiet)
 709                         ieee_retval = NaN;
 710                 else
 711                         ieee_retval = setexception(3, 1.0);
 712 
 713                 if (lib_version == c_issue_4)
 714                         exc.retval = x;
 715                 else
 716                         exc.retval = ieee_retval;
 717 
 718                 if (lib_version == strict_ansi) {
 719                         errno = EDOM;
 720                 } else if (!matherr(&exc)) {
 721                         if (lib_version == c_issue_4)
 722                                 (void) write(2, "fmod:  DOMAIN error\n", 20);
 723 
 724                         errno = EDOM;
 725                 }
 726 
 727                 break;
 728         case 28:
 729                 /* remainder(x,0) */
 730                 exc.type = DOMAIN;
 731                 exc.name = "remainder";
 732 
 733                 if (fp_class(x) == fp_quiet)
 734                         ieee_retval = NaN;
 735                 else
 736                         ieee_retval = setexception(3, 1.0);
 737 
 738                 exc.retval = NaN;
 739 
 740                 if (lib_version == strict_ansi) {
 741                         errno = EDOM;
 742                 } else if (!matherr(&exc)) {
 743                         if (lib_version == c_issue_4)
 744                                 (void) write(2, "remainder: DOMAIN error\n",
 745                                     24);
 746 
 747                         errno = EDOM;
 748                 }
 749 
 750                 break;
 751         case 29:
 752                 /* acosh(x<1) */
 753                 exc.type = DOMAIN;
 754                 exc.name = "acosh";
 755                 ieee_retval = setexception(3, 1.0);
 756                 exc.retval = NaN;
 757 
 758                 if (lib_version == strict_ansi) {
 759                         errno = EDOM;
 760                 } else if (!matherr(&exc)) {
 761                         if (lib_version == c_issue_4)
 762                                 (void) write(2, "acosh: DOMAIN error\n", 20);
 763 
 764                         errno = EDOM;
 765                 }
 766 
 767                 break;
 768         case 30:
 769                 /* atanh(|x|>1) */
 770                 exc.type = DOMAIN;
 771                 exc.name = "atanh";
 772                 ieee_retval = setexception(3, 1.0);
 773                 exc.retval = NaN;
 774 
 775                 if (lib_version == strict_ansi) {
 776                         errno = EDOM;
 777                 } else if (!matherr(&exc)) {
 778                         if (lib_version == c_issue_4)
 779                                 (void) write(2, "atanh: DOMAIN error\n", 20);
 780 
 781                         errno = EDOM;
 782                 }
 783 
 784                 break;
 785         case 31:
 786                 /* atanh(|x|=1) */
 787                 exc.type = SING;
 788                 exc.name = "atanh";
 789                 ieee_retval = setexception(0, x);
 790                 exc.retval = ieee_retval;
 791 
 792                 if (lib_version == strict_ansi) {
 793                         errno = ERANGE;
 794                 } else if (!matherr(&exc)) {
 795                         if (lib_version == c_issue_4) {
 796                                 (void) write(2, "atanh: SING error\n", 18);
 797                                 errno = EDOM;
 798                         } else {
 799                                 errno = ERANGE;
 800                         }
 801                 }
 802 
 803                 break;
 804         case 32:
 805                 /* scalb overflow; SVID also returns +-HUGE_VAL */
 806                 exc.type = OVERFLOW;
 807                 exc.name = "scalb";
 808                 ieee_retval = setexception(2, x);
 809                 exc.retval = x > 0.0 ? HUGE_VAL : -HUGE_VAL;
 810 
 811                 if (lib_version == strict_ansi)
 812                         errno = ERANGE;
 813                 else if (!matherr(&exc))
 814                         errno = ERANGE;
 815 
 816                 break;
 817         case 33:
 818                 /* scalb underflow */
 819                 exc.type = UNDERFLOW;
 820                 exc.name = "scalb";
 821                 ieee_retval = setexception(1, x);
 822                 exc.retval = ieee_retval;       /* +-0.0 */
 823 
 824                 if (lib_version == strict_ansi)
 825                         errno = ERANGE;
 826                 else if (!matherr(&exc))
 827                         errno = ERANGE;
 828 
 829                 break;
 830         case 34:
 831                 /* j0(|x|>X_TLOSS) */
 832                 exc.type = TLOSS;
 833                 exc.name = "j0";
 834                 exc.retval = 0.0;
 835                 ieee_retval = y;
 836 
 837                 if (lib_version == strict_ansi) {
 838                         errno = ERANGE;
 839                 } else if (!matherr(&exc)) {
 840                         if (lib_version == c_issue_4) {
 841                                 (void) write(2, exc.name, 2);
 842                                 (void) write(2, ": TLOSS error\n", 14);
 843                         }
 844 
 845                         errno = ERANGE;
 846                 }
 847 
 848                 break;
 849         case 35:
 850                 /* y0(x>X_TLOSS) */
 851                 exc.type = TLOSS;
 852                 exc.name = "y0";
 853                 exc.retval = 0.0;
 854                 ieee_retval = y;
 855 
 856                 if (lib_version == strict_ansi) {
 857                         errno = ERANGE;
 858                 } else if (!matherr(&exc)) {
 859                         if (lib_version == c_issue_4) {
 860                                 (void) write(2, exc.name, 2);
 861                                 (void) write(2, ": TLOSS error\n", 14);
 862                         }
 863 
 864                         errno = ERANGE;
 865                 }
 866 
 867                 break;
 868         case 36:
 869                 /* j1(|x|>X_TLOSS) */
 870                 exc.type = TLOSS;
 871                 exc.name = "j1";
 872                 exc.retval = 0.0;
 873                 ieee_retval = y;
 874 
 875                 if (lib_version == strict_ansi) {
 876                         errno = ERANGE;
 877                 } else if (!matherr(&exc)) {
 878                         if (lib_version == c_issue_4) {
 879                                 (void) write(2, exc.name, 2);
 880                                 (void) write(2, ": TLOSS error\n", 14);
 881                         }
 882 
 883                         errno = ERANGE;
 884                 }
 885 
 886                 break;
 887         case 37:
 888                 /* y1(x>X_TLOSS) */
 889                 exc.type = TLOSS;
 890                 exc.name = "y1";
 891                 exc.retval = 0.0;
 892                 ieee_retval = y;
 893 
 894                 if (lib_version == strict_ansi) {
 895                         errno = ERANGE;
 896                 } else if (!matherr(&exc)) {
 897                         if (lib_version == c_issue_4) {
 898                                 (void) write(2, exc.name, 2);
 899                                 (void) write(2, ": TLOSS error\n", 14);
 900                         }
 901 
 902                         errno = ERANGE;
 903                 }
 904 
 905                 break;
 906         case 38:
 907 
 908                 /*
 909                  * jn(|x|>X_TLOSS)
 910                  * incorrect ieee value: ieee should never be here
 911                  */
 912                 exc.type = TLOSS;
 913                 exc.name = "jn";
 914                 exc.retval = 0.0;
 915                 ieee_retval = 0.0;      /* shall not be used */
 916 
 917                 if (lib_version == strict_ansi) {
 918                         errno = ERANGE;
 919                 } else if (!matherr(&exc)) {
 920                         if (lib_version == c_issue_4) {
 921                                 (void) write(2, exc.name, 2);
 922                                 (void) write(2, ": TLOSS error\n", 14);
 923                         }
 924 
 925                         errno = ERANGE;
 926                 }
 927 
 928                 break;
 929         case 39:
 930 
 931                 /*
 932                  * yn(x>X_TLOSS)
 933                  * incorrect ieee value: ieee should never be here
 934                  */
 935                 exc.type = TLOSS;
 936                 exc.name = "yn";
 937                 exc.retval = 0.0;
 938                 ieee_retval = 0.0;      /* shall not be used */
 939 
 940                 if (lib_version == strict_ansi) {
 941                         errno = ERANGE;
 942                 } else if (!matherr(&exc)) {
 943                         if (lib_version == c_issue_4) {
 944                                 (void) write(2, exc.name, 2);
 945                                 (void) write(2, ": TLOSS error\n", 14);
 946                         }
 947 
 948                         errno = ERANGE;
 949                 }
 950 
 951                 break;
 952         case 40:
 953                 /* gamma(finite) overflow */
 954                 exc.type = OVERFLOW;
 955                 exc.name = "gamma";
 956                 ieee_retval = setexception(2, 1.0);
 957 
 958                 if (lib_version == c_issue_4)
 959                         exc.retval = HUGE;
 960                 else
 961                         exc.retval = HUGE_VAL;
 962 
 963                 if (lib_version == strict_ansi)
 964                         errno = ERANGE;
 965                 else if (!matherr(&exc))
 966                         errno = ERANGE;
 967 
 968                 break;
 969         case 41:
 970                 /* gamma(-integer) or gamma(0) */
 971                 exc.type = SING;
 972                 exc.name = "gamma";
 973                 ieee_retval = setexception(0, 1.0);
 974 
 975                 if (lib_version == c_issue_4)
 976                         exc.retval = HUGE;
 977                 else
 978                         exc.retval = HUGE_VAL;
 979 
 980                 if (lib_version == strict_ansi) {
 981                         errno = EDOM;
 982                 } else if (!matherr(&exc)) {
 983                         if (lib_version == c_issue_4)
 984                                 (void) write(2, "gamma: SING error\n", 18);
 985 
 986                         errno = EDOM;
 987                 }
 988 
 989                 break;
 990         case 42:
 991 
 992                 /*
 993                  * pow(NaN,0.0)
 994                  * error if lib_version == c_issue_4 or ansi_1
 995                  */
 996                 exc.type = DOMAIN;
 997                 exc.name = "pow";
 998                 exc.retval = x;
 999                 ieee_retval = 1.0;
1000 
1001                 if (lib_version == strict_ansi) {
1002                         exc.retval = 1.0;
1003                 } else if (!matherr(&exc)) {
1004                         if ((lib_version == c_issue_4) || (lib_version ==
1005                             ansi_1))
1006                                 errno = EDOM;
1007                 }
1008 
1009                 break;
1010         case 43:
1011                 /* log1p(-1) */
1012                 exc.type = SING;
1013                 exc.name = "log1p";
1014                 ieee_retval = setexception(0, -1.0);
1015 
1016                 if (lib_version == c_issue_4)
1017                         exc.retval = -HUGE;
1018                 else
1019                         exc.retval = -HUGE_VAL;
1020 
1021                 if (lib_version == strict_ansi) {
1022                         errno = ERANGE;
1023                 } else if (!matherr(&exc)) {
1024                         if (lib_version == c_issue_4) {
1025                                 (void) write(2, "log1p: SING error\n", 18);
1026                                 errno = EDOM;
1027                         } else {
1028                                 errno = ERANGE;
1029                         }
1030                 }
1031 
1032                 break;
1033         case 44:
1034                 /* log1p(x<-1) */
1035                 exc.type = DOMAIN;
1036                 exc.name = "log1p";
1037                 ieee_retval = setexception(3, 1.0);
1038                 exc.retval = ieee_retval;
1039 
1040                 if (lib_version == strict_ansi) {
1041                         errno = EDOM;
1042                 } else if (!matherr(&exc)) {
1043                         if (lib_version == c_issue_4)
1044                                 (void) write(2, "log1p: DOMAIN error\n", 20);
1045 
1046                         errno = EDOM;
1047                 }
1048 
1049                 break;
1050         case 45:
1051                 /* logb(0) */
1052                 exc.type = DOMAIN;
1053                 exc.name = "logb";
1054                 ieee_retval = setexception(0, -1.0);
1055                 exc.retval = -HUGE_VAL;
1056 
1057                 if (lib_version == strict_ansi)
1058                         errno = EDOM;
1059                 else if (!matherr(&exc))
1060                         errno = EDOM;
1061 
1062                 break;
1063         case 46:
1064                 /* nextafter overflow */
1065                 exc.type = OVERFLOW;
1066                 exc.name = "nextafter";
1067 
1068                 /*
1069                  * The value as returned by setexception is +/-DBL_MAX in
1070                  * round-to-{zero,-/+Inf} mode respectively, which is not
1071                  * usable.
1072                  */
1073                 (void) setexception(2, x);
1074                 ieee_retval = x > 0 ? Inf : -Inf;
1075                 exc.retval = x > 0 ? HUGE_VAL : -HUGE_VAL;
1076 
1077                 if (lib_version == strict_ansi)
1078                         errno = ERANGE;
1079                 else if (!matherr(&exc))
1080                         errno = ERANGE;
1081 
1082                 break;
1083         case 47:
1084                 /* scalb(x,inf) */
1085                 iy = ((int *)&y)[HIWORD];
1086 
1087                 if (lib_version == c_issue_4)
1088                         /* SVID3: ERANGE in all cases */
1089                         errno = ERANGE;
1090                 else if ((x == 0.0 && iy > 0) || (!finite(x) && iy < 0))
1091                         /* EDOM for scalb(0,+inf) or scalb(inf,-inf) */
1092                         errno = EDOM;
1093 
1094                 exc.retval = ieee_retval = ((iy < 0) ? x / -y : x * y);
1095                 break;
1096         }
1097 
1098         switch (lib_version) {
1099         case c_issue_4:
1100         case ansi_1:
1101         case strict_ansi:
1102                 return (exc.retval);
1103         /* NOTREACHED */
1104         default:
1105                 return (ieee_retval);
1106         }
1107 
1108         /* NOTREACHED */
1109 }
1110 
1111 static double
1112 setexception(int n, double x)
1113 {
1114         /*
1115          * n =
1116          * 0    division by zero
1117          * 1    underflow
1118          * 2    overflow
1119          * 3    invalid
1120          */
1121         volatile double one = 1.0, zero = 0.0, retv;
1122 
1123         switch (n) {
1124         case 0:                         /* division by zero */
1125                 retv = copysign(one / zero, x);
1126                 break;
1127         case 1:                         /* underflow */
1128                 retv = DBL_MIN * copysign(DBL_MIN, x);
1129                 break;
1130         case 2:                         /* overflow */
1131                 retv = DBL_MAX * copysign(DBL_MAX, x);
1132                 break;
1133         case 3:                         /* invalid */
1134                 retv = zero * Inf;      /* for Cheetah */
1135                 break;
1136         }
1137 
1138         return (retv);
1139 }