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

Annotation of OpenXM_contrib2/asir2000/engine/dalg.c, Revision 1.6

1.1       noro        1: /*
1.6     ! noro        2:  * $OpenXM: OpenXM_contrib2/asir2000/engine/dalg.c,v 1.5 2004/12/04 09:39:27 noro Exp $
1.1       noro        3: */
                      4:
                      5: #include "ca.h"
                      6: #include "base.h"
                      7:
                      8: static NumberField current_numberfield;
                      9: extern struct order_spec *dp_current_spec;
1.2       noro       10: void simpdalg(DAlg da,DAlg *r);
1.3       noro       11: void invdalg(DAlg a,DAlg *c);
                     12: void rmcontdalg(DAlg a, DAlg *c);
1.6     ! noro       13: void algtodalg(Alg a,DAlg *r);
        !            14: void dalgtoalg(DAlg da,Alg *r);
1.1       noro       15:
1.5       noro       16: NumberField get_numberfield()
                     17: {
                     18:        return current_numberfield;
                     19: }
                     20:
1.1       noro       21: void setfield_dalg(NODE alist)
                     22: {
                     23:        NumberField nf;
                     24:        VL vl,vl1,vl2;
                     25:        int n,i,dim;
                     26:        Alg *gen;
                     27:        P *defpoly;
                     28:        P p;
                     29:        Q c,iq,two;
                     30:        DP *ps,*mb;
1.3       noro       31:        DP one;
1.1       noro       32:        NODE t,b,b1,b2,hlist,mblist;
                     33:        struct order_spec *current_spec;
                     34:
                     35:        nf = (NumberField)MALLOC(sizeof(struct oNumberField));
                     36:        current_numberfield = nf;
                     37:        vl = 0;
                     38:        for ( t = alist; t; t = NEXT(t) ) {
                     39:                clctalg(BDY((Alg)BDY(t)),&vl1);
                     40:                mergev(ALG,vl,vl1,&vl2); vl = vl2;
                     41:        }
                     42:        for ( n = 0, vl1 = vl; vl1; vl1 = NEXT(vl1), n++ );
                     43:        nf->n = n;
                     44:        nf->vl = vl;
                     45:        nf->defpoly = defpoly = (P *)MALLOC(n*sizeof(P));
                     46:        nf->ps = ps = (DP *)MALLOC(n*sizeof(DP));
                     47:        current_spec = dp_current_spec;
                     48:        STOQ(2,two);
                     49:        create_order_spec(0,(Obj)two,&nf->spec);
                     50:        initd(nf->spec);
                     51:        for ( b = hlist = 0, i = 0, vl1 = vl; i < n; vl1 = NEXT(vl1), i++ ) {
                     52:                ptozp(vl1->v->attr,1,&c,&defpoly[i]);
                     53:                ptod(ALG,vl,defpoly[i],&ps[i]);
                     54:                STOQ(i,iq); MKNODE(b1,(pointer)iq,b); b = b1;
1.3       noro       55:                MKNODE(b2,(pointer)ps[i],hlist); hlist = b2;
1.1       noro       56:        }
1.3       noro       57:        ptod(ALG,vl,(P)ONE,&one);
                     58:        MKDAlg(one,ONE,nf->one);
                     59:        nf->ind = b;
                     60:        dp_mbase(hlist,&mblist);
1.1       noro       61:        initd(current_spec);
                     62:        nf->dim = dim = length(mblist);
                     63:        nf->mb = mb = (DP *)MALLOC(dim*sizeof(DP));
                     64:        for ( i = 0, t = mblist; t; t = NEXT(t), i++ )
1.3       noro       65:                mb[i] = (DP)BDY(t);
1.1       noro       66: }
                     67:
