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

Annotation of OpenXM_contrib2/asir2000/builtin/math.c, Revision 1.8

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.8     ! noro       48:  * $OpenXM: OpenXM_contrib2/asir2000/builtin/math.c,v 1.7 2003/12/24 08:00:38 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();
                     55: void Pabs(),Pdfloor(),Pdceil(),Pdrint();
                     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},
                     74:        {0,0,0},
                     75: };
                     76:
1.8     ! noro       77: void get_ri(Num z,double *r,double *i)
        !            78: {
        !            79:        if ( !z ) {
        !            80:                *r = 0; *i = 0; return;
        !            81:        }
        !            82:        if ( OID(z) != O_N )
        !            83:                error("get_ri : invalid argument");
        !            84:        switch ( NID(z) ) {
        !            85:                case N_Q: case N_R: case N_B:
        !            86:                        *r = ToReal(z); *i = 0;
        !            87:                        break;
        !            88:                case N_C:
        !            89:                        *r = ToReal(((C)z)->r);
        !            90:                        *i = ToReal(((C)z)->i);
        !            91:                        break;
        !            92:                default:
        !            93:                        error("get_ri : invalid argument");
        !            94:                        break;
        !            95:        }
        !            96: }
        !            97:
1.1       noro       98: void Pabs(arg,rp)
                     99: NODE arg;
                    100: Real *rp;
                    101: {
1.8     ! noro      102:        double s,r,i;
1.1       noro      103:
1.8     ! noro      104:        if ( !ARG0(arg) ) {
        !           105:                *rp = 0; return;
        !           106:        }
        !           107:        get_ri((Num)ARG0(arg),&r,&i);
        !           108:        if ( i == 0 )
        !           109:                s = fabs(r);
        !           110:        else if ( r == 0 )
        !           111:                s = fabs(i);
        !           112:        else
        !           113:                s = sqrt(r*r+i*i);
1.1       noro      114:        MKReal(s,*rp);
                    115: }
                    116:
                    117: void Pdsqrt(arg,rp)
                    118: NODE arg;
