Annotation of OpenXM_contrib2/asir2000/builtin/math.c, Revision 1.11
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.11 ! fujimoto 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/math.c,v 1.10 2011/12/15 16:53:26 ohara Exp $
1.2 noro 49: */
1.1 noro 50: #include "ca.h"
51: #include <math.h>
52: #include "parse.h"
1.11 ! fujimoto 53: #if defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__)
1.10 ohara 54: #include <float.h>
55: #endif
1.1 noro 56:
57: void Pdsqrt(),Pdsin(),Pdcos(),Pdtan(),Pdasin(),Pdacos(),Pdatan(),Pdlog(),Pdexp();
1.9 noro 58: void Pabs(),Pdfloor(),Pdceil(),Pdrint(),Pdisnan();
1.1 noro 59:
60: struct ftab math_tab[] = {
61: {"dsqrt",Pdsqrt,1},
62: {"dabs",Pabs,1},
63: {"dsin",Pdsin,1},
64: {"dcos",Pdcos,1},
65: {"dtan",Pdtan,1},
66: {"dlog",Pdlog,1},
67: {"dexp",Pdexp,1},
68: {"dasin",Pdasin,1},
69: {"dacos",Pdacos,1},
70: {"datan",Pdatan,1},
1.5 noro 71: {"floor",Pdfloor,1},
1.1 noro 72: {"dfloor",Pdfloor,1},
1.5 noro 73: {"ceil",Pdceil,1},
1.1 noro 74: {"dceil",Pdceil,1},
1.5 noro 75: {"rint",Pdrint,1},
1.1 noro 76: {"drint",Pdrint,1},
1.9 noro 77: {"disnan",Pdisnan,1},
1.1 noro 78: {0,0,0},
79: };
80:
1.8 noro 81: void get_ri(Num z,double *r,double *i)
82: {
83: if ( !z ) {
84: *r = 0; *i = 0; return;
85: }
86: if ( OID(z) != O_N )
87: error("get_ri : invalid argument");
88: switch ( NID(z) ) {
89: case N_Q: case N_R: case N_B:
90: *r = ToReal(z); *i = 0;
91: break;
92: case N_C:
93: *r = ToReal(((C)z)->r);
94: *i = ToReal(((C)z)->i);
95: break;
96: default:
97: error("get_ri : invalid argument");
98: break;
99: }
100: }
101:
1.1 noro 102: void Pabs(arg,rp)
103: NODE arg;
104: Real *rp;
105: {
1.8 noro 106: double s,r,i;
1.1 noro 107:
1.8 noro 108: if ( !ARG0(arg) ) {
109: *rp = 0; return;
110: }
111: get_ri((Num)ARG0(arg),&r,&i);
112: if ( i == 0 )
113: s = fabs(r);
114: else if ( r == 0 )
115: s = fabs(i);
116: else
117: s = sqrt(r*r+i*i);
1.1 noro 118: MKReal(s,*rp);
119: }
120:
121: void Pdsqrt(arg,rp)
122: NODE arg;
1.8 noro 123: Num *rp;
1.1 noro 124: {
1.8 noro 125: double s,r,i,a;
126: C z;
127: Real real;
128:
129: if ( !ARG0(arg) ) {
130: *rp = 0; return;
131: }
132: get_ri((Num)ARG0(arg),&r,&i);
133: if ( i == 0 )
134: if ( r > 0 ) {
135: s = sqrt(r);
136: MKReal(s,real);
137: *rp = (Num)real;
138: } else {
139: NEWC(z);
140: z->r = 0;
141: s = sqrt(-r); MKReal(s,real); z->i = (Num)real;
142: *rp = (Num)z;
143: }
144: else {
145: a = sqrt(r*r+i*i);
146: NEWC(z);
147: s = sqrt((r+a)/2); MKReal(s,real); z->r = (Num)real;
148: s = i>0?sqrt((-r+a)/2):-sqrt((-r+a)/2);
149: MKReal(s,real); z->i = (Num)real;
150: *rp = (Num)z;
151: }
1.1 noro 152: }
153:
154: void Pdsin(arg,rp)
155: NODE arg;
156: Real *rp;
157: {
158: double s;
159:
160: s = sin(ToReal(ARG0(arg)));
161: MKReal(s,*rp);
162: }
163:
164: void Pdcos(arg,rp)
165: NODE arg;
166: Real *rp;
167: {
168: double s;
169:
170: s = cos(ToReal(ARG0(arg)));
171: MKReal(s,*rp);
172: }
173:
174: void Pdtan(arg,rp)
175: NODE arg;
176: Real *rp;
177: {
178: double s;
179:
180: s = tan(ToReal(ARG0(arg)));
181: MKReal(s,*rp);
182: }
183:
184: void Pdasin(arg,rp)
185: NODE arg;
186: Real *rp;
187: {
188: double s;
189:
190: s = asin(ToReal(ARG0(arg)));
191: MKReal(s,*rp);
192: }
193:
194: void Pdacos(arg,rp)
195: NODE arg;
196: Real *rp;
197: {
198: double s;
199:
200: s = acos(ToReal(ARG0(arg)));
201: MKReal(s,*rp);
202: }
203:
204: void Pdatan(arg,rp)
205: NODE arg;
206: Real *rp;
207: {
208: double s;
209:
210: s = atan(ToReal(ARG0(arg)));
211: MKReal(s,*rp);
212: }
213:
214: void Pdlog(arg,rp)
215: NODE arg;
216: Real *rp;
217: {
218: double s;
219:
220: s = log(ToReal(ARG0(arg)));
221: MKReal(s,*rp);
222: }
223:
224: void Pdexp(arg,rp)
225: NODE arg;
226: Real *rp;
227: {
228: double s;
229:
230: s = exp(ToReal(ARG0(arg)));
231: MKReal(s,*rp);
232: }
233:
234: void Pdfloor(arg,rp)
235: NODE arg;
236: Q *rp;
237: {
238: L a;
239: unsigned int au,al;
240: int sgn;
241: Q q;
242: double d;
243:
244: if ( !ARG0(arg) ) {
245: *rp = 0;
246: return;
247: }
248: d = floor(ToReal(ARG0(arg)));
249: if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
250: error("dfloor : OverFlow");
1.4 noro 251: if ( !d ) {
252: *rp = 0;
253: return;
254: }
1.1 noro 255: a = (L)d;
256: if ( a < 0 ) {
257: sgn = -1;
258: a = -a;
259: } else
260: sgn = 1;
1.11 ! fujimoto 261: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__) || defined(__x86_64)
1.1 noro 262: au = ((unsigned int *)&a)[1];
263: al = ((unsigned int *)&a)[0];
264: #else
265: al = ((unsigned int *)&a)[1];
266: au = ((unsigned int *)&a)[0];
267: #endif
268: if ( au ) {
269: NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
270: PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
271: } else {
272: UTOQ(al,q); SGN(q) = sgn;
273: }
274: *rp = q;
275: }
276:
277: void Pdceil(arg,rp)
278: NODE arg;
279: Q *rp;
280: {
281: L a;
282: unsigned int au,al;
283: int sgn;
284: Q q;
285: double d;
286:
287: if ( !ARG0(arg) ) {
288: *rp = 0;
289: return;
290: }
291: d = ceil(ToReal(ARG0(arg)));
292: if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
293: error("dceil : OverFlow");
1.4 noro 294: if ( !d ) {
295: *rp = 0;
296: return;
297: }
1.1 noro 298: a = (L)d;
299: if ( a < 0 ) {
300: sgn = -1;
301: a = -a;
302: } else
303: sgn = 1;
1.11 ! fujimoto 304: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__) || defined(__x86_64)
1.1 noro 305: au = ((unsigned int *)&a)[1];
306: al = ((unsigned int *)&a)[0];
307: #else
308: al = ((unsigned int *)&a)[1];
309: au = ((unsigned int *)&a)[0];
310: #endif
311: if ( au ) {
312: NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
313: PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
314: } else {
315: UTOQ(al,q); SGN(q) = sgn;
316: }
317: *rp = q;
318: }
319:
320: void Pdrint(arg,rp)
321: NODE arg;
322: Q *rp;
323: {
324: L a;
325: unsigned int au,al;
326: int sgn;
327: Q q;
328: double d;
329:
330: if ( !ARG0(arg) ) {
331: *rp = 0;
332: return;
333: }
1.11 ! fujimoto 334: #if defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__)
1.1 noro 335: d = ToReal(ARG0(arg));
336: if ( d > 0 )
337: d = floor(d+0.5);
338: else
339: d = ceil(d-0.5);
340: #else
341: d = rint(ToReal(ARG0(arg)));
342: #endif
343: if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
344: error("drint : OverFlow");
345: a = (L)d;
346: if ( a < 0 ) {
347: sgn = -1;
348: a = -a;
349: } else
350: sgn = 1;
1.11 ! fujimoto 351: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__) || defined(__x86_64)
1.1 noro 352: au = ((unsigned int *)&a)[1];
353: al = ((unsigned int *)&a)[0];
354: #else
355: al = ((unsigned int *)&a)[1];
356: au = ((unsigned int *)&a)[0];
357: #endif
358: if ( au ) {
359: NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
360: PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
1.6 noro 361: } else if ( al ) {
1.1 noro 362: UTOQ(al,q); SGN(q) = sgn;
1.6 noro 363: } else
364: q = 0;
1.1 noro 365: *rp = q;
366: }
1.9 noro 367:
368: void Pdisnan(NODE arg,Q *rp)
369: {
370: Real r;
371: double d;
1.11 ! fujimoto 372: #if defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__)
1.10 ohara 373: int c;
374: #endif
1.9 noro 375:
376: r = (Real)ARG0(arg);
377: if ( !r || !NUM(r) || !REAL(r) ) {
378: *rp = 0;
379: return;
380: }
381: d = ToReal(r);
1.11 ! fujimoto 382: #if defined(VISUAL) || defined(__MINGW32__) || defined(__MINGW64__)
1.10 ohara 383: c = _fpclass(d);
384: if ( c == _FPCLASS_SNAN || c == _FPCLASS_QNAN ) *rp = ONE;
385: else if ( c == _FPCLASS_PINF || c == _FPCLASS_NINF ) STOQ(2,*rp);
386: #else
1.9 noro 387: if ( isnan(d) ) *rp = ONE;
388: else if ( isinf(d) ) STOQ(2,*rp);
1.10 ohara 389: #endif
1.9 noro 390: else *rp = 0;
391: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to math.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / builtin