1.4       noro       68: void qtodalg(Q q,DAlg *r)
                     69: {
                     70:        NumberField nf;
                     71:        Q t;
                     72:        DP nm;
                     73:
                     74:        if ( !(nf=current_numberfield) )
                     75:                error("qtodalg : current_numberfield is not set");
                     76:        if ( !q )
                     77:                *r = 0;
                     78:        else if ( NID(q) == N_DA )
                     79:                *r = (DAlg)q;
                     80:        else if ( NID(q) == N_Q ) {
                     81:                if ( INT(q) ) {
                     82:                        muldc(CO,nf->one->nm,(P)q,&nm);
                     83:                        MKDAlg(nm,ONE,*r);
                     84:                } else {
                     85:                        NTOQ(NM(q),SGN(q),t);
                     86:                        muldc(CO,nf->one->nm,(P)t,&nm);
                     87:                        NTOQ(DN(q),1,t);
                     88:                        MKDAlg(nm,t,*r);
                     89:                }
                     90:        } else
                     91:                error("qtodalg : invalid argument");
1.6     ! noro       92: }
        !            93:
        !            94: void obj_algtodalg(Obj obj,Obj *r)
        !            95: {
        !            96:        DAlg d;
        !            97:        DCP dc,dcr0,dcr;
        !            98:        P c,p;
        !            99:        Obj t;
        !           100:        Obj nm,dn;
        !           101:        NODE b,s,s0;
        !           102:        R rat;
        !           103:        VECT v;
        !           104:        MAT mat;
        !           105:        LIST list;
        !           106:        pointer *a;
        !           107:        pointer **m;
        !           108:        int len,row,col,i,j,l;
        !           109:
        !           110:        if ( !obj ) {
        !           111:                *r = 0;
        !           112:                return;
        !           113:        }
        !           114:        switch ( OID(obj) ) {
        !           115:                case O_N:
        !           116:                        algtodalg((Alg)obj,&d); *r = (Obj)d;
        !           117:                        break;
        !           118:                case O_P:
        !           119:                        for ( dcr0 = 0, dc = DC((P)obj); dc; dc = NEXT(dc) ) {
        !           120:                                obj_algtodalg((Obj)COEF(dc),&t);
        !           121:                                if ( t ) {
        !           122:                                        NEXTDC(dcr0,dcr);
        !           123:                                        COEF(dcr) = (P)t;
        !           124:                                        DEG(dcr) = DEG(dc);
        !           125:                                }
        !           126:                        }
        !           127:                        if ( dcr0 ) {
        !           128:                                MKP(VR((P)obj),dcr0,p);
        !           129:                                *r = (Obj)p;
        !           130:                        } else
        !           131:                                *r = 0;
        !           132:                        break;
        !           133:                case O_R:
        !           134:                        obj_algtodalg((Obj)NM((R)obj),&nm);
        !           135:                        obj_algtodalg((Obj)DN((R)obj),&dn);
        !           136:                        if ( !dn )
        !           137:                                error("obj_algtodalg : division by 0");
        !           138:                        if ( !nm )
        !           139:                                *r = 0;
        !           140:                        else {
        !           141:                                MKRAT((P)nm,(P)dn,0,rat); *r = (Obj)rat;
        !           142:                        }
        !           143:                        break;
        !           144:                case O_LIST:
        !           145:                        s0 = 0;
        !           146:                        for ( b = BDY((LIST)obj); b; b = NEXT(b) ) {
        !           147:                                NEXTNODE(s0,s);
        !           148:                                obj_algtodalg((Obj)BDY(b),&t);
        !           149:                                BDY(s) = (pointer)t;
        !           150:                        }
        !           151:                        NEXT(s) = 0;
        !           152:                        MKLIST(list,s0);
        !           153:                        *r = (Obj)list;
        !           154:                        break;
        !           155:                case O_VECT:
        !           156:                        l = ((VECT)obj)->len;
        !           157:                        a = BDY((VECT)obj);
        !           158:                        MKVECT(v,l);
        !           159:                        for ( i = 0; i < l; i++ ) {
        !           160:                                obj_algtodalg((Obj)a[i],&t);
        !           161:                                BDY(v)[i] = (pointer)t;
        !           162:                        }
        !           163:                        *r = (Obj)v;
        !           164:                        break;
        !           165:                case O_MAT:
        !           166:                        row = ((MAT)obj)->row; col = ((MAT)obj)->col;
        !           167:                        m = BDY((MAT)obj);
        !           168:                        MKMAT(mat,row,col);
        !           169:                        for ( i = 0; i < row; i++ )
        !           170:                                for ( j = 0; j < col; j++ ) {
        !           171:                                        obj_algtodalg((Obj)m[i][j],&t);
        !           172:                                        BDY(mat)[i][j] = (pointer)t;
        !           173:                                }
        !           174:                        *r = (Obj)mat;
        !           175:                        break;
        !           176:                default:
        !           177:                        *r = obj;
        !           178:                        break;
        !           179:        }
        !           180: }
        !           181:
        !           182: void obj_dalgtoalg(Obj obj,Obj *r)
        !           183: {
        !           184:        Alg d;
        !           185:        DCP dc,dcr0,dcr;
        !           186:        P c,p;
        !           187:        Obj t;
        !           188:        Obj nm,dn;
        !           189:        NODE b,s,s0;
        !           190:        R rat;
        !           191:        VECT v;
        !           192:        MAT mat;
        !           193:        LIST list;
        !           194:        pointer *a;
        !           195:        pointer **m;
        !           196:        int len,row,col,i,j,l;
        !           197:
        !           198:        if ( !obj ) {
        !           199:                *r = 0;
        !           200:                return;
        !           201:        }
        !           202:        switch ( OID(obj) ) {
        !           203:                case O_N:
        !           204:                        dalgtoalg((DAlg)obj,&d); *r = (Obj)d;
        !           205:                        break;
        !           206:                case O_P:
        !           207:                        for ( dcr0 = 0, dc = DC((P)obj); dc; dc = NEXT(dc) ) {
        !           208:                                obj_dalgtoalg((Obj)COEF(dc),&t);
        !           209:                                if ( t ) {
        !           210:                                        NEXTDC(dcr0,dcr);
        !           211:                                        COEF(dcr) = (P)t;
        !           212:                                        DEG(dcr) = DEG(dc);
        !           213:                                }
        !           214:                        }
        !           215:                        if ( dcr0 ) {
        !           216:                                MKP(VR((P)obj),dcr0,p);
        !           217:                                *r = (Obj)p;
        !           218:                        } else
        !           219:                                *r = 0;
        !           220:                        break;
        !           221:                case O_R:
        !           222:                        obj_dalgtoalg((Obj)NM((R)obj),&nm);
        !           223:                        obj_dalgtoalg((Obj)DN((R)obj),&dn);
        !           224:                        if ( !dn )
        !           225:                                error("obj_dalgtoalg : division by 0");
        !           226:                        if ( !nm )
        !           227:                                *r = 0;
        !           228:                        else {
        !           229:                                MKRAT((P)nm,(P)dn,0,rat); *r = (Obj)rat;
        !           230:                        }
        !           231:                        break;
        !           232:                case O_LIST:
        !           233:                        s0 = 0;
        !           234:                        for ( b = BDY((LIST)obj); b; b = NEXT(b) ) {
        !           235:                                NEXTNODE(s0,s);
        !           236:                                obj_dalgtoalg((Obj)BDY(b),&t);
        !           237:                                BDY(s) = (pointer)t;
        !           238:                        }
        !           239:                        NEXT(s) = 0;
        !           240:                        MKLIST(list,s0);
        !           241:                        *r = (Obj)list;
        !           242:                        break;
        !           243:                case O_VECT:
        !           244:                        l = ((VECT)obj)->len;
        !           245:                        a = BDY((VECT)obj);
        !           246:                        MKVECT(v,l);
        !           247:                        for ( i = 0; i < l; i++ ) {
        !           248:                                obj_dalgtoalg((Obj)a[i],&t);
        !           249:                                BDY(v)[i] = (pointer)t;
        !           250:                        }
        !           251:                        *r = (Obj)v;
        !           252:                        break;
        !           253:                case O_MAT:
        !           254:                        row = ((MAT)obj)->row; col = ((MAT)obj)->col;
        !           255:                        m = BDY((MAT)obj);
        !           256:                        MKMAT(mat,row,col);
        !           257:                        for ( i = 0; i < row; i++ )
        !           258:                                for ( j = 0; j < col; j++ ) {
        !           259:                                        obj_dalgtoalg((Obj)m[i][j],&t);
        !           260:                                        BDY(mat)[i][j] = (pointer)t;
        !           261:                                }
        !           262:                        *r = (Obj)mat;
        !           263:                        break;
        !           264:                default:
        !           265:                        *r = obj;
        !           266:                        break;
        !           267:        }
1.4       noro      268: }
                    269:
