Return to polydiv.c CVS log

Up to [local] / OpenXM_contrib2 / asir2000 / engine

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

1.1       noro        1: void D_DIVSRP(vl,p1,p2,q,r)
                      2: VL vl;
                      3: P p1,p2;
                      4: P *q,*r;
                      5: {
1.2     ! noro        6:   register int i,j;
        !             7:   register DCP dc1,dc2,dc;
        !             8:   P m,s,dvr;
        !             9:   P *pq,*pr,*pd;
        !            10:   V v1,v2;
        !            11:   Q deg1,deg2;
        !            12:   int d1,d2,sgn;
        !            13:
        !            14:   if ( !p1 ) {
        !            15:     *q = 0; *r = 0;
        !            16:   } else if ( NUM(p2) )
        !            17:     if ( NUM(p1) ) {
        !            18:       DIVNUM(p1,p2,q); *r = 0;
        !            19:     } else {
        !            20:       DIVSDCP(vl,p1,p2,q); *r = 0;
        !            21:     }
        !            22:   else if ( NUM(p1) ) {
        !            23:     *q = 0; *r = p1;
        !            24:   } else if ( ( v1 = VR(p1) ) == ( v2 = VR(p2) ) ) {
        !            25:     dc1 = DC(p1); dc2 = DC(p2);
        !            26:     deg1 = DEG(dc1); deg2 = DEG(dc2);
        !            27:     sgn = cmpq(deg1,deg2);
        !            28:     if ( sgn == 0 ) {
        !            29:       DIVSP(vl,COEF(dc1),COEF(dc2),q);
        !            30:       MULP(vl,p2,*q,&m); SUBP(vl,p1,m,r);
        !            31:     } else if ( sgn < 0 ) {
        !            32:       *q = 0; *r = p1;
        !            33:     } else {
        !            34:       if ( (PL(NM(deg1)) > 1) )
        !            35:         error("divsrp : invalid input");
        !            36:       d1 = QTOS(deg1); d2 = QTOS(deg2);
        !            37:       W_CALLOC(d1-d2,P,pq); W_CALLOC(d1,P,pr); W_CALLOC(d2,P,pd);
        !            38:       for ( dc = dc1; dc; dc = NEXT(dc) )
        !            39:         pr[QTOS(DEG(dc))] = COEF(dc);
        !            40:       for ( dc = dc2; dc; dc = NEXT(dc) )
        !            41:         pd[QTOS(DEG(dc))] = COEF(dc);
        !            42:       for ( dvr = COEF(dc2), i = d1 - d2; i >= 0; i-- ) {
        !            43:         if ( !pr[i+d2] )
        !            44:           continue;
        !            45:         DIVSP(vl,pr[i+d2],dvr,&pq[i]);
        !            46:         for ( j = d2; j >= 0; j-- ) {
        !            47:           MULP(vl,pq[i],pd[j],&m);
        !            48:           SUBP(vl,pr[i + j],m,&s); pr[i + j] = s;
        !            49:         }
        !            50:       }
        !            51:       plisttop(pq,v1,d1 - d2,q); plisttop(pr,v1,d1 - 1,r);
        !            52:     }
        !            53:   } else {
        !            54:     for ( ; (v1 != vl->v) && (v2 != vl->v); vl = NEXT(vl) );
        !            55:     if ( v2 == vl->v ) {
        !            56:       *q = 0; *r = p1;
        !            57:     } else
        !            58:       DIVSRDCP(vl,p1,p2,q,r);
        !            59:   }
1.1       noro       60: }
                     61:
                     62: void D_DIVSRDCP(vl,p1,p2,q,r)
                     63: VL vl;
                     64: P p1,p2,*q,*r;
                     65: {
                     66:
1.2     ! noro       67:   P qc,rc;
        !            68:   DCP dc,dcq,dcq0,dcr,dcr0;
1.1       noro       69:
1.2     ! noro       70:   for ( dc = DC(p1), dcq0 = 0, dcr0 = 0; dc; dc = NEXT(dc) ) {
        !            71:     DIVSRP(vl,COEF(dc),p2,&qc,&rc);
        !            72:     if ( qc ) {
        !            73:       NEXTDC(dcq0,dcq); DEG(dcq) = DEG(dc); COEF(dcq) = qc;
        !            74:     }
        !            75:     if ( rc ) {
        !            76:       NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = rc;
        !            77:     }
        !            78:   }
        !            79:   if ( dcq0 ) {
        !            80:     NEXT(dcq) = 0; MKP(VR(p1),dcq0,*q);
        !            81:   } else
        !            82:     *q = 0;
        !            83:   if ( dcr0 ) {
        !            84:     NEXT(dcr) = 0; MKP(VR(p1),dcr0,*r);
        !            85:   } else
        !            86:     *r = 0;
1.1       noro       87: }
                     88:
                     89: void D_DIVSP(vl,p1,p2,q)
                     90: VL vl;
                     91: P p1,p2,*q;
                     92: {
1.2     ! noro       93:   P t;
1.1       noro       94:
1.2     ! noro       95:   DIVSRP(vl,p1,p2,q,&t);
        !            96:   if ( t )
        !            97:     error("divsp: cannot happen");
1.1       noro       98: }
                     99:
                    100: void D_DIVSDCP(vl,p1,p2,q)
                    101: VL vl;
                    102: P p1,p2,*q;
                    103: {
                    104:
1.2     ! noro      105:   P m;
        !           106:   register DCP dc,dcr,dcr0;
1.1       noro      107:
1.2     ! noro      108:   for ( dc = DC(p1), dcr0 = 0; dc; dc = NEXT(dc) ) {
        !           109:     DIVSP(vl,COEF(dc),p2,&m);
        !           110:     NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = m; NEXT(dcr) = 0;
        !           111:   }
        !           112:   MKP(VR(p1),dcr0,*q);
1.1       noro      113: }
                    114: