![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to ao.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [CSRG BSD Unix] / 3BSD / cmd / apl|
Return to ao.c CVS log | Up to [CSRG BSD Unix] / 3BSD / cmd / apl |
1.1 root 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: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.