1 /* $OpenBSD: tgamma.c,v 1.3 2016/10/24 15:34:59 otto Exp $ */
2
3 /* Written by Martynas Venckus, 2008, Public domain. */
4
5 #include <err.h>
6 #include <errno.h>
7 #include <math.h>
8 #include <float.h>
9
10 extern int errno;
11
12 #if defined(__vax__)
13 #define _IEEE 0
14 #else
15 #define _IEEE 1
16 #endif
17
18 double
infnan(int iarg)19 infnan(int iarg)
20 {
21 switch (iarg) {
22 case ERANGE:
23 errno = ERANGE;
24 return (HUGE);
25 case -ERANGE:
26 errno = EDOM;
27 return (-HUGE);
28 default:
29 errno = EDOM;
30 return (0);
31 }
32 }
33
34 int
_isinf(double x)35 _isinf(double x)
36 {
37 if (_IEEE) {
38 return isinf(x);
39 }
40 else {
41 return errno == ERANGE;
42 }
43 }
44
45 int
_isnan(double x)46 _isnan(double x)
47 {
48 if (_IEEE) {
49 return isnan(x);
50 }
51 else {
52 return errno == ERANGE;
53 }
54 }
55
56 int
main(void)57 main(void)
58 {
59 double x;
60
61 /* Random values, approx. -177.79..171.63 */
62 x = tgamma(11.0); /* (11 - 1)! */
63 if (floor(x) != 3628800.0)
64 errx(1, "tgamma(11.0) = %f", x);
65
66 x = tgamma(3.5); /* 15/8 * sqrt(pi) */
67 if (floor(x * 100) != 332.0)
68 errx(1, "tgamma(3.5) = %f", x);
69
70 x = tgamma(-0.5); /* -2 * sqrt(pi) */
71 if (floor(x * 100) != -355.0)
72 errx(1, "tgamma(-0.5) = %f", x);
73
74 /* Special cases */
75 x = tgamma(-1); /* Negative integers */
76 if (!_isnan(x))
77 errx(1, "tgamma(-1) = %f", x);
78
79 x = tgamma(-177.8); /* x ~< -177.79 */
80 if (fabs(x) > DBL_EPSILON)
81 errx(1, "tgamma(-177.8) = %f", x);
82
83 x = tgamma(171.64); /* x ~> 171.63 */
84 if (!_isinf(x))
85 errx(1, "tgamma(171.64) = %f", x);
86
87 x = tgamma(0.0);
88 if (!_isinf(x) || x < 0)
89 errx(1, "tgamma(0) = %f", x);
90
91 x = tgamma(-0.0);
92 if (!_isinf(x) || x > 0)
93 errx(1, "tgamma(0) = %f", x);
94
95 x = tgamma(-HUGE_VAL);
96 if (!_isnan(x))
97 errx(1, "tgamma(-HUGE_VAL) = %f", x);
98
99 x = tgamma(HUGE_VAL);
100 if (!_isinf(x))
101 errx(1, "tgamma(HUGE_VAL) = %f", x);
102
103 #ifdef NAN
104 x = tgamma(NAN);
105 if (!_isnan(x))
106 errx(1, "tgamma(NaN) = %f", x);
107 #endif
108
109 return 0;
110 }
111