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