Annotation of OpenXM_contrib2/asir2000/builtin/bfaux.c, Revision 1.6
1.6 ! takayama 1: /* $OpenXM: OpenXM_contrib2/asir2000/builtin/bfaux.c,v 1.5 2015/08/07 05:30:35 noro Exp $ */
1.1 noro 2: #include "ca.h"
3: #include "parse.h"
4:
1.4 noro 5: void Peval(), Psetprec(), Psetbprec(), Ptodouble(), Psetround();
1.6 ! takayama 6: void Pmpfr_gamma(), Pmpfr_floor();
1.1 noro 7:
8: struct ftab bf_tab[] = {
9: {"eval",Peval,-2},
10: {"setprec",Psetprec,-1},
1.3 noro 11: {"setbprec",Psetbprec,-1},
1.4 noro 12: {"setround",Psetround,-1},
1.1 noro 13: {"todouble",Ptodouble,1},
1.5 noro 14: {"mpfr_gamma",Pmpfr_gamma,-2},
1.6 ! takayama 15: {"mpfr_floor",Pmpfr_floor,-1},
1.1 noro 16: {0,0,0},
17: };
18:
1.4 noro 19: int mpfr_roundmode = MPFR_RNDN;
20:
1.1 noro 21: void Ptodouble(NODE arg,Num *rp)
22: {
23: double r,i;
24: Real real,imag;
25: Num num;
26:
27: asir_assert(ARG0(arg),O_N,"todouble");
28: num = (Num)ARG0(arg);
29: if ( !num ) {
30: *rp = 0;
31: return;
32: }
33: switch ( NID(num) ) {
34: case N_R: case N_Q: case N_B:
35: r = ToReal(num);
36: MKReal(r,real);
37: *rp = (Num)real;
38: break;
39: case N_C:
40: r = ToReal(((C)num)->r);
41: i = ToReal(((C)num)->i);
42: MKReal(r,real);
43: MKReal(i,imag);
44: reimtocplx((Num)real,(Num)imag,rp);
45: break;
46: default:
47: *rp = num;
48: break;
49: }
50: }
51:
1.4 noro 52: void Peval(NODE arg,Obj *rp)
1.1 noro 53: {
54: int prec;
55:
56: asir_assert(ARG0(arg),O_R,"eval");
57: if ( argc(arg) == 2 ) {
1.3 noro 58: prec = QTOS((Q)ARG1(arg))*3.32193;
1.1 noro 59: if ( prec < MPFR_PREC_MIN ) prec = MPFR_PREC_MIN;
60: else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX;
61: } else
62: prec = 0;
1.2 noro 63: evalr(CO,(Obj)ARG0(arg),prec,rp);
1.1 noro 64: }
65:
1.3 noro 66: /* set/get decimal precision */
1.1 noro 67:
68: void Psetprec(NODE arg,Obj *rp)
69: {
70: int p;
71: Q q;
1.3 noro 72: int prec,dprec;
73:
74: prec = mpfr_get_default_prec();
75: /* decimal precision */
76: dprec = prec*0.30103;
77: STOQ(dprec,q); *rp = (Obj)q;
78: if ( arg ) {
79: asir_assert(ARG0(arg),O_N,"setprec");
80: prec = QTOS((Q)ARG0(arg))*3.32193;
81: if ( p > 0 )
82: prec = p;
83: }
84: if ( prec < MPFR_PREC_MIN ) prec = MPFR_PREC_MIN;
85: else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX;
86: mpfr_set_default_prec(prec);
87: }
1.1 noro 88:
1.3 noro 89: /* set/get bit precision */
1.1 noro 90:
1.3 noro 91: void Psetbprec(NODE arg,Obj *rp)
92: {
93: int p;
94: Q q;
95: int prec;
96:
97: prec = mpfr_get_default_prec();
98: STOQ(prec,q); *rp = (Obj)q;
1.1 noro 99: if ( arg ) {
1.3 noro 100: asir_assert(ARG0(arg),O_N,"setbprec");
101: prec = QTOS((Q)ARG0(arg));
1.1 noro 102: if ( p > 0 )
103: prec = p;
104: }
105: if ( prec < MPFR_PREC_MIN ) prec = MPFR_PREC_MIN;
106: else if ( prec > MPFR_PREC_MAX ) prec = MPFR_PREC_MAX;
107: mpfr_set_default_prec(prec);
108: }
109:
1.4 noro 110: void Psetround(NODE arg,Q *rp)
111: {
112: int round;
113:
114: STOQ(mpfr_roundmode,*rp);
115: if ( arg ) {
116: asir_assert(ARG0(arg),O_N,"setround");
117: round = QTOS((Q)ARG0(arg));
118: switch ( round ) {
119: case 0:
120: mpfr_roundmode = MPFR_RNDN;
121: break;
122: case 1:
123: mpfr_roundmode = MPFR_RNDZ;
124: break;
125: case 2:
126: mpfr_roundmode = MPFR_RNDU;
127: break;
128: case 3:
129: mpfr_roundmode = MPFR_RNDD;
130: break;
131: case 4:
132: mpfr_roundmode = MPFR_RNDA;
133: break;
134: case 5:
135: mpfr_roundmode = MPFR_RNDF;
136: break;
137: case 6:
138: mpfr_roundmode = MPFR_RNDNA;
139: break;
140: default:
141: error("setround : invalid rounding mode");
142: break;
143: }
144: }
145: }
146:
1.1 noro 147: Num tobf(Num a,int prec);
148:
149: void mp_pi(NODE arg,BF *rp)
150: {
151: int prec;
152: BF r;
153:
154: prec = arg ? QTOS((Q)ARG0(arg)) : 0;
155: NEWBF(r);
156: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 157: mpfr_const_pi(r->body,mpfr_roundmode);
1.1 noro 158: *rp = r;
159: }
160:
161: void mp_e(NODE arg,BF *rp)
162: {
163: int prec;
164: mpfr_t one;
165: BF r;
166:
167: prec = arg ? QTOS((Q)ARG0(arg)) : 0;
168: NEWBF(r);
169: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
170: mpfr_init(one);
1.4 noro 171: mpfr_set_ui(one,1,mpfr_roundmode);
172: mpfr_exp(r->body,one,mpfr_roundmode);
1.1 noro 173: *rp = r;
174: }
175:
176: void mp_sin(NODE arg,BF *rp)
177: {
178: Num a;
179: int prec;
180: BF r;
181:
182: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
183: a = tobf(ARG0(arg),prec);
184: NEWBF(r);
185: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 186: mpfr_sin(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 187: *rp = r;
188: }
189:
190: void mp_cos(NODE arg,BF *rp)
191: {
192: Num a;
193: int prec;
194: BF r;
195:
196: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
197: a = tobf(ARG0(arg),prec);
198: NEWBF(r);
199: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 200: mpfr_cos(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 201: *rp = r;
202: }
203:
204: void mp_tan(NODE arg,BF *rp)
205: {
206: Num a;
207: int prec;
208: BF r;
209:
210: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
211: a = tobf(ARG0(arg),prec);
212: NEWBF(r);
213: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 214: mpfr_tan(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 215: *rp = r;
216: }
217:
218: void mp_asin(NODE arg,BF *rp)
219: {
220: Num a;
221: int prec;
222: BF r;
223:
224: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
225: a = tobf(ARG0(arg),prec);
226: NEWBF(r);
227: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 228: mpfr_asin(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 229: *rp = r;
230: }
231: void mp_acos(NODE arg,BF *rp)
232: {
233: Num a;
234: int prec;
235: BF r;
236:
237: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
238: a = tobf(ARG0(arg),prec);
239: NEWBF(r);
240: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 241: mpfr_acos(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 242: *rp = r;
243: }
244: void mp_atan(NODE arg,BF *rp)
245: {
246: Num a;
247: int prec;
248: BF r;
249:
250: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
251: a = tobf(ARG0(arg),prec);
252: NEWBF(r);
253: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 254: mpfr_atan(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 255: *rp = r;
256: }
257:
258: void mp_sinh(NODE arg,BF *rp)
259: {
260: Num a;
261: int prec;
262: BF r;
263:
264: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
265: a = tobf(ARG0(arg),prec);
266: NEWBF(r);
267: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 268: mpfr_sinh(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 269: *rp = r;
270: }
271:
272: void mp_cosh(NODE arg,BF *rp)
273: {
274: Num a;
275: int prec;
276: BF r;
277:
278: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
279: a = tobf(ARG0(arg),prec);
280: NEWBF(r);
281: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 282: mpfr_cosh(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 283: *rp = r;
284: }
285:
286: void mp_tanh(NODE arg,BF *rp)
287: {
288: Num a;
289: int prec;
290: BF r;
291:
292: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
293: a = tobf(ARG0(arg),prec);
294: NEWBF(r);
295: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 296: mpfr_tanh(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 297: *rp = r;
298: }
299:
300: void mp_asinh(NODE arg,BF *rp)
301: {
302: Num a;
303: int prec;
304: BF r;
305:
306: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
307: a = tobf(ARG0(arg),prec);
308: NEWBF(r);
309: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 310: mpfr_asinh(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 311: *rp = r;
312: }
313: void mp_acosh(NODE arg,BF *rp)
314: {
315: Num a;
316: int prec;
317: BF r;
318:
319: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
320: a = tobf(ARG0(arg),prec);
321: NEWBF(r);
322: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 323: mpfr_acosh(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 324: *rp = r;
325: }
326: void mp_atanh(NODE arg,BF *rp)
327: {
328: Num a;
329: int prec;
330: BF r;
331:
332: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
333: a = tobf(ARG0(arg),prec);
334: NEWBF(r);
335: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 336: mpfr_atanh(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 337: *rp = r;
338: }
339:
340: void mp_exp(NODE arg,BF *rp)
341: {
342: Num a;
343: int prec;
344: BF r;
345:
346: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
347: a = tobf(ARG0(arg),prec);
348: NEWBF(r);
349: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 350: mpfr_exp(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 351: *rp = r;
352: }
353:
354: void mp_log(NODE arg,BF *rp)
355: {
356: Num a;
357: int prec;
358: BF r;
359:
360: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
361: a = tobf(ARG0(arg),prec);
362: NEWBF(r);
363: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 364: mpfr_log(r->body,((BF)a)->body,mpfr_roundmode);
1.1 noro 365: *rp = r;
366: }
367:
368: void mp_pow(NODE arg,BF *rp)
369: {
370: Num a,e;
371: int prec;
372: BF r;
373:
374: prec = NEXT(NEXT(arg)) ? QTOS((Q)ARG2(arg)) : 0;
375: a = tobf(ARG0(arg),prec);
376: e = tobf(ARG1(arg),prec);
377: NEWBF(r);
378: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
1.4 noro 379: mpfr_pow(r->body,((BF)a)->body,((BF)e)->body,mpfr_roundmode);
1.1 noro 380: *rp = r;
381: }
1.5 noro 382:
383: void Pmpfr_gamma(NODE arg,BF *rp)
384: {
385: Num a;
386: int prec;
387: BF r;
388:
389: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
390: prec *= 3.32193;
391: a = tobf(ARG0(arg),prec);
392: NEWBF(r);
393: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
394: mpfr_gamma(r->body,((BF)a)->body,mpfr_roundmode);
395: *rp = r;
396: }
1.6 ! takayama 397:
! 398: void Pmpfr_floor(NODE arg,BF *rp)
! 399: {
! 400: Num a;
! 401: int prec;
! 402: BF r;
! 403:
! 404: prec = NEXT(arg) ? QTOS((Q)ARG1(arg)) : 0;
! 405: prec *= 3.32193;
! 406: a = tobf(ARG0(arg),prec);
! 407: NEWBF(r);
! 408: prec ? mpfr_init2(r->body,prec) : mpfr_init(r->body);
! 409: mpfr_floor(r->body,((BF)a)->body);
! 410: *rp = r;
! 411: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>