1.8     ! noro      119: Num *rp;
1.1       noro      120: {
1.8     ! noro      121:        double s,r,i,a;
        !           122:        C z;
        !           123:        Real real;
        !           124:
        !           125:        if ( !ARG0(arg) ) {
        !           126:                *rp = 0; return;
        !           127:        }
        !           128:        get_ri((Num)ARG0(arg),&r,&i);
        !           129:        if ( i == 0 )
        !           130:                if ( r > 0 ) {
        !           131:                        s = sqrt(r);
        !           132:                        MKReal(s,real);
        !           133:                        *rp = (Num)real;
        !           134:                } else {
        !           135:                        NEWC(z);
        !           136:                        z->r = 0;
        !           137:                        s = sqrt(-r); MKReal(s,real); z->i = (Num)real;
        !           138:                        *rp = (Num)z;
        !           139:                }
        !           140:        else {
        !           141:                a = sqrt(r*r+i*i);
        !           142:                NEWC(z);
        !           143:                s = sqrt((r+a)/2); MKReal(s,real); z->r = (Num)real;
        !           144:                s = i>0?sqrt((-r+a)/2):-sqrt((-r+a)/2);
        !           145:                MKReal(s,real); z->i = (Num)real;
        !           146:                *rp = (Num)z;
        !           147:        }
1.1       noro      148: }
                    149:
                    150: void Pdsin(arg,rp)
                    151: NODE arg;
                    152: Real *rp;
                    153: {
                    154:        double s;
                    155:
                    156:        s = sin(ToReal(ARG0(arg)));
                    157:        MKReal(s,*rp);
                    158: }
                    159:
                    160: void Pdcos(arg,rp)
                    161: NODE arg;
                    162: Real *rp;
                    163: {
                    164:        double s;
                    165:
                    166:        s = cos(ToReal(ARG0(arg)));
                    167:        MKReal(s,*rp);
                    168: }
                    169:
                    170: void Pdtan(arg,rp)
                    171: NODE arg;
                    172: Real *rp;
                    173: {
                    174:        double s;
                    175:
                    176:        s = tan(ToReal(ARG0(arg)));
                    177:        MKReal(s,*rp);
                    178: }
                    179:
                    180: void Pdasin(arg,rp)
                    181: NODE arg;
                    182: Real *rp;
                    183: {
                    184:        double s;
                    185:
                    186:        s = asin(ToReal(ARG0(arg)));
                    187:        MKReal(s,*rp);
                    188: }
                    189:
                    190: void Pdacos(arg,rp)
                    191: NODE arg;
                    192: Real *rp;
                    193: {
                    194:        double s;
                    195:
                    196:        s = acos(ToReal(ARG0(arg)));
                    197:        MKReal(s,*rp);
                    198: }
                    199:
                    200: void Pdatan(arg,rp)
                    201: NODE arg;
                    202: Real *rp;
                    203: {
                    204:        double s;
                    205:
                    206:        s = atan(ToReal(ARG0(arg)));
                    207:        MKReal(s,*rp);
                    208: }
                    209:
                    210: void Pdlog(arg,rp)
                    211: NODE arg;
                    212: Real *rp;
                    213: {
                    214:        double s;
                    215:
                    216:        s = log(ToReal(ARG0(arg)));
                    217:        MKReal(s,*rp);
                    218: }
                    219:
                    220: void Pdexp(arg,rp)
                    221: NODE arg;
                    222: Real *rp;
                    223: {
                    224:        double s;
                    225:
                    226:        s = exp(ToReal(ARG0(arg)));
                    227:        MKReal(s,*rp);
                    228: }
                    229:
                    230: void Pdfloor(arg,rp)
                    231: NODE arg;
                    232: Q *rp;
                    233: {
                    234:        L a;
                    235:        unsigned int au,al;
                    236:        int sgn;
                    237:        Q q;
                    238:        double d;
                    239:
                    240:        if ( !ARG0(arg) ) {
                    241:                *rp = 0;
                    242:                return;
                    243:        }
                    244:        d = floor(ToReal(ARG0(arg)));
                    245:        if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
                    246:                error("dfloor : OverFlow");
1.4       noro      247:        if ( !d ) {
                    248:                *rp = 0;
                    249:                return;
                    250:        }
1.1       noro      251:        a = (L)d;
                    252:        if ( a < 0 ) {
                    253:                sgn = -1;
                    254:                a = -a;
                    255:        } else
                    256:                sgn = 1;
1.7       noro      257: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__x86_64)
1.1       noro      258:        au = ((unsigned int *)&a)[1];
                    259:        al = ((unsigned int *)&a)[0];
                    260: #else
                    261:        al = ((unsigned int *)&a)[1];
                    262:        au = ((unsigned int *)&a)[0];
                    263: #endif
                    264:        if ( au ) {
                    265:                NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
                    266:                PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
                    267:        } else {
                    268:                UTOQ(al,q); SGN(q) = sgn;
                    269:        }
                    270:        *rp = q;
                    271: }
                    272:
                    273: void Pdceil(arg,rp)
                    274: NODE arg;
                    275: Q *rp;
                    276: {
                    277:        L a;
                    278:        unsigned int au,al;
                    279:        int sgn;
                    280:        Q q;
                    281:        double d;
                    282:
                    283:        if ( !ARG0(arg) ) {
                    284:                *rp = 0;
                    285:                return;
                    286:        }
                    287:        d = ceil(ToReal(ARG0(arg)));
                    288:        if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
                    289:                error("dceil : OverFlow");
1.4       noro      290:        if ( !d ) {
                    291:                *rp = 0;
                    292:                return;
                    293:        }
1.1       noro      294:        a = (L)d;
                    295:        if ( a < 0 ) {
                    296:                sgn = -1;
                    297:                a = -a;
                    298:        } else
                    299:                sgn = 1;
1.7       noro      300: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__x86_64)
1.1       noro      301:        au = ((unsigned int *)&a)[1];
                    302:        al = ((unsigned int *)&a)[0];
                    303: #else
                    304:        al = ((unsigned int *)&a)[1];
                    305:        au = ((unsigned int *)&a)[0];
                    306: #endif
                    307:        if ( au ) {
                    308:                NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
                    309:                PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
                    310:        } else {
                    311:                UTOQ(al,q); SGN(q) = sgn;
                    312:        }
                    313:        *rp = q;
                    314: }
                    315:
                    316: void Pdrint(arg,rp)
                    317: NODE arg;
                    318: Q *rp;
                    319: {
                    320:        L a;
                    321:        unsigned int au,al;
                    322:        int sgn;
                    323:        Q q;
                    324:        double d;
                    325:
                    326:        if ( !ARG0(arg) ) {
                    327:                *rp = 0;
                    328:                return;
                    329:        }
                    330: #if defined(VISUAL)
                    331:        d = ToReal(ARG0(arg));
                    332:        if ( d > 0 )
                    333:                d = floor(d+0.5);
                    334:        else
                    335:                d = ceil(d-0.5);
                    336: #else
                    337:        d = rint(ToReal(ARG0(arg)));
                    338: #endif
                    339:        if ( d < -9.223372036854775808e18 || d >= 9.223372036854775808e18 )
                    340:                error("drint : OverFlow");
                    341:        a = (L)d;
                    342:        if ( a < 0 ) {
                    343:                sgn = -1;
                    344:                a = -a;
                    345:        } else
                    346:                sgn = 1;
1.7       noro      347: #if defined(i386) || defined(__alpha) || defined(VISUAL) || defined(__x86_64)
1.1       noro      348:        au = ((unsigned int *)&a)[1];
                    349:        al = ((unsigned int *)&a)[0];
                    350: #else
                    351:        al = ((unsigned int *)&a)[1];
                    352:        au = ((unsigned int *)&a)[0];
                    353: #endif
                    354:        if ( au ) {
                    355:                NEWQ(q); SGN(q) = sgn; NM(q)=NALLOC(2); DN(q)=0;
                    356:                PL(NM(q))=2; BD(NM(q))[0]=al; BD(NM(q))[1] = au;
1.6       noro      357:        } else if ( al ) {
1.1       noro      358:                UTOQ(al,q); SGN(q) = sgn;
1.6       noro      359:        } else
                    360:                q = 0;
1.1       noro      361:        *rp = q;
                    362: }

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