Annotation of OpenXM_contrib2/asir2000/builtin/math.c, Revision 1.9
1.2 noro 1: /*
2: * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
3: * All rights reserved.
4: *
5: * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
6: * non-exclusive and royalty-free license to use, copy, modify and
7: * redistribute, solely for non-commercial and non-profit purposes, the
8: * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
9: * conditions of this Agreement. For the avoidance of doubt, you acquire
10: * only a limited right to use the SOFTWARE hereunder, and FLL or any
11: * third party developer retains all rights, including but not limited to
12: * copyrights, in and to the SOFTWARE.
13: *
14: * (1) FLL does not grant you a license in any way for commercial
15: * purposes. You may use the SOFTWARE only for non-commercial and
16: * non-profit purposes only, such as academic, research and internal
17: * business use.
18: * (2) The SOFTWARE is protected by the Copyright Law of Japan and
19: * international copyright treaties. If you make copies of the SOFTWARE,
20: * with or without modification, as permitted hereunder, you shall affix
21: * to all such copies of the SOFTWARE the above copyright notice.
22: * (3) An explicit reference to this SOFTWARE and its copyright owner
23: * shall be made on your publication or presentation in any form of the
24: * results obtained by use of the SOFTWARE.
25: * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
1.3 noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.2 noro 27: * for such modification or the source code of the modified part of the
28: * SOFTWARE.
29: *
30: * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
31: * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
32: * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
33: * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
34: * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
35: * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
36: * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
37: * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
38: * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
39: * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
40: * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
41: * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
42: * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
43: * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
44: * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
45: * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
46: * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
47: *
1.9 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/math.c,v 1.8 2003/12/26 05:47:37 noro Exp $
1.2 noro 49: */
1.1 noro 50: #include "ca.h"
51: #include <math.h>
52: #include "parse.h"
53:
54: void Pdsqrt(),Pdsin(),Pdcos(),Pdtan(),Pdasin(),Pdacos(),Pdatan(),Pdlog(),Pdexp();
1.9 ! noro 55: void Pabs(),Pdfloor(),Pdceil(),Pdrint(),Pdisnan();
1.1 noro 56:
57: struct ftab math_tab[] = {
58: {"dsqrt",Pdsqrt,1},
59: {"dabs",Pabs,1},
60: {"dsin",Pdsin,1},
61: {"dcos",Pdcos,1},
62: {"dtan",Pdtan,1},
63: {"dlog",Pdlog,1},
64: {"dexp",Pdexp,1},
65: {"dasin",Pdasin,1},
66: {"dacos",Pdacos,1},
67: {"datan",Pdatan,1},
1.5 noro 68: {"floor",Pdfloor,1},
1.1 noro 69: {"dfloor",Pdfloor,1},
1.5 noro 70: {"ceil",Pdceil,1},
1.1 noro 71: {"dceil",Pdceil,1},
1.5 noro 72: {"rint",Pdrint,1},
1.1 noro 73: {"drint",Pdrint,1},
1.9 ! noro 74: {"disnan",Pdisnan,1},
1.1 noro 75: {0,0,0},
76: };
77:
1.8 noro 78: void get_ri(Num z,double *r,double *i)
79: {
80: if ( !z ) {
81: *r = 0; *i = 0; return;
82: }
83: if ( OID(z) != O_N )
84: error("get_ri : invalid argument");
85: switch ( NID(z) ) {
86: case N_Q: case N_R: case N_B:
87: *r = ToReal(z); *i = 0;
88: break;
89: case N_C:
90: *r = ToReal(((C)z)->r);
91: *i = ToReal(((C)z)->i);
92: break;
93: default:
94: error("get_ri : invalid argument");
95: break;
96: }
97: }
98:
1.1 noro 99: void Pabs(arg,rp)
100: NODE arg;
101: Real *rp;
102: {
1.8 noro 103: double s,r,i;
1.1 noro 104:
1.8 noro 105: if ( !ARG0(arg) ) {
106: *rp = 0; return;
107: }
108: get_ri((Num)ARG0(arg),&r,&i);
109: if ( i == 0 )
110: s = fabs(r);
111: else if ( r == 0 )
112: s = fabs(i);
113: else
114: s = sqrt(r*r+i*i);
1.1 noro 115: MKReal(s,*rp);
116: }
117:
118: void Pdsqrt(arg,rp)
119: NODE arg;
1.8 noro 120: Num *rp;
1.1 noro 121: {
1.8 noro 122: double s,r,i,a;
123: C z;
124: Real real;
125:
126: if ( !ARG0(arg) ) {
127: *rp = 0; return;
128: }
129: get_ri((Num)ARG0(arg),&r,&i);
130: if ( i == 0 )
131: if ( r > 0 ) {
132: s = sqrt(r);
133: MKReal(s,real);
134: *rp = (Num)real;
135: } else {
136: NEWC(z);
137: z->r = 0;
138: s = sqrt(-r); MKReal(s,real); z->i = (Num)real;
139: *rp = (Num)z;
140: }
141: else {
142: a = sqrt(r*r+i*i);
143: NEWC(z);
144: s = sqrt((r+a)/2); MKReal(s,real); z->r = (Num)real;
145: s = i>0?sqrt((-r+a)/2):-sqrt((-r+a)/2);
146: MKReal(s,real); z->i = (Num)real;
147: *rp = (Num)z;
148: }
1.1 noro 149: }
150:
151: void Pdsin(arg,rp)
152: NODE arg;
153: Real *rp;
154: {
155: double s;
156:
157: s = sin(ToReal(ARG0(arg)));
158: MKReal(s,*rp);
159: }
160:
161: void Pdcos(arg,rp)
162: NODE arg;
163: Real *rp;
164: {
165: double s;
166:
167: s = cos(ToReal(ARG0(arg)));
168: MKReal(s,*rp);
169: }
170:
171: void Pdtan(arg,rp)
172: NODE arg;
173: Real *rp;
174: {
175: double s;
176:
177: s = tan(ToReal(ARG0(arg)));
178: MKReal(s,*rp);
179: }
180:
181: void Pdasin(arg,rp)
182: NODE arg;
183: Real *rp;
184: {
185: double s;
186:
187: s = asin(ToReal(ARG0(arg)));
188: MKReal(s,*rp);
189: }
190:
191: void Pdacos(arg,rp)
192: NODE arg;
193: Real *rp;
194: {
195: double s;
196:
197: s = acos(ToReal(ARG0(arg)));
198: MKReal(s,*rp);
199: }
200:
201: void Pdatan(arg,rp)
202: NODE arg;
203: Real *rp;
204: {
205: double s;
206:
207: s = atan(ToReal(ARG0(arg)));
208: MKReal(s,*rp);
209: }
210:
211: void Pdlog(arg,rp)
212: NODE arg;
213: Real *rp;
214: {
215: double s;
216:
217: s = log(ToReal(ARG0(arg)));
218: MKReal(s,*rp);
219: }
220:
221: void Pdexp(arg,rp)
222: NODE arg;
223: Real *rp;
224: {
225: double s;
226:
227: s = exp(ToReal(ARG0(arg)));
228: MKReal(s,*rp);
229: }
230:
231: void Pdfloor(arg,rp)
232: NODE arg;
233: Q *rp;
234: {
235: L a;
236: unsigned int au,al;
237: int sgn;
238: Q q;
239: double d;
240:
241: if ( !ARG0(arg) ) {
242: *rp = 0;
243: return;
244: }
245: d = floor(ToReal(ARG0(arg)));
246: if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
247: error("dfloor : OverFlow");
1.4 noro 248: if ( !d ) {
249: *rp = 0;
250: return;
251: }
1.1 noro 252: a = (L)d;
253: if ( a < 0 ) {
254: sgn = -1;
255: a = -a;
256: } else
257: sgn = 1;
1.7 noro 258: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__x86_64)
1.1 noro 259: au = ((unsigned int *)&a)[1];
260: al = ((unsigned int *)&a)[0];
261: #else
262: al = ((unsigned int *)&a)[1];
263: au = ((unsigned int *)&a)[0];
264: #endif
265: if ( au ) {
266: NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
267: PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
268: } else {
269: UTOQ(al,q); SGN(q) = sgn;
270: }
271: *rp = q;
272: }
273:
274: void Pdceil(arg,rp)
275: NODE arg;
276: Q *rp;
277: {
278: L a;
279: unsigned int au,al;
280: int sgn;
281: Q q;
282: double d;
283:
284: if ( !ARG0(arg) ) {
285: *rp = 0;
286: return;
287: }
288: d = ceil(ToReal(ARG0(arg)));
289: if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
290: error("dceil : OverFlow");
1.4 noro 291: if ( !d ) {
292: *rp = 0;
293: return;
294: }
1.1 noro 295: a = (L)d;
296: if ( a < 0 ) {
297: sgn = -1;
298: a = -a;
299: } else
300: sgn = 1;
1.7 noro 301: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__x86_64)
1.1 noro 302: au = ((unsigned int *)&a)[1];
303: al = ((unsigned int *)&a)[0];
304: #else
305: al = ((unsigned int *)&a)[1];
306: au = ((unsigned int *)&a)[0];
307: #endif
308: if ( au ) {
309: NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
310: PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
311: } else {
312: UTOQ(al,q); SGN(q) = sgn;
313: }
314: *rp = q;
315: }
316:
317: void Pdrint(arg,rp)
318: NODE arg;
319: Q *rp;
320: {
321: L a;
322: unsigned int au,al;
323: int sgn;
324: Q q;
325: double d;
326:
327: if ( !ARG0(arg) ) {
328: *rp = 0;
329: return;
330: }
331: #if defined(VISUAL)
332: d = ToReal(ARG0(arg));
333: if ( d > 0 )
334: d = floor(d+0.5);
335: else
336: d = ceil(d-0.5);
337: #else
338: d = rint(ToReal(ARG0(arg)));
339: #endif
340: if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
341: error("drint : OverFlow");
342: a = (L)d;
343: if ( a < 0 ) {
344: sgn = -1;
345: a = -a;
346: } else
347: sgn = 1;
1.7 noro 348: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__x86_64)
1.1 noro 349: au = ((unsigned int *)&a)[1];
350: al = ((unsigned int *)&a)[0];
351: #else
352: al = ((unsigned int *)&a)[1];
353: au = ((unsigned int *)&a)[0];
354: #endif
355: if ( au ) {
356: NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
357: PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
1.6 noro 358: } else if ( al ) {
1.1 noro 359: UTOQ(al,q); SGN(q) = sgn;
1.6 noro 360: } else
361: q = 0;
1.1 noro 362: *rp = q;
363: }
1.9 ! noro 364:
! 365: void Pdisnan(NODE arg,Q *rp)
! 366: {
! 367: Real r;
! 368: double d;
! 369:
! 370: r = (Real)ARG0(arg);
! 371: if ( !r || !NUM(r) || !REAL(r) ) {
! 372: *rp = 0;
! 373: return;
! 374: }
! 375: d = ToReal(r);
! 376: if ( isnan(d) ) *rp = ONE;
! 377: else if ( isinf(d) ) STOQ(2,*rp);
! 378: else *rp = 0;
! 379: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>