[BACK]Return to bfaux.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / builtin

Annotation of OpenXM_contrib2/asir2000/builtin/bfaux.c, Revision 1.5

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

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>