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

Annotation of OpenXM_contrib2/asir2000/engine/gfspn.c, Revision 1.2

1.2     ! noro        1: /* OpenXM */
        !             2:
1.1       noro        3: #include "ca.h"
                      4: #include "base.h"
                      5:
                      6: UM current_mod_gfspn;
                      7:
                      8: void setmod_gfspn(p)
                      9: UM p;
                     10: {
                     11:        current_mod_gfspn = p;
                     12: }
                     13:
                     14: void getmod_gfspn(up)
                     15: UM *up;
                     16: {
                     17:        *up = current_mod_gfspn;
                     18: }
                     19:
                     20: void simpgfspn(n,r)
                     21: GFSPN n;
                     22: GFSPN *r;
                     23: {
                     24:        UM t,q;
                     25:
                     26:        if ( !n )
                     27:                *r = 0;
                     28:        else if ( NID(n) != N_GFSPN )
                     29:                *r = n;
                     30:        else {
                     31:                t = UMALLOC(DEG(BDY(n)));
                     32:                q = W_UMALLOC(DEG(BDY(n)));
                     33:                cpyum(BDY(n),t);
                     34:                divsfum(t,current_mod_gfspn,q);
                     35:                MKGFSPN(t,*r);
                     36:        }
                     37: }
                     38:
                     39: #define NZGFSPN(a) ((a)&&(OID(a)==O_N)&&(NID(a)==N_GFSPN))
                     40:
                     41: void ntogfspn(a,b)
                     42: Obj a;
                     43: GFSPN *b;
                     44: {
                     45:        UM t;
                     46:        GFS c;
                     47:
                     48:        if ( !a || (OID(a)==O_N && NID(a) == N_GFSPN) )
                     49:                *b = (GFSPN)a;
                     50:        else if ( OID(a) == O_N && (NID(a) == N_GFS || NID(a) == N_Q) ) {
                     51:                ntogfs((Num)a,&c);
                     52:                if ( !b )
                     53:                        *b = 0;
                     54:                else {
                     55:                        t = UMALLOC(0); ptosfum((P)c,t);
                     56:                        MKGFSPN(t,*b);
                     57:                }
                     58:        } else
                     59:                error("ntogfspn : invalid argument");
                     60: }
                     61:
                     62: void addgfspn(a,b,c)
                     63: GFSPN a,b;
                     64: GFSPN *c;
                     65: {
                     66:        UM t,q;
                     67:        GFSPN z;
                     68:        int d;
                     69:
                     70:        ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
                     71:        if ( !a )
                     72:                *c = b;
                     73:        else if ( !b )
                     74:                *c = a;
                     75:        else {
                     76:                d = MAX(DEG(BDY(a)),DEG(BDY(b)));
                     77:                t = UMALLOC(d);
                     78:                q = W_UMALLOC(d);
                     79:                addsfum(BDY(a),BDY(b),t);
                     80:                divsfum(t,current_mod_gfspn,q);
                     81:                MKGFSPN(t,*c);
                     82:        }
                     83: }
                     84:
                     85: void subgfspn(a,b,c)
                     86: GFSPN a,b;
                     87: GFSPN *c;
                     88: {
                     89:        UM t,q;
                     90:        GFSPN z;
                     91:        int d;
                     92:
                     93:        ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
                     94:        if ( !a )
                     95:                chsgngfspn(b,c);
                     96:        else if ( !b )
                     97:                *c = a;
                     98:        else {
                     99:                d = MAX(DEG(BDY(a)),DEG(BDY(b)));
                    100:                t = UMALLOC(d);
                    101:                q = W_UMALLOC(d);
                    102:                subsfum(BDY(a),BDY(b),t);
                    103:                divsfum(t,current_mod_gfspn,q);
                    104:                MKGFSPN(t,*c);
                    105:        }
                    106: }
                    107:
                    108: extern int up_lazy;
                    109:
                    110: void mulgfspn(a,b,c)
                    111: GFSPN a,b;
                    112: GFSPN *c;
                    113: {
                    114:        UM t,q;
                    115:        GFSPN z;
                    116:        int d;
                    117:
                    118:        ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
                    119:        if ( !a || !b )
                    120:                *c = 0;
                    121:        else {
                    122:                d = DEG(BDY(a))+DEG(BDY(b));
                    123:                t = UMALLOC(d);
                    124:                q = W_UMALLOC(d);
                    125:                mulsfum(BDY(a),BDY(b),t);
                    126:                divsfum(t,current_mod_gfspn,q);
                    127:                MKGFSPN(t,*c);
                    128:        }
                    129: }
                    130:
                    131: void divgfspn(a,b,c)
                    132: GFSPN a,b;
                    133: GFSPN *c;
                    134: {
                    135:        GFSPN z;
                    136:        int d;
                    137:        UM wb,wc,wd,we,t,q;
                    138:
                    139:        ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
                    140:        if ( !b )
                    141:                error("divgfspn: division by 0");
                    142:        else if ( !a )
                    143:                *c = 0;
                    144:        else {
                    145:                wb = W_UMALLOC(DEG(BDY(b))); cpyum(BDY(b),wb);
                    146:                d = DEG(current_mod_gfspn);
                    147:                wc = W_UMALLOC(d); cpyum(current_mod_gfspn,wc);
                    148:                wd = W_UMALLOC(2*d); we = W_UMALLOC(2*d);
                    149:                /* wd*wb+we*wc=1 */
                    150:                eucsfum(wb,wc,wd,we);
                    151:                d = DEG(BDY(a))+DEG(wd);
                    152:                t = UMALLOC(d);
                    153:                q = W_UMALLOC(d);
                    154:                mulsfum(a,wd,t);
                    155:                divsfum(t,current_mod_gfspn,q);
                    156:                MKGFSPN(t,*c);
                    157:        }
                    158: }
                    159:
                    160: void invgfspn(b,c)
                    161: GFSPN b;
                    162: GFSPN *c;
                    163: {
                    164:        GFSPN z;
                    165:        int d;
                    166:        UM wb,wc,wd,we,t;
                    167:
                    168:        ntogfspn((Obj)b,&z); b = z;
                    169:        if ( !b )
                    170:                error("divgfspn: division by 0");
                    171:        else {
                    172:                wb = W_UMALLOC(DEG(BDY(b))); cpyum(BDY(b),wb);
                    173:                d = DEG(current_mod_gfspn);
                    174:                wc = W_UMALLOC(d); cpyum(current_mod_gfspn,wc);
                    175:                wd = W_UMALLOC(2*d); we = W_UMALLOC(2*d);
                    176:                /* wd*wb+we*wc=1 */
                    177:                eucsfum(wb,wc,wd,we);
                    178:                d = DEG(wd);
                    179:                t = UMALLOC(d);
                    180:                cpyum(wd,t);
                    181:                MKGFSPN(t,*c);
                    182:        }
                    183: }
                    184:
                    185: void chsgngfspn(a,c)
                    186: GFSPN a,*c;
                    187: {
                    188:        GFSPN z;
                    189:        int d;
                    190:        struct oUM zero;
                    191:        UM t;
                    192:
                    193:        ntogfspn((Obj)a,&z); a = z;
                    194:        if ( !a )
                    195:                *c = 0;
                    196:        else {
                    197:                d = DEG(BDY(a));
                    198:                t = UMALLOC(d);
                    199:                DEG(&zero) = -1;
                    200:                subsfum(&zero,BDY(a),t);
                    201:                MKGFSPN(t,*c);
                    202:        }
                    203: }
                    204:
                    205: void pwrgfspn(a,b,c)
                    206: GFSPN a;
                    207: Q b;
                    208: GFSPN *c;
                    209: {
                    210:        GFSPN z;
                    211:        UM t,x,y,q;
                    212:        int d,k;
                    213:        N e;
                    214:
                    215:        ntogfspn((Obj)a,&z); a = z;
                    216:        if ( !b ) {
                    217:                t = UMALLOC(0); DEG(t) = 0; COEF(t)[0] = _onesf();
                    218:                MKGFSPN(t,*c);
                    219:        } else if ( !a )
                    220:                *c = 0;
                    221:        else {
                    222:                d = DEG(current_mod_gfspn);
                    223:
                    224:                /* y = 1 */
                    225:                y = UMALLOC(d); DEG(y) = 0; COEF(y)[0] = _onesf();
                    226:
                    227:                t = W_UMALLOC(2*d); q = W_UMALLOC(2*d);
                    228:
                    229:                /* x = simplify(a) */
                    230:                x = W_UMALLOC(DEG(BDY(a))); cpyum(BDY(a),x);
                    231:                divsfum(x,current_mod_gfspn,q);
                    232:                if ( DEG(x) < 0 ) {
                    233:                        *c = 0;
                    234:                } else {
                    235:                        e = NM(b);
                    236:                        for ( k = n_bits(e)-1; k >= 0; k-- ) {
                    237:                                mulsfum(y,y,t);
                    238:                                divsfum(t,current_mod_gfspn,q);
                    239:                                cpyum(t,y);
                    240:                                if ( e->b[k/32] & (1<<(k%32)) ) {
                    241:                                        mulsfum(y,x,t);
                    242:                                        divsfum(t,current_mod_gfspn,q);
                    243:                                        cpyum(t,y);
                    244:                                }
                    245:                        }
                    246:                        MKGFSPN(y,*c);
                    247:                }
                    248:        }
                    249: }
                    250:
                    251: int cmpgfspn(a,b)
                    252: GFSPN a,b;
                    253: {
                    254:        GFSPN z;
                    255:
                    256:        ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
                    257:        if ( !a )
                    258:                if ( !b )
                    259:                        return 0;
                    260:                else
                    261:                        return -1;
                    262:        else if ( !b )
                    263:                        return 1;
                    264:        else
                    265:                return compsfum(BDY(a),BDY(b));
                    266: }
                    267:
                    268: void randomgfspn(r)
                    269: GFSPN *r;
                    270: {
                    271:        int i,d;
                    272:        UM t;
                    273:
                    274:        if ( !current_mod_gfspn )
                    275:                error("randomgfspn : current_mod_gfspn is not set");
                    276:        d = DEG(current_mod_gfspn);
                    277:        t = UMALLOC(d-1);
                    278:        randsfum(d,t);
                    279:        degum(t,d-1);
                    280:        if ( DEG(t) < 0 )
                    281:                *r = 0;
                    282:        else {
                    283:                MKGFSPN(t,*r);
                    284:        }
                    285: }

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