| 571 |
ihaka |
1 |
#ifndef MATHLIB_H
|
|
|
2 |
#define MATHLIB_H
|
| 2 |
r |
3 |
|
|
|
4 |
#include "Arith.h"
|
|
|
5 |
|
| 571 |
ihaka |
6 |
#include <errno.h>
|
|
|
7 |
#include <float.h>
|
|
|
8 |
#include <limits.h>
|
| 2 |
r |
9 |
#include <math.h>
|
|
|
10 |
#include <stdlib.h>
|
|
|
11 |
|
|
|
12 |
/* 30 Decimal-place constants computed with bc (scale=32; proper round) */
|
|
|
13 |
|
|
|
14 |
#ifndef M_1_SQRT_2
|
| 571 |
ihaka |
15 |
#define M_1_SQRT_2 0.707106781186547524400844362105
|
| 2 |
r |
16 |
#endif
|
|
|
17 |
|
|
|
18 |
#ifndef M_PI
|
| 571 |
ihaka |
19 |
#define M_PI 3.141592653589793238462643383279
|
| 2 |
r |
20 |
#endif
|
|
|
21 |
|
|
|
22 |
#ifndef M_PI_half
|
| 571 |
ihaka |
23 |
#define M_PI_half 1.570796326794896619231321691640
|
| 2 |
r |
24 |
#endif
|
|
|
25 |
|
|
|
26 |
#ifndef M_SQRT_PI
|
| 571 |
ihaka |
27 |
#define M_SQRT_PI 1.772453850905516027298167483341
|
| 2 |
r |
28 |
#endif
|
|
|
29 |
|
|
|
30 |
#ifndef M_1_SQRT_2PI
|
| 571 |
ihaka |
31 |
#define M_1_SQRT_2PI 0.398942280401432677939946059934
|
| 2 |
r |
32 |
#endif
|
|
|
33 |
|
|
|
34 |
/* log(sqrt(2*pi)) = log(2*pi)/2 : */
|
|
|
35 |
#ifndef M_LN_SQRT_2PI
|
| 571 |
ihaka |
36 |
#define M_LN_SQRT_2PI 0.918938533204672741780329736406
|
| 2 |
r |
37 |
#endif
|
|
|
38 |
|
| 571 |
ihaka |
39 |
#define MATHLIB_ERROR(x) { printf("%s\n",x); exit(1); }
|
|
|
40 |
#define MATHLIB_WARNING(x) { printf("%s\n",x); }
|
| 2 |
r |
41 |
|
| 571 |
ihaka |
42 |
#define ME_NONE 0
|
|
|
43 |
#define ME_DOMAIN 1
|
|
|
44 |
#define ME_RANGE 2
|
|
|
45 |
#define ME_NOCONV 3
|
|
|
46 |
#define ME_PRECISION 4
|
|
|
47 |
#define ME_UNDERFLOW 5
|
|
|
48 |
|
|
|
49 |
#undef ML_PRECISION_WARNINGS
|
|
|
50 |
|
|
|
51 |
#ifdef IEEE_754
|
|
|
52 |
extern double m_zero;
|
|
|
53 |
extern double m_one;
|
|
|
54 |
extern double m_tiny;
|
|
|
55 |
#define ML_ERROR(x) /* nothing */
|
|
|
56 |
#define ML_POSINF (m_one / m_zero)
|
|
|
57 |
#define ML_NEGINF ((-m_one) / m_zero)
|
|
|
58 |
#define ML_NAN (m_zero / m_zero)
|
|
|
59 |
#define ML_UNDERFLOW (m_tiny * m_tiny)
|
|
|
60 |
#define ML_VALID(x) (!isnan(x))
|
|
|
61 |
#else
|
|
|
62 |
#define ML_ERROR(x) ml_error(x)
|
|
|
63 |
#define ML_POSINF DBL_MAX
|
|
|
64 |
#define ML_NEGINF (-DBL_MAX)
|
|
|
65 |
#define ML_NAN (-DBL_MAX)
|
|
|
66 |
#define ML_UNDERFLOW 0
|
|
|
67 |
#define ML_VALID(x) (errno == 0)
|
| 2 |
r |
68 |
#endif
|
|
|
69 |
|
| 571 |
ihaka |
70 |
/* Splus Compatibility */
|
| 2 |
r |
71 |
|
| 571 |
ihaka |
72 |
#define snorm norm_rand
|
|
|
73 |
#define sunif unif_rand
|
|
|
74 |
#define sexp exp_rand
|
| 2 |
r |
75 |
|
| 571 |
ihaka |
76 |
/* Name Hiding to Avoid Clashes with Fortran */
|
| 2 |
r |
77 |
|
| 571 |
ihaka |
78 |
#ifdef HIDE_NAMES
|
|
|
79 |
#define d1mach c_d1mach
|
|
|
80 |
#define i1mach c_i1mach
|
|
|
81 |
#endif
|
| 2 |
r |
82 |
|
| 571 |
ihaka |
83 |
#define rround fround
|
|
|
84 |
#define prec fprec
|
|
|
85 |
#define trunc ftrunc
|
|
|
86 |
#define sign fsign
|
| 2 |
r |
87 |
|
| 571 |
ihaka |
88 |
/* Machine Characteristics */
|
| 2 |
r |
89 |
|
| 571 |
ihaka |
90 |
double d1mach(int);
|
|
|
91 |
double d1mach_(int*);
|
|
|
92 |
int i1mach(int);
|
|
|
93 |
int i1mach_(int*);
|
| 2 |
r |
94 |
|
| 571 |
ihaka |
95 |
/* General Support Functions */
|
| 2 |
r |
96 |
|
| 571 |
ihaka |
97 |
int imax2(int, int);
|
|
|
98 |
int imin2(int, int);
|
|
|
99 |
double fmax2(double, double);
|
|
|
100 |
double fmin2(double, double);
|
|
|
101 |
double fmod(double, double);
|
|
|
102 |
double fprec(double, double);
|
|
|
103 |
double fround(double, double);
|
|
|
104 |
double ftrunc(double);
|
|
|
105 |
double fsign(double, double);
|
|
|
106 |
double fsquare(double);
|
|
|
107 |
double fcube(double);
|
| 2 |
r |
108 |
|
| 571 |
ihaka |
109 |
/* Random Number Generators */
|
| 2 |
r |
110 |
|
| 571 |
ihaka |
111 |
double snorm(void);
|
|
|
112 |
double sunif(void);
|
|
|
113 |
double sexp(void);
|
| 2 |
r |
114 |
|
| 571 |
ihaka |
115 |
/* Chebyshev Series */
|
| 2 |
r |
116 |
|
| 571 |
ihaka |
117 |
int chebyshev_init(double*, int, double);
|
|
|
118 |
double chebyshev_eval(double, double *, int);
|
| 2 |
r |
119 |
|
| 571 |
ihaka |
120 |
/* Gamma and Related Functions */
|
| 2 |
r |
121 |
|
| 571 |
ihaka |
122 |
double logrelerr(double);
|
|
|
123 |
void gammalims(double*, double*);
|
|
|
124 |
double lgammacor(double);
|
|
|
125 |
double gamma(double);
|
|
|
126 |
double lgamma(double);
|
|
|
127 |
void dpsifn(double, int, int, int, double*, int*, int*);
|
|
|
128 |
double digamma(double);
|
|
|
129 |
double trigamma(double);
|
|
|
130 |
double tetragamma(double);
|
|
|
131 |
double pentagamma(double);
|
| 2 |
r |
132 |
|
| 571 |
ihaka |
133 |
double choose(double, double);
|
|
|
134 |
double lchoose(double, double);
|
|
|
135 |
double fastchoose(double, double);
|
|
|
136 |
double lfastchoose(double, double);
|
| 2 |
r |
137 |
|
| 571 |
ihaka |
138 |
/* Beta and Related Functions */
|
| 2 |
r |
139 |
|
| 571 |
ihaka |
140 |
double beta(double, double);
|
|
|
141 |
double lbeta(double, double);
|
| 2 |
r |
142 |
|
| 571 |
ihaka |
143 |
/* Normal Distribution */
|
| 2 |
r |
144 |
|
| 571 |
ihaka |
145 |
double dnorm(double, double, double);
|
|
|
146 |
double pnorm(double, double, double);
|
|
|
147 |
double qnorm(double, double, double);
|
|
|
148 |
double rnorm(double, double);
|
| 2 |
r |
149 |
|
| 571 |
ihaka |
150 |
/* Uniform Distribution */
|
| 2 |
r |
151 |
|
| 571 |
ihaka |
152 |
double dunif(double, double, double);
|
|
|
153 |
double punif(double, double, double);
|
|
|
154 |
double qunif(double, double, double);
|
|
|
155 |
double runif(double, double);
|
| 2 |
r |
156 |
|
| 571 |
ihaka |
157 |
/* Gamma Distribution */
|
| 2 |
r |
158 |
|
| 571 |
ihaka |
159 |
double dgamma(double, double, double);
|
|
|
160 |
double pgamma(double, double, double);
|
|
|
161 |
double qgamma(double, double, double);
|
|
|
162 |
double rgamma(double, double);
|
| 2 |
r |
163 |
|
| 571 |
ihaka |
164 |
/* Beta Distribution */
|
| 2 |
r |
165 |
|
| 571 |
ihaka |
166 |
double dbeta(double, double, double);
|
|
|
167 |
double pbeta(double, double, double);
|
|
|
168 |
double pbeta_raw(double, double, double);
|
|
|
169 |
double qbeta(double, double, double);
|
|
|
170 |
double rbeta(double, double);
|
| 2 |
r |
171 |
|
| 571 |
ihaka |
172 |
/* Lognormal Distribution */
|
| 2 |
r |
173 |
|
| 571 |
ihaka |
174 |
double dlnorm(double, double, double);
|
|
|
175 |
double plnorm(double, double, double);
|
|
|
176 |
double qlnorm(double, double, double);
|
|
|
177 |
double rlnorm(double, double);
|
| 2 |
r |
178 |
|
| 571 |
ihaka |
179 |
/* Chi-squared Distribution */
|
| 2 |
r |
180 |
|
| 571 |
ihaka |
181 |
double dchisq(double, double);
|
|
|
182 |
double pchisq(double, double);
|
|
|
183 |
double qchisq(double, double);
|
|
|
184 |
double rchisq(double);
|
| 2 |
r |
185 |
|
| 571 |
ihaka |
186 |
/* Non-central Chi-squared Distribution */
|
| 2 |
r |
187 |
|
| 605 |
ihaka |
188 |
double dnchisq(double, double, double);
|
| 571 |
ihaka |
189 |
double pnchisq(double, double, double);
|
|
|
190 |
double qnchisq(double, double, double);
|
| 605 |
ihaka |
191 |
double rnchisq(double, double);
|
| 571 |
ihaka |
192 |
|
|
|
193 |
/* F Distibution */
|
|
|
194 |
|
|
|
195 |
double df(double, double, double);
|
|
|
196 |
double pf(double, double, double);
|
|
|
197 |
double qf(double, double, double);
|
|
|
198 |
double rf(double, double);
|
|
|
199 |
|
|
|
200 |
/* Student t Distibution */
|
|
|
201 |
|
|
|
202 |
double dt(double, double);
|
|
|
203 |
double pt(double, double);
|
|
|
204 |
double qt(double, double);
|
|
|
205 |
double rt(double);
|
|
|
206 |
|
|
|
207 |
/* Binomial Distribution */
|
|
|
208 |
|
|
|
209 |
double dbinom(double, double, double);
|
|
|
210 |
double pbinom(double, double, double);
|
|
|
211 |
double qbinom(double, double, double);
|
|
|
212 |
double rbinom(double, double);
|
|
|
213 |
|
|
|
214 |
/* Cauchy Distribution */
|
|
|
215 |
|
|
|
216 |
double dcauchy(double, double, double);
|
|
|
217 |
double pcauchy(double, double, double);
|
|
|
218 |
double qcauchy(double, double, double);
|
|
|
219 |
double rcauchy(double, double);
|
|
|
220 |
|
|
|
221 |
/* Exponential Distribution */
|
|
|
222 |
|
|
|
223 |
double dexp(double, double);
|
|
|
224 |
double pexp(double, double);
|
|
|
225 |
double qexp(double, double);
|
|
|
226 |
double rexp(double);
|
|
|
227 |
|
|
|
228 |
/* Geometric Distribution */
|
|
|
229 |
|
|
|
230 |
double dgeom(double, double);
|
|
|
231 |
double pgeom(double, double);
|
|
|
232 |
double qgeom(double, double);
|
|
|
233 |
double rgeom(double);
|
|
|
234 |
|
|
|
235 |
/* Hypergeometric Distibution */
|
|
|
236 |
|
|
|
237 |
double dhyper(double, double, double, double);
|
|
|
238 |
double phyper(double, double, double, double);
|
|
|
239 |
double qhyper(double, double, double, double);
|
|
|
240 |
double rhyper(double, double, double);
|
|
|
241 |
|
|
|
242 |
/* Negative Binomial Distribution */
|
|
|
243 |
|
|
|
244 |
double dnbinom(double, double, double);
|
|
|
245 |
double pnbinom(double, double, double);
|
|
|
246 |
double qnbinom(double, double, double);
|
|
|
247 |
double rnbinom(double, double);
|
|
|
248 |
|
|
|
249 |
/* Poisson Distribution */
|
|
|
250 |
|
|
|
251 |
double dpois(double, double);
|
|
|
252 |
double ppois(double, double);
|
|
|
253 |
double qpois(double, double);
|
|
|
254 |
double rpois(double);
|
|
|
255 |
|
|
|
256 |
/* Weibull Distribution */
|
|
|
257 |
|
|
|
258 |
double dweibull(double, double, double);
|
|
|
259 |
double pweibull(double, double, double);
|
|
|
260 |
double qweibull(double, double, double);
|
|
|
261 |
double rweibull(double, double);
|
|
|
262 |
|
|
|
263 |
/* Logistic Distribution */
|
|
|
264 |
|
|
|
265 |
double dlogis(double, double, double);
|
|
|
266 |
double plogis(double, double, double);
|
|
|
267 |
double qlogis(double, double, double);
|
|
|
268 |
double rlogis(double, double);
|
|
|
269 |
|
| 602 |
ihaka |
270 |
/* Non-central Beta Distribution */
|
|
|
271 |
|
|
|
272 |
double dnbeta(double, double, double, double);
|
|
|
273 |
double pnbeta(double, double, double, double);
|
|
|
274 |
double qnbeta(double, double, double, double);
|
|
|
275 |
double rnbeta(double, double, double);
|
|
|
276 |
|
|
|
277 |
/* Non-central F Distribution */
|
|
|
278 |
|
|
|
279 |
double dnf(double, double, double, double);
|
|
|
280 |
double pnf(double, double, double, double);
|
|
|
281 |
double qnf(double, double, double, double);
|
|
|
282 |
double rnf(double, double, double);
|
|
|
283 |
|
|
|
284 |
/* Non-central Student t Distribution */
|
|
|
285 |
|
|
|
286 |
double dnt(double, double, double);
|
|
|
287 |
double pnt(double, double, double);
|
|
|
288 |
double qnt(double, double, double);
|
|
|
289 |
double rnt(double, double);
|
|
|
290 |
|
| 624 |
ihaka |
291 |
/* Studentized Range Distribution */
|
|
|
292 |
|
|
|
293 |
double dtukey(double, double, double, double);
|
|
|
294 |
double ptukey(double, double, double, double);
|
|
|
295 |
double qtukey(double, double, double, double);
|
|
|
296 |
double rtukey(double, double, double);
|
|
|
297 |
|
| 2 |
r |
298 |
#endif
|