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

Annotation of OpenXM_contrib2/asir2018/engine/poly.c, Revision 1.2

1.1       noro        1: #ifdef MODULAR
                      2: void addmp(VL vl,int mod,P p1,P p2,P *pr)
                      3: #else
                      4: void addp(VL vl,P p1,P p2,P *pr)
                      5: #endif
                      6: {
                      7:   DCP dc1,dc2,dcr0,dcr;
                      8:   V v1,v2;
                      9:   P t;
                     10:
                     11:   if ( !p1 )
                     12:     *pr = p2;
                     13:   else if ( !p2 )
                     14:     *pr = p1;
                     15:   else if ( NUM(p1) )
                     16:     if ( NUM(p2) )
                     17:       ADDNUM(p1,p2,pr);
                     18:     else
                     19:       ADDPQ(p2,p1,pr);
                     20:   else if ( NUM(p2) )
                     21:     ADDPQ(p1,p2,pr);
                     22:   else if ( ( v1 = VR(p1) ) ==  ( v2 = VR(p2) ) ) {
                     23:     for ( dc1 = DC(p1), dc2 = DC(p2), dcr0 = 0; dc1 && dc2; )
                     24:       switch ( cmpz(DEG(dc1),DEG(dc2)) ) {
                     25:         case 0:
                     26:           ADDP(vl,COEF(dc1),COEF(dc2),&t);
                     27:           if ( t )  {
                     28:             NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc1); COEF(dcr) = t;
                     29:           }
                     30:           dc1 = NEXT(dc1); dc2 = NEXT(dc2); break;
                     31:         case 1:
                     32:           NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc1); COEF(dcr) = COEF(dc1);
                     33:           dc1 = NEXT(dc1); break;
                     34:         case -1:
                     35:           NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc2); COEF(dcr) = COEF(dc2);
                     36:           dc2 = NEXT(dc2); break;
                     37:       }
                     38:     if ( !dcr0 )
                     39:       if ( dc1 )
                     40:         dcr0 = dc1;
                     41:       else if ( dc2 )
                     42:         dcr0 = dc2;
                     43:       else {
                     44:         *pr = 0;
                     45:         return;
                     46:       }
                     47:     else
                     48:       if ( dc1 )
                     49:         NEXT(dcr) = dc1;
                     50:       else if ( dc2 )
                     51:         NEXT(dcr) = dc2;
                     52:       else
                     53:         NEXT(dcr) = 0;
                     54:     MKP(v1,dcr0,*pr);
                     55:   } else {
                     56:     while ( v1 != VR(vl) && v2 != VR(vl) )
                     57:       vl = NEXT(vl);
                     58:     if ( v1 == VR(vl) )
                     59:       ADDPTOC(vl,p1,p2,pr);
                     60:     else
                     61:       ADDPTOC(vl,p2,p1,pr);
                     62:   }
                     63: }
                     64:
                     65: #ifdef MODULAR
                     66: void addmpq(int mod,P p,P q,P *pr)
                     67: #else
                     68: void addpq(P p,P q,P *pr)
                     69: #endif
                     70: {
                     71:   DCP dc,dcr,dcr0;
                     72:   P t;
                     73:
                     74:   if ( NUM(p) )
                     75:     ADDNUM(p,q,pr);
                     76:   else {
                     77:     dc = DC(p);
                     78:     for ( dcr0 = 0; dc && cmpz(DEG(dc),0) > 0; dc = NEXT(dc) ) {
                     79:       NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = COEF(dc);
                     80:     }
                     81:     if ( !dc ) {
                     82:       NEXTDC(dcr0,dcr);
                     83:       DEG(dcr) = 0; COEF(dcr) = q; NEXT(dcr) = 0;
                     84:     } else if ( !DEG(dc) ) {
                     85:       ADDPQ(COEF(dc),q,&t);
                     86:       if ( t ) {
                     87:         NEXTDC(dcr0,dcr);
                     88:         DEG(dcr) = 0; COEF(dcr) = t;
                     89:       }
                     90:       NEXT(dcr) = NEXT(dc);
                     91:     } else {
                     92:       NEXTDC(dcr0,dcr);
                     93:       DEG(dcr) = 0; COEF(dcr) = q; NEXT(dcr) = dc;
                     94:     }
                     95:     MKP(VR(p),dcr0,*pr);
                     96:   }
                     97: }
                     98:
                     99: #ifdef MODULAR
                    100: void addmptoc(VL vl,int mod,P p,P c,P *pr)
                    101: #else
                    102: void addptoc(VL vl,P p,P c,P *pr)
                    103: #endif
                    104: {
                    105:   DCP dc,dcr,dcr0;
                    106:   P t;
                    107:
                    108:   for ( dcr0 = 0, dc = DC(p); dc && cmpz(DEG(dc),0) > 0; dc = NEXT(dc) ) {
                    109:     NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = COEF(dc);
                    110:   }
                    111:   if ( !dc ) {
                    112:     NEXTDC(dcr0,dcr);
                    113:     DEG(dcr) = 0; COEF(dcr) = c; NEXT(dcr) = 0;
                    114:   } else if ( !DEG(dc) ) {
                    115:     ADDP(vl,COEF(dc),c,&t);
                    116:     if ( t ) {
                    117:       NEXTDC(dcr0,dcr);
                    118:       DEG(dcr) = 0; COEF(dcr) = t;
                    119:     }
                    120:     NEXT(dcr) = NEXT(dc);
                    121:   } else {
                    122:     NEXTDC(dcr0,dcr);
                    123:     DEG(dcr) = 0; COEF(dcr) = c; NEXT(dcr) = dc;
                    124:   }
                    125:   MKP(VR(p),dcr0,*pr);
                    126: }
                    127:
                    128: #ifdef MODULAR
                    129: void submp(VL vl,int mod,P p1,P p2,P *pr)
                    130: #else
                    131: void subp(VL vl,P p1,P p2,P *pr)
                    132: #endif
                    133: {
                    134:   P t;
                    135:
                    136:   if ( !p2 )
                    137:     *pr = p1;
                    138:   else {
                    139:     CHSGNP(p2,&t); ADDP(vl,p1,t,pr);
                    140:   }
                    141: }
                    142:
                    143: #ifdef MODULAR
                    144: void mulmp(VL vl,int mod,P p1,P p2,P *pr)
                    145: #else
                    146: void mulp(VL vl,P p1,P p2,P *pr)
                    147: #endif
                    148: {
                    149:   DCP dc,dct,dcr,dcr0;
                    150:   V v1,v2;
                    151:   P t,s,u;
                    152:   int n1,n2;
                    153:
                    154:   if ( !p1 || !p2 ) *pr = 0;
                    155:   else if ( NUM(p1) )
                    156:     MULPQ(p2,p1,pr);
                    157:   else if ( NUM(p2) )
                    158:     MULPQ(p1,p2,pr);
                    159:   else if ( ( v1 = VR(p1) ) ==  ( v2 = VR(p2) ) ) {
                    160:     for ( dc = DC(p1), n1 = 0; dc; dc = NEXT(dc), n1++ );
                    161:     for ( dc = DC(p2), n2 = 0; dc; dc = NEXT(dc), n2++ );
                    162:     if ( n1 > n2 )
                    163:       for ( dc = DC(p2), s = 0; dc; dc = NEXT(dc) ) {
                    164:         for ( dcr0 = 0, dct = DC(p1); dct; dct = NEXT(dct) ) {
                    165:           NEXTDC(dcr0,dcr); MULP(vl,COEF(dct),COEF(dc),&COEF(dcr));
                    166:           addz(DEG(dct),DEG(dc),&DEG(dcr));
                    167:         }
                    168:         NEXT(dcr) = 0; MKP(v1,dcr0,t);
                    169:         ADDP(vl,s,t,&u); s = u; t = u = 0;
                    170:       }
                    171:     else
                    172:       for ( dc = DC(p1), s = 0; dc; dc = NEXT(dc) ) {
                    173:         for ( dcr0 = 0, dct = DC(p2); dct; dct = NEXT(dct) ) {
                    174:           NEXTDC(dcr0,dcr); MULP(vl,COEF(dct),COEF(dc),&COEF(dcr));
                    175:           addz(DEG(dct),DEG(dc),&DEG(dcr));
                    176:         }
                    177:         NEXT(dcr) = 0; MKP(v1,dcr0,t);
                    178:         ADDP(vl,s,t,&u); s = u; t = u = 0;
                    179:       }
                    180:     *pr = s;
                    181:   } else {
                    182:     while ( v1 != VR(vl) && v2 != VR(vl) )
                    183:       vl = NEXT(vl);
                    184:     if ( v1 == VR(vl) )
                    185:       MULPC(vl,p1,p2,pr);
                    186:     else
                    187:       MULPC(vl,p2,p1,pr);
                    188:   }
                    189: }
                    190:
                    191: #ifdef MODULAR
                    192: void mulmpq(int mod,P p,P q,P *pr)
                    193: #else
                    194: void mulpq(P p,P q,P *pr)
                    195: #endif
                    196: {
                    197:   DCP dc,dcr,dcr0;
                    198:   P t;
                    199:
                    200:   if (!p || !q)
                    201:     *pr = 0;
                    202:   else if ( Uniq(q) )
                    203:     *pr = p;
                    204:   else if ( NUM(p) )
                    205:     MULNUM(p,q,pr);
                    206:   else {
                    207:     for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
                    208:       MULPQ(COEF(dc),q,&t);
                    209:       if ( t ) {
                    210:         NEXTDC(dcr0,dcr); COEF(dcr) = t; DEG(dcr) = DEG(dc);
                    211:       }
                    212:     }
                    213:     if ( dcr0 ) {
                    214:       NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
                    215:     } else
                    216:       *pr = 0;
                    217:   }
                    218: }
                    219:
                    220: #ifdef MODULAR
                    221: void mulmpc(VL vl,int mod,P p,P c,P *pr)
                    222: #else
                    223: void mulpc(VL vl,P p,P c,P *pr)
                    224: #endif
                    225: {
                    226:   DCP dc,dcr,dcr0;
                    227:   P t;
                    228:
                    229:   if ( NUM(c) )
                    230:     MULPQ(p,c,pr);
                    231:   else {
                    232:     for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
                    233:       MULP(vl,COEF(dc),c,&t);
                    234:       if ( t ) {
                    235:         NEXTDC(dcr0,dcr); COEF(dcr) = t; DEG(dcr) = DEG(dc);
                    236:       }
                    237:     }
                    238:     if ( dcr0 ) {
                    239:       NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
                    240:     } else
                    241:       *pr = 0;
                    242:   }
                    243: }
                    244:
                    245: #ifdef MODULAR
                    246: void pwrmp(VL vl,int mod,P p,Z q,P *pr)
                    247: #else
                    248: void pwrp(VL vl,P p,Z q,P *pr)
                    249: #endif
                    250: {
                    251:   DCP dc,dcr;
                    252:   int n,i;
                    253:   P *x,*y;
                    254:   P t,s,u;
                    255:   DCP dct;
                    256:   P *pt;
                    257:
                    258:   if ( !q ) {
                    259:     *pr = (P)One;
                    260:   } else if ( !p )
                    261:     *pr = 0;
                    262:   else if ( UNIQ(q) )
                    263:     *pr = p;
                    264:   else if ( NUM(p) )
                    265:     PWRNUM(p,(Q)q,pr);
                    266:   else {
                    267:     dc = DC(p);
                    268:     if ( !NEXT(dc) ) {
                    269:       NEWDC(dcr);
                    270:       PWRP(vl,COEF(dc),q,&COEF(dcr)); mulz(DEG(dc),(Z)q,&DEG(dcr));
                    271:       NEXT(dcr) = 0; MKP(VR(p),dcr,*pr);
                    272:     } else if ( !INT(q) ) {
                    273:       error("pwrp: can't calculate fractional power."); *pr = 0;
                    274:     } else if ( smallz(q) ) {
1.2     ! noro      275:       n = ZTOS(q); x = (P *)ALLOCA((n+1)*sizeof(pointer));
1.1       noro      276:       NEWDC(dct); DEG(dct) = DEG(dc); COEF(dct) = COEF(dc);
                    277:       NEXT(dct) = 0; MKP(VR(p),dct,t);
                    278:       for ( i = 0, u = (P)One; i < n; i++ ) {
                    279:         x[i] = u; MULP(vl,u,t,&s); u = s;
                    280:       }
                    281:       x[n] = u; y = (P *)ALLOCA((n+1)*sizeof(pointer));
                    282:
                    283:       MKP(VR(p),NEXT(dc),t);
                    284:       for ( i = 0, u = (P)One; i < n; i++ ) {
                    285:         y[i] = u; MULP(vl,u,t,&s); u = s;
                    286:       }
                    287:       y[n] = u;
                    288:       pt = (P *)ALLOCA((n+1)*sizeof(pointer)); MKBC(n,pt);
                    289:       for ( i = 0, u = 0; i <= n; i++ ) {
                    290:         MULP(vl,x[i],y[n-i],&t); MULP(vl,t,pt[i],&s);
                    291:         ADDP(vl,u,s,&t); u = t;
                    292:       }
                    293:       *pr = u;
                    294:     } else {
                    295:       error("exponent too big");
                    296:       *pr = 0;
                    297:     }
                    298:   }
                    299: }
                    300:
                    301: #ifdef MODULAR
                    302: void chsgnmp(int mod,P p,P *pr)
                    303: #else
                    304: void chsgnp(P p,P *pr)
                    305: #endif
                    306: {
                    307:   register DCP dc,dcr,dcr0;
                    308:
                    309:   if ( !p )
                    310:     *pr = NULL;
                    311:   else if ( NUM(p) ) {
                    312: #if defined(_PA_RISC1_1) || defined(__alpha) || defined(mips) || defined(_IBMR2)
                    313: #ifdef FBASE
                    314:     chsgnnum((Num)p,(Num *)pr);
                    315: #else
                    316:     MQ mq;
                    317:
                    318:     NEWMQ(mq); CONT(mq)=mod-CONT((MQ)p); *pr = (P)mq;
                    319: #endif
                    320: #else
                    321:     CHSGNNUM(p,*pr);
                    322: #endif
                    323:   } else {
                    324:     for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
                    325:       NEXTDC(dcr0,dcr); CHSGNP(COEF(dc),&COEF(dcr)); DEG(dcr) = DEG(dc);
                    326:     }
                    327:     NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
                    328:   }
                    329: }
                    330:

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