Commit | Line | Data |
---|---|---|
9432a5c3 KT |
1 | #include "apl.h" |
2 | ||
3 | data | |
4 | ex_add(d1, d2) | |
5 | data d1, d2; | |
6 | { | |
7 | ||
8 | if(fuzz(d1, -d2) == 0) | |
9 | return(zero); | |
10 | return(d1 + d2); | |
11 | } | |
12 | ||
13 | data | |
14 | ex_sub(d1, d2) | |
15 | data d1, d2; | |
16 | { | |
17 | ||
18 | if(fuzz(d1, d2) == 0) | |
19 | return(zero); | |
20 | return(d1 - d2); | |
21 | } | |
22 | ||
23 | data | |
24 | ex_mul(d1, d2) | |
25 | data d1, d2; | |
26 | { | |
27 | return(d1 * d2); | |
28 | } | |
29 | ||
30 | data | |
31 | ex_div(d1, d2) | |
32 | data d1, d2; | |
33 | { | |
34 | ||
35 | /* 0 div 0 is NOT 1 */ | |
36 | if(d2 == zero) | |
37 | error("div D"); | |
38 | return(d1 / d2); | |
39 | } | |
40 | ||
41 | data | |
42 | ex_mod(d1, d2) | |
43 | data d1, d2; | |
44 | { | |
45 | double f; | |
46 | ||
47 | /* see 0 div 0 comment */ | |
48 | if(d1 == zero) | |
49 | error("mod D"); | |
50 | if(d1 < zero) | |
51 | d1 = -d1; | |
52 | f = d2; | |
53 | d2 = d2 - d1 * floor(f/d1); | |
54 | return(d2); | |
55 | } | |
56 | ||
57 | data | |
58 | ex_min(d1, d2) | |
59 | data d1, d2; | |
60 | { | |
61 | ||
62 | if(d1 < d2) | |
63 | return(d1); | |
64 | return(d2); | |
65 | } | |
66 | ||
67 | data | |
68 | ex_max(d1, d2) | |
69 | data d1, d2; | |
70 | { | |
71 | ||
72 | if(d1 > d2) | |
73 | return(d1); | |
74 | return(d2); | |
75 | } | |
76 | ||
77 | data | |
78 | ex_pwr(d1, d2) | |
79 | data d1, d2; | |
80 | { | |
81 | register s; | |
82 | double f1; | |
83 | ||
84 | s = 0; | |
85 | f1 = d1; | |
86 | if(f1 > 0.) { | |
87 | f1 = d2 * log(f1); | |
88 | goto chk; | |
89 | } | |
90 | if(f1 == 0.) { | |
91 | if(d2 == zero) | |
92 | goto bad; | |
93 | return(zero); | |
94 | } | |
95 | /* should check rational d2 here */ | |
96 | goto bad; | |
97 | ||
98 | chk: | |
99 | if(f1 < maxexp) { | |
100 | d1 = exp(f1); | |
101 | if(s) | |
102 | d1 = -d1; | |
103 | return(d1); | |
104 | } | |
105 | bad: | |
106 | error("pwr D"); | |
107 | } | |
108 | ||
109 | data | |
110 | ex_log(d1, d2) | |
111 | data d1, d2; | |
112 | { | |
113 | double f1, f2; | |
114 | ||
115 | f1 = d1; | |
116 | f2 = d2; | |
117 | if(f1 <= 0. || f2 <= 0.) | |
118 | error("log D"); | |
119 | d1 = log(f2)/log(f1); | |
120 | return(d1); | |
121 | } | |
122 | ||
123 | data | |
124 | ex_cir(d1, d2) | |
125 | data d1, d2; | |
126 | { | |
127 | double f, f1; | |
128 | ||
129 | switch(fix(d1)) { | |
130 | ||
131 | default: | |
132 | error("crc D"); | |
133 | ||
134 | case -7: | |
135 | /* arc tanh */ | |
136 | f = d2; | |
137 | if(f <= -1. || f >= 1.) | |
138 | error("tanh D"); | |
139 | f = log((1.+f)/(1.-f))/2.; | |
140 | goto ret; | |
141 | ||
142 | case -6: | |
143 | /* arc cosh */ | |
144 | f = d2; | |
145 | if(f < 1.) | |
146 | error("cosh D"); | |
147 | f = log(f+sqrt(f*f-1.)); | |
148 | goto ret; | |
149 | ||
150 | case -5: | |
151 | /* arc sinh */ | |
152 | f = d2; | |
153 | f = log(f+sqrt(f*f+1.)); | |
154 | goto ret; | |
155 | ||
156 | case -4: | |
157 | /* sqrt(-1 + x*x) */ | |
158 | f = -one + d2*d2; | |
159 | goto sq; | |
160 | ||
161 | case -3: | |
162 | /* arc tan */ | |
163 | f = d2; | |
164 | f = atan(f); | |
165 | goto ret; | |
166 | ||
167 | case -2: | |
168 | /* arc cos */ | |
169 | f = d2; | |
170 | if(f < -1. || f > 1.) | |
171 | error("arc cos D"); | |
172 | f = atan2(sqrt(1.-f*f), f); | |
173 | goto ret; | |
174 | ||
175 | case -1: | |
176 | /* arc sin */ | |
177 | f = d2; | |
178 | if(f < -1. || f > 1.) | |
179 | error("arc sin D"); | |
180 | f = atan2(f, sqrt(1.-f*f)); | |
181 | goto ret; | |
182 | ||
183 | case 0: | |
184 | /* sqrt(1 - x*x) */ | |
185 | f = one - d2*d2; | |
186 | goto sq; | |
187 | ||
188 | case 1: | |
189 | /* sin */ | |
190 | f = d2; | |
191 | f = sin(f); | |
192 | goto ret; | |
193 | ||
194 | case 2: | |
195 | /* cos */ | |
196 | f = d2; | |
197 | f = cos(f); | |
198 | goto ret; | |
199 | ||
200 | case 3: | |
201 | /* tan */ | |
202 | f = d2; | |
203 | f1 = cos(f); | |
204 | if(f1 == 0.) | |
205 | f = exp(maxexp); else | |
206 | f = sin(f)/f1; | |
207 | goto ret; | |
208 | ||
209 | case 4: | |
210 | /* sqrt(1 + x*x) */ | |
211 | f = one + d2*d2; | |
212 | goto sq; | |
213 | ||
214 | case 5: | |
215 | /* sinh */ | |
216 | f = d2; | |
217 | if(f < -maxexp || f > maxexp) | |
218 | error("sinh D"); | |
219 | f = exp(f); | |
220 | f = (f-1./f)/2.; | |
221 | goto ret; | |
222 | ||
223 | case 6: | |
224 | /* cosh */ | |
225 | f = d2; | |
226 | if(f < -maxexp || f > maxexp) | |
227 | error("cosh D"); | |
228 | f = exp(f); | |
229 | f = (f+1./f)/2.; | |
230 | goto ret; | |
231 | ||
232 | case 7: | |
233 | /* tanh */ | |
234 | f = d2; | |
235 | if(f > maxexp) { | |
236 | f = 1.; | |
237 | goto ret; | |
238 | } | |
239 | if(f < -maxexp) { | |
240 | f = -1.; | |
241 | goto ret; | |
242 | } | |
243 | f = exp(f); | |
244 | f = (f-1./f)/(f+1./f); | |
245 | goto ret; | |
246 | } | |
247 | ||
248 | sq: | |
249 | if(f < 0.) | |
250 | error("sqrt D"); | |
251 | f = sqrt(f); | |
252 | ||
253 | ret: | |
254 | d1 = f; | |
255 | return(d1); | |
256 | } | |
257 | ||
258 | data | |
259 | ex_comb(d1, d2) | |
260 | data d1, d2; | |
261 | { | |
262 | double f; | |
263 | ||
264 | if(d1 > d2) | |
265 | return(zero); | |
266 | f = gamma(d2+1.) - gamma(d1+1.) - gamma(d2-d1+1.); | |
267 | if(f > maxexp) | |
268 | error("comb D"); | |
269 | d1 = exp(f); | |
270 | return(d1); | |
271 | } | |
272 | ||
273 | data | |
274 | ex_and(d1, d2) | |
275 | data d1, d2; | |
276 | { | |
277 | ||
278 | if(d1!=zero && d2!=zero) | |
279 | return(one); | |
280 | return(zero); | |
281 | } | |
282 | ||
283 | data | |
284 | ex_or(d1, d2) | |
285 | data d1, d2; | |
286 | { | |
287 | ||
288 | if(d1!=zero || d2!=zero) | |
289 | return(one); | |
290 | return(zero); | |
291 | } | |
292 | ||
293 | data | |
294 | ex_nand(d1, d2) | |
295 | data d1, d2; | |
296 | { | |
297 | ||
298 | if(d1!=zero && d2!=zero) | |
299 | return(zero); | |
300 | return(one); | |
301 | } | |
302 | ||
303 | data | |
304 | ex_nor(d1, d2) | |
305 | data d1, d2; | |
306 | { | |
307 | ||
308 | if(d1!=zero || d2!=zero) | |
309 | return(zero); | |
310 | return(one); | |
311 | } | |
312 | ||
313 | data | |
314 | ex_lt(d1, d2) | |
315 | data d1, d2; | |
316 | { | |
317 | ||
318 | if(fuzz(d1, d2) < 0) | |
319 | return(one); | |
320 | return(zero); | |
321 | } | |
322 | ||
323 | data | |
324 | ex_le(d1, d2) | |
325 | data d1, d2; | |
326 | { | |
327 | ||
328 | if(fuzz(d1, d2) <= 0) | |
329 | return(one); | |
330 | return(zero); | |
331 | } | |
332 | ||
333 | data | |
334 | ex_eq(d1, d2) | |
335 | data d1, d2; | |
336 | { | |
337 | ||
338 | if(fuzz(d1, d2) == 0) | |
339 | return(one); | |
340 | return(zero); | |
341 | } | |
342 | ||
343 | data | |
344 | ex_ge(d1, d2) | |
345 | data d1, d2; | |
346 | { | |
347 | ||
348 | if(fuzz(d1, d2) >= 0) | |
349 | return(one); | |
350 | return(zero); | |
351 | } | |
352 | ||
353 | data | |
354 | ex_gt(d1, d2) | |
355 | data d1, d2; | |
356 | { | |
357 | ||
358 | if(fuzz(d1, d2) > 0) | |
359 | return(one); | |
360 | return(zero); | |
361 | } | |
362 | ||
363 | data | |
364 | ex_ne(d1, d2) | |
365 | data d1, d2; | |
366 | { | |
367 | ||
368 | if(fuzz(d1, d2) != 0) | |
369 | return(one); | |
370 | return(zero); | |
371 | } | |
372 | ||
373 | data | |
374 | ex_plus(d) | |
375 | data d; | |
376 | { | |
377 | ||
378 | return(d); | |
379 | } | |
380 | ||
381 | data | |
382 | ex_minus(d) | |
383 | data d; | |
384 | { | |
385 | ||
386 | return(-d); | |
387 | } | |
388 | ||
389 | data | |
390 | ex_sgn(d) | |
391 | data d; | |
392 | { | |
393 | ||
394 | if(d == zero) | |
395 | return(zero); | |
396 | if(d < zero) | |
397 | return(-one); | |
398 | return(one); | |
399 | } | |
400 | ||
401 | data | |
402 | ex_recip(d) | |
403 | data d; | |
404 | { | |
405 | ||
406 | if(d == zero) | |
407 | error("recip D"); | |
408 | return(one/d); | |
409 | } | |
410 | ||
411 | data | |
412 | ex_abs(d) | |
413 | data d; | |
414 | { | |
415 | ||
416 | if(d < zero) | |
417 | return(-d); | |
418 | return(d); | |
419 | } | |
420 | ||
421 | data | |
422 | ex_floor(d) | |
423 | data d; | |
424 | { | |
425 | ||
426 | d = floor(d + thread.fuzz); | |
427 | return(d); | |
428 | } | |
429 | ||
430 | data | |
431 | ex_ceil(d) | |
432 | data d; | |
433 | { | |
434 | ||
435 | d = ceil(d - thread.fuzz); | |
436 | return(d); | |
437 | } | |
438 | ||
439 | data | |
440 | ex_exp(d) | |
441 | data d; | |
442 | { | |
443 | double f; | |
444 | ||
445 | f = d; | |
446 | if(f > maxexp) | |
447 | error ("exp D"); | |
448 | d = exp(f); | |
449 | return(d); | |
450 | } | |
451 | ||
452 | data | |
453 | ex_loge(d) | |
454 | data d; | |
455 | { | |
456 | double f; | |
457 | ||
458 | f = d; | |
459 | if(f <= 0.) | |
460 | error("log D"); | |
461 | d = log(f); | |
462 | return(d); | |
463 | } | |
464 | ||
465 | data | |
466 | ex_pi(d) | |
467 | data d; | |
468 | { | |
469 | ||
470 | d = pi * d; | |
471 | return(d); | |
472 | } | |
473 | ||
474 | data | |
475 | ex_rand(d) | |
476 | data d; | |
477 | { | |
478 | double f; | |
479 | ||
480 | f = (rand()/(32768.*32768.*2.)) * d; | |
481 | d = floor(f) + thread.iorg; | |
482 | return(d); | |
483 | } | |
484 | ||
485 | data | |
486 | ex_fac(d) | |
487 | data d; | |
488 | { | |
489 | double f; | |
490 | ||
491 | f = gamma(d+1.); | |
492 | if(f > maxexp) | |
493 | error("fac D"); | |
494 | d = exp(f); | |
495 | if(signgam == -1) | |
496 | d = -d; | |
497 | return(d); | |
498 | } | |
499 | ||
500 | data | |
501 | ex_not(d) | |
502 | data d; | |
503 | { | |
504 | ||
505 | if(d == zero) | |
506 | return(one); | |
507 | return(zero); | |
508 | } |