1.1       noro      270: void algtodalg(Alg a,DAlg *r)
                    271: {
                    272:        P ap,p,p1;
1.3       noro      273:        Q c,c1,d1,dn,nm;
1.1       noro      274:        DP dp;
                    275:        DAlg da;
                    276:        NumberField nf;
                    277:        struct order_spec *current_spec;
1.3       noro      278:        VL vl,tvl,svl;
                    279:        V v;
1.1       noro      280:
                    281:        if ( !(nf=current_numberfield) )
                    282:                error("algtodalg : current_numberfield is not set");
1.3       noro      283:        if ( !a ) {
                    284:                *r = 0;
                    285:                return;
                    286:        }
                    287:        switch (NID((Num)a) ) {
                    288:                case N_Q:
                    289:                        c = (Q)a;
                    290:                        if ( INT(c) ) {
                    291:                                muldc(CO,nf->one->nm,(P)c,&dp);
                    292:                                MKDAlg(dp,ONE,*r);
                    293:                        } else {
                    294:                                NTOQ(NM(c),SGN(c),c1);
                    295:                                NTOQ(DN(c),1,d1);
                    296:                                muldc(CO,nf->one->nm,(P)c1,&dp);
                    297:                                MKDAlg(dp,c1,*r);
                    298:                        }
                    299:                        break;
                    300:                case N_A:
                    301:                        ap = (P)BDY(a);
                    302:                        ptozp(ap,1,&c,&p);
                    303:                        if ( INT(c) ) {
                    304:                                p = ap;
                    305:                                dn = ONE;
                    306:                        } else {
                    307:                                NTOQ(NM(c),SGN(c),nm);
                    308:                                NTOQ(DN(c),1,dn);
                    309:                                mulpq(p,(P)nm,&p1); p = p1;
                    310:                        }
                    311:                        current_spec = dp_current_spec; initd(nf->spec);
                    312:                        get_vars(p,&vl);
                    313:                        for ( tvl = vl; tvl; tvl = NEXT(tvl) ) {
                    314:                                v = tvl->v;
                    315:                                for ( svl = nf->vl; svl; svl = NEXT(svl) )
                    316:                                        if ( v == svl->v )
                    317:                                                break;
                    318:                                if ( !svl )
                    319:                                        error("algtodalg : incompatible numberfield");
                    320:                        }
                    321:                        ptod(ALG,nf->vl,p,&dp);
                    322:                        MKDAlg(dp,dn,da);
                    323:                        simpdalg(da,r);
                    324:                        break;
                    325:                default:
                    326:                        error("algtodalg : invalid argument");
                    327:                        break;
1.1       noro      328:        }
                    329: }
                    330:
1.2       noro      331: void dalgtoalg(DAlg da,Alg *r)
1.1       noro      332: {
1.2       noro      333:        NumberField nf;
                    334:        P p,p1;
                    335:        Q inv;
                    336:
                    337:        if ( !(nf=current_numberfield) )
1.3       noro      338:                error("dalgtoalg : current_numberfield is not set");
1.2       noro      339:        dtop(ALG,nf->vl,da->nm,&p);
                    340:        invq(da->dn,&inv);
                    341:        mulpq(p,(P)inv,&p1);
                    342:        MKAlg(p1,*r);
1.1       noro      343: }
                    344:
                    345: void simpdalg(DAlg da,DAlg *r)
                    346: {
1.2       noro      347:        NumberField nf;
                    348:        DP nm;
1.3       noro      349:        DAlg d;
1.2       noro      350:        Q dn,dn1;
1.5       noro      351:        struct order_spec *current_spec;
1.2       noro      352:
                    353:        if ( !(nf=current_numberfield) )
1.3       noro      354:                error("simpdalg : current_numberfield is not set");
                    355:        if ( !da ) {
                    356:                *r = 0;
                    357:                return;
                    358:        }
1.5       noro      359:        current_spec = dp_current_spec; initd(nf->spec);
1.2       noro      360:        dp_true_nf(nf->ind,da->nm,nf->ps,1,&nm,&dn);
1.5       noro      361:        initd(current_spec);
1.2       noro      362:        mulq(da->dn,dn,&dn1);
1.3       noro      363:        MKDAlg(nm,dn1,d);
                    364:        rmcontdalg(d,r);
1.1       noro      365: }
                    366:
                    367: void adddalg(DAlg a,DAlg b,DAlg *c)
                    368: {
1.3       noro      369:        NumberField nf;
                    370:        Q dna,dnb,a1,b1,dn,g;
                    371:        N an,bn,gn;
1.4       noro      372:        DAlg t;
1.3       noro      373:        DP ta,tb,nm;
                    374:        struct order_spec *current_spec;
                    375:
                    376:        if ( !(nf=current_numberfield) )
                    377:                error("adddalg : current_numberfield is not set");
                    378:        if ( !a )
                    379:                *c = b;
                    380:        else if ( !b )
                    381:                *c = a;
                    382:        else {
1.4       noro      383:                qtodalg((Q)a,&t); a = t; qtodalg((Q)b,&t); b = t;
1.3       noro      384:                dna = a->dn;
                    385:                dnb = b->dn;
                    386:                gcdn(NM(dna),NM(dnb),&gn);
                    387:                divsn(NM(dna),gn,&an); divsn(NM(dnb),gn,&bn);
                    388:                NTOQ(an,SGN(dna),a1); NTOQ(bn,SGN(dnb),b1);
                    389:                /* nma/dna+nmb/dnb = (nma*b1+nmb*a1)/(dna*b1) */
                    390:                muldc(CO,a->nm,(P)b1,&ta); muldc(CO,b->nm,(P)a1,&tb);
                    391:                current_spec = dp_current_spec; initd(nf->spec);
                    392:                addd(CO,ta,tb,&nm);
                    393:                initd(current_spec);
                    394:                if ( !nm )
                    395:                        *c = 0;
                    396:                else {
                    397:                        mulq(dna,b1,&dn);
                    398:                        MKDAlg(nm,dn,*c);
                    399:                }
                    400:        }
1.1       noro      401: }
                    402:
                    403: void subdalg(DAlg a,DAlg b,DAlg *c)
                    404: {
1.3       noro      405:        NumberField nf;
                    406:        Q dna,dnb,a1,b1,dn,g;
                    407:        N an,bn,gn;
                    408:        DP ta,tb,nm;
1.4       noro      409:        DAlg t;
1.3       noro      410:        struct order_spec *current_spec;
                    411:
                    412:        if ( !(nf=current_numberfield) )
                    413:                error("subdalg : current_numberfield is not set");
                    414:        if ( !a )
                    415:                *c = b;
                    416:        else if ( !b )
                    417:                *c = a;
                    418:        else {
1.4       noro      419:                qtodalg((Q)a,&t); a = t; qtodalg((Q)b,&t); b = t;
1.3       noro      420:                dna = a->dn;
                    421:                dnb = b->dn;
                    422:                gcdn(NM(dna),NM(dnb),&gn);
                    423:                divsn(NM(dna),gn,&an); divsn(NM(dnb),gn,&bn);
                    424:                NTOQ(an,SGN(dna),a1); NTOQ(bn,SGN(dnb),b1);
                    425:                /* nma/dna-nmb/dnb = (nma*b1-nmb*a1)/(dna*b1) */
                    426:                muldc(CO,a->nm,(P)b1,&ta); muldc(CO,b->nm,(P)a1,&tb);
                    427:                current_spec = dp_current_spec; initd(nf->spec);
                    428:                subd(CO,ta,tb,&nm);
                    429:                initd(current_spec);
                    430:                if ( !nm )
                    431:                        *c = 0;
                    432:                else {
                    433:                        mulq(dna,b1,&dn);
                    434:                        MKDAlg(nm,dn,*c);
                    435:                }
                    436:        }
1.1       noro      437: }
                    438:
1.2       noro      439: void muldalg(DAlg a,DAlg b,DAlg *c)
                    440: {
1.3       noro      441:        NumberField nf;
                    442:        DP nm;
                    443:        Q dn;
                    444:        DAlg t;
                    445:        struct order_spec *current_spec;
                    446:
                    447:        if ( !(nf=current_numberfield) )
                    448:                error("muldalg : current_numberfield is not set");
                    449:        if ( !a || !b )
                    450:                *c = 0;
                    451:        else {
1.4       noro      452:                qtodalg((Q)a,&t); a = t; qtodalg((Q)b,&t); b = t;
1.3       noro      453:                current_spec = dp_current_spec; initd(nf->spec);
                    454:                muld(CO,a->nm,b->nm,&nm);
                    455:                initd(current_spec);
                    456:                mulq(a->dn,b->dn,&dn);
                    457:                MKDAlg(nm,dn,t);
                    458:                simpdalg(t,c);
                    459:        }
1.2       noro      460: }
                    461:
                    462:
1.1       noro      463: void divdalg(DAlg a,DAlg b,DAlg *c)
                    464: {
1.4       noro      465:        DAlg inv,t;
1.3       noro      466:
1.1       noro      467:        if ( !current_numberfield )
1.3       noro      468:                error("divdalg : current_numberfield is not set");
                    469:        if ( !b )
                    470:                error("divdalg : division by 0");
                    471:        if ( !a )
                    472:                c = 0;
                    473:        else {
1.4       noro      474:                qtodalg((Q)a,&t); a = t; qtodalg((Q)b,&t); b = t;
1.3       noro      475:                invdalg(b,&inv);
                    476:                muldalg(a,inv,c);
                    477:        }
                    478: }
                    479:
                    480: void rmcontdalg(DAlg a, DAlg *r)
                    481: {
                    482:        DP u,u1;
                    483:        Q cont,c,d;
                    484:        N gn,cn,dn;
                    485:
                    486:        if ( !a )
                    487:                *r = a;
                    488:        else {
                    489:                dp_ptozp(a->nm,&u);
                    490:                divq((Q)BDY(a->nm)->c,(Q)BDY(u)->c,&cont);
                    491:                gcdn(NM(cont),NM(a->dn),&gn);
                    492:                divsn(NM(cont),gn,&cn); NTOQ(cn,SGN(cont),c);
                    493:                divsn(NM(a->dn),gn,&dn); NTOQ(dn,SGN(a->dn),d);
                    494:                muldc(CO,u,(P)c,&u1);
                    495:                MKDAlg(u1,d,*r);
                    496:        }
1.1       noro      497: }
                    498:
                    499: void invdalg(DAlg a,DAlg *c)
                    500: {
1.3       noro      501:        NumberField nf;
                    502:        int dim,n,i,j;
                    503:        DP *mb;
                    504:        DP m,d,u;
                    505:        N ln,gn,qn;
                    506:        DAlg *simp;
                    507:        DAlg t,a0,r;
                    508:        Q dn,dnsol,mul;
                    509:        MAT mobj,sol;
                    510:        Q **mat,**solmat;
                    511:        MP mp0,mp;
                    512:        int *rinfo,*cinfo;
                    513:        struct order_spec *current_spec;
1.5       noro      514:        struct oEGT eg0,eg1;
                    515:        extern struct oEGT eg_le;
1.3       noro      516:
                    517:        if ( !(nf=current_numberfield) )
                    518:                error("invdalg : current_numberfield is not set");
                    519:        if ( !a )
                    520:                error("invdalg : division by 0");
1.4       noro      521:        else if ( NID(a) == N_Q ) {
                    522:                invq((Q)a,&dn); *c = (DAlg)dn;
                    523:                return;
                    524:        }
1.3       noro      525:        dim = nf->dim;
                    526:        mb = nf->mb;
                    527:        n = nf->n;
                    528:        ln = ONEN;
                    529:        MKDAlg(a->nm,ONE,a0);
                    530:        simp = (DAlg *)ALLOCA(dim*sizeof(DAlg));
                    531:        current_spec = dp_current_spec; initd(nf->spec);
                    532:        for ( i = dim-1; i >= 0; i-- ) {
                    533:                m = mb[i];
                    534:                for ( j = i+1; j < dim; j++ )
                    535:                        if ( dp_redble(m,mb[j]) )
                    536:                                break;
                    537:                if ( j < dim ) {
                    538:                        dp_subd(m,mb[j],&d);
                    539:                        muld(CO,d,simp[j]->nm,&u);
                    540:                        MKDAlg(u,simp[j]->dn,t);
                    541:                        simpdalg(t,&simp[i]);
                    542:                } else {
                    543:                        MKDAlg(m,ONE,t);
                    544:                        muldalg(t,a0,&simp[i]);
                    545:                }
                    546:                gcdn(NM(simp[i]->dn),ln,&gn); divsn(ln,gn,&qn);
                    547:                muln(NM(simp[i]->dn),qn,&ln);
                    548:        }
                    549:        initd(current_spec);
                    550:        NTOQ(ln,1,dn);
                    551:        MKMAT(mobj,dim,dim+1);
                    552:        mat = (Q **)BDY(mobj);
                    553:        mulq(dn,a->dn,&mat[dim-1][dim]);
                    554:        for ( j = 0; j < dim; j++ ) {
                    555:                divq(dn,simp[j]->dn,&mul);
                    556:                for ( i = 0, mp = BDY(simp[j]->nm); mp && i < dim; i++ )
                    557:                        if ( dl_equal(n,BDY(mb[i])->dl,mp->dl) ) {
                    558:                                mulq(mul,(Q)mp->c,&mat[i][j]);
                    559:                                mp = NEXT(mp);
                    560:                        }
                    561:        }
1.5       noro      562:        get_eg(&eg0);
1.3       noro      563:        generic_gauss_elim_hensel(mobj,&sol,&dnsol,&rinfo,&cinfo);
1.5       noro      564:        get_eg(&eg1); add_eg(&eg_le,&eg0,&eg1);
1.3       noro      565:        solmat = (Q **)BDY(sol);
                    566:        for ( i = 0, mp0 = 0; i < dim; i++ )
                    567:                if ( solmat[i][0] ) {
                    568:                        NEXTMP(mp0,mp);
                    569:                        mp->c = (P)solmat[i][0];
                    570:                        mp->dl = BDY(mb[i])->dl;
                    571:                }
                    572:        NEXT(mp) = 0; MKDP(n,mp0,u);
                    573:        MKDAlg(u,dnsol,r);
                    574:        rmcontdalg(r,c);
1.1       noro      575: }
                    576:
                    577: void chsgndalg(DAlg a,DAlg *c)
                    578: {
1.3       noro      579:        DP nm;
1.4       noro      580:        Q t;
1.3       noro      581:
                    582:        if ( !a ) *c = 0;
1.4       noro      583:        else if ( NID(a) == N_Q ) {
                    584:                chsgnq((Q)a,&t); *c = (DAlg)t;
                    585:        } else {
1.3       noro      586:                chsgnd(a->nm,&nm);
                    587:                MKDAlg(nm,a->dn,*c);
                    588:        }
1.1       noro      589: }
                    590:
1.3       noro      591: void pwrdalg(DAlg a,Q e,DAlg *c)
1.1       noro      592: {
1.3       noro      593:        NumberField nf;
                    594:        DAlg t,z,y;
1.4       noro      595:        Q q;
1.3       noro      596:        N en,qn;
                    597:        int r;
                    598:
                    599:        if ( !(nf=current_numberfield) )
                    600:                error("pwrdalg : current_numberfield is not set");
1.4       noro      601:        if ( !a )
                    602:                *c = !e ? (DAlg)ONE : 0;
                    603:        else if ( NID(a) == N_Q ) {
                    604:                pwrq((Q)a,e,&q); *c = (DAlg)q;
                    605:        } else if ( !e )
1.3       noro      606:                *c = nf->one;
                    607:        else if ( UNIQ(e) )
                    608:                *c = a;
                    609:        else {
                    610:                if ( SGN(e) < 0 ) {
                    611:                        invdalg(a,&t); a = t;
                    612:                }
                    613:                en = NM(e);
                    614:                y = nf->one;
                    615:                z = a;
                    616:                while ( 1 ) {
                    617:                        r = divin(en,2,&qn); en = qn;
                    618:                        if ( r ) {
                    619:                                muldalg(z,y,&t); y = t;
                    620:                                if ( !en ) {
                    621:                                        *c = y;
                    622:                                        return;
                    623:                                }
                    624:                        }
                    625:                        muldalg(z,z,&t); z = t;
                    626:                }
                    627:        }
1.1       noro      628: }
                    629:
1.3       noro      630: int cmpdalg(DAlg a,DAlg b)
1.1       noro      631: {
1.3       noro      632:        DAlg c;
                    633:
                    634:        subdalg(a,b,&c);
                    635:        if ( !c ) return 0;
                    636:        else
                    637:                return SGN((Q)BDY(c->nm)->c);
1.1       noro      638: }

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