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 }
|