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

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

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
        !            26:  * e-mail at risa-admin@flab.fujitsu.co.jp of the detailed specification
        !            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:  *
        !            48:  * $OpenXM: OpenXM_contrib2/asir2000/engine/_distm.c,v 1.1 2000/07/13 05:50:39 noro Exp $
        !            49: */
1.1       noro       50: #include "ca.h"
                     51: #include "inline.h"
                     52:
                     53: #define NV(p) ((p)->nv)
                     54: #define C(p) ((p)->c)
                     55: #if 0
                     56: #define ITOS(p) (((unsigned int)(p))>>1)
                     57: #define STOI(i) ((P)((((unsigned int)(i))<<1)|1))
                     58: #else
                     59: #define ITOS(p) (((unsigned int)(p))&0x7fffffff)
                     60: #define STOI(i) ((P)((unsigned int)(i)|0x80000000))
                     61: #endif
                     62:
                     63: extern int (*cmpdl)();
                     64: extern int do_weyl;
                     65:
                     66: void _dptodp();
                     67: void adddl_dup();
                     68: void _mulmdm_dup();
                     69: void _free_dp();
                     70: void _comm_mulmd_dup();
                     71: void _weyl_mulmd_dup();
                     72: void _weyl_mulmmm_dup();
                     73: void _weyl_mulmdm_dup();
                     74:
                     75: MP _mp_free_list;
                     76: DP _dp_free_list;
                     77: DL _dl_free_list;
                     78: int current_dl_length;
                     79:
                     80: #define _NEWDL_NOINIT(d,n) if ((n)!= current_dl_length){_dl_free_list=0; current_dl_length=(n);} if(!_dl_free_list)_DL_alloc(); (d)=_dl_free_list; _dl_free_list = *((DL *)_dl_free_list)
                     81: #define _NEWDL(d,n) if ((n)!= current_dl_length){_dl_free_list=0; current_dl_length=(n);} if(!_dl_free_list)_DL_alloc(); (d)=_dl_free_list; _dl_free_list = *((DL *)_dl_free_list); bzero((d),(((n)+1)*sizeof(int)))
                     82: #define _NEWMP(m) if(!_mp_free_list)_MP_alloc(); (m)=_mp_free_list; _mp_free_list = NEXT(_mp_free_list)
                     83: #define _MKDP(n,m,d) if(!_dp_free_list)_DP_alloc(); (d)=_dp_free_list; _dp_free_list = (DP)BDY(_dp_free_list); (d)->nv=(n); BDY(d)=(m)
                     84:
                     85: #define _NEXTMP(r,c) \
                     86: if(!(r)){_NEWMP(r);(c)=(r);}else{_NEWMP(NEXT(c));(c)=NEXT(c);}
                     87:
                     88: #define _NEXTMP2(r,c,s) \
                     89: if(!(r)){(c)=(r)=(s);}else{NEXT(c)=(s);(c)=(s);}
                     90:
                     91: #define FREEDL(m) *((DL *)m)=_dl_free_list; _dl_free_list=(m)
                     92: #define FREEMP(m) NEXT(m)=_mp_free_list; _mp_free_list=(m)
                     93: #define FREEDP(m) BDY(m)=(MP)_dp_free_list; _dp_free_list=(m)
                     94:
                     95:
                     96: void _DL_alloc()
                     97: {
                     98:        int *p;
                     99:        int i,dl_len;
                    100:        static int DL_alloc_count;
                    101:
                    102: /*     fprintf(stderr,"DL_alloc : %d \n",++DL_alloc_count); */
                    103:        dl_len = (current_dl_length+1);
                    104:        p = (int *)GC_malloc(128*dl_len*sizeof(int));
                    105:        for ( i = 0; i < 128; i++, p += dl_len ) {
                    106:                *(DL *)p = _dl_free_list;
                    107:                _dl_free_list = (DL)p;
                    108:        }
                    109: }
                    110:
                    111: void _MP_alloc()
                    112: {
                    113:        MP p;
                    114:        int i;
                    115:        static int MP_alloc_count;
                    116:
                    117: /*     fprintf(stderr,"MP_alloc : %d \n",++MP_alloc_count); */
                    118:        p = (MP)GC_malloc(1024*sizeof(struct oMP));
                    119:        for ( i = 0; i < 1024; i++ ) {
                    120:                p[i].next = _mp_free_list; _mp_free_list = &p[i];
                    121:        }
                    122: }
                    123:
                    124: void _DP_alloc()
                    125: {
                    126:        DP p;
                    127:        int i;
                    128:        static int DP_alloc_count;
                    129:
                    130: /*     fprintf(stderr,"DP_alloc : %d \n",++DP_alloc_count); */
                    131:        p = (DP)GC_malloc(1024*sizeof(struct oDP));
                    132:        for ( i = 0; i < 1024; i++ ) {
                    133:                p[i].body = (MP)_dp_free_list; _dp_free_list = &p[i];
                    134:        }
                    135: }
                    136:
                    137: /* merge p1 and p2 into pr */
                    138:
                    139: void _addmd_destructive(mod,p1,p2,pr)
                    140: int mod;
                    141: DP p1,p2,*pr;
                    142: {
                    143:        int n;
                    144:        MP m1,m2,mr,mr0,s;
                    145:        int t;
                    146:
                    147:        if ( !p1 )
                    148:                *pr = p2;
                    149:        else if ( !p2 )
                    150:                *pr = p1;
                    151:        else {
                    152:                for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; )
                    153:                        switch ( (*cmpdl)(n,m1->dl,m2->dl) ) {
                    154:                                case 0:
                    155:                                        t = (ITOS(C(m1))+ITOS(C(m2))) - mod;
                    156:                                        if ( t < 0 )
                    157:                                                t += mod;
                    158:                                        s = m1; m1 = NEXT(m1);
                    159:                                        if ( t ) {
                    160:                                                _NEXTMP2(mr0,mr,s); C(mr) = STOI(t);
                    161:                                        } else {
                    162:                                                FREEDL(s->dl); FREEMP(s);
                    163:                                        }
                    164:                                        s = m2; m2 = NEXT(m2); FREEDL(s->dl); FREEMP(s);
                    165:                                        break;
                    166:                                case 1:
                    167:                                        s = m1; m1 = NEXT(m1); _NEXTMP2(mr0,mr,s);
                    168:                                        break;
                    169:                                case -1:
                    170:                                        s = m2; m2 = NEXT(m2); _NEXTMP2(mr0,mr,s);
                    171:                                        break;
                    172:                        }
                    173:                if ( !mr0 )
                    174:                        if ( m1 )
                    175:                                mr0 = m1;
                    176:                        else if ( m2 )
                    177:                                mr0 = m2;
                    178:                        else {
                    179:                                *pr = 0;
                    180:                                return;
                    181:                        }
                    182:                else if ( m1 )
                    183:                        NEXT(mr) = m1;
                    184:                else if ( m2 )
                    185:                        NEXT(mr) = m2;
                    186:                else
                    187:                        NEXT(mr) = 0;
                    188:                _MKDP(NV(p1),mr0,*pr);
                    189:                if ( *pr )
                    190:                        (*pr)->sugar = MAX(p1->sugar,p2->sugar);
                    191:                FREEDP(p1); FREEDP(p2);
                    192:        }
                    193: }
                    194:
                    195: void _mulmd_dup(mod,p1,p2,pr)
                    196: int mod;
                    197: DP p1,p2,*pr;
                    198: {
                    199:        if ( !do_weyl )
                    200:                _comm_mulmd_dup(mod,p1,p2,pr);
                    201:        else
                    202:                _weyl_mulmd_dup(mod,p1,p2,pr);
                    203: }
                    204:
                    205: void _comm_mulmd_dup(mod,p1,p2,pr)
                    206: int mod;
                    207: DP p1,p2,*pr;
                    208: {
                    209:        MP m;
                    210:        DP s,t,u;
                    211:        int i,l,l1;
                    212:        static MP *w;
                    213:        static int wlen;
                    214:
                    215:        if ( !p1 || !p2 )
                    216:                *pr = 0;
                    217:        else {
                    218:                for ( m = BDY(p1), l1 = 0; m; m = NEXT(m), l1++ );
                    219:                for ( m = BDY(p2), l = 0; m; m = NEXT(m), l++ );
                    220:                if ( l1 < l ) {
                    221:                        t = p1; p1 = p2; p2 = t;
                    222:                        l = l1;
                    223:                }
                    224:                if ( l > wlen ) {
                    225:                        if ( w ) GC_free(w);
                    226:                        w = (MP *)MALLOC(l*sizeof(MP));
                    227:                        wlen = l;
                    228:                }
                    229:                for ( m = BDY(p2), i = 0; i < l; m = NEXT(m), i++ )
                    230:                        w[i] = m;
                    231:                for ( s = 0, i = l-1; i >= 0; i-- ) {
                    232:                        _mulmdm_dup(mod,p1,w[i],&t); _addmd_destructive(mod,s,t,&u); s = u;
                    233:                }
                    234:                bzero(w,l*sizeof(MP));
                    235:                *pr = s;
                    236:        }
                    237: }
                    238:
                    239: void _weyl_mulmd_dup(mod,p1,p2,pr)
                    240: int mod;
                    241: DP p1,p2,*pr;
                    242: {
                    243:        MP m;
                    244:        DP s,t,u;
                    245:        int i,l,l1;
                    246:        static MP *w;
                    247:        static int wlen;
                    248:
                    249:        if ( !p1 || !p2 )
                    250:                *pr = 0;
                    251:        else {
                    252:                for ( m = BDY(p2), l = 0; m; m = NEXT(m), l++ );
                    253:                if ( l > wlen ) {
                    254:                        if ( w ) GC_free(w);
                    255:                        w = (MP *)MALLOC(l*sizeof(MP));
                    256:                        wlen = l;
                    257:                }
                    258:                for ( m = BDY(p2), i = 0; i < l; m = NEXT(m), i++ )
                    259:                        w[i] = m;
                    260:                for ( s = 0, i = l-1; i >= 0; i-- ) {
                    261:                        _weyl_mulmdm_dup(mod,p1,w[i],&t); _addmd_destructive(mod,s,t,&u); s = u;
                    262:                }
                    263:                bzero(w,l*sizeof(MP));
                    264:                *pr = s;
                    265:        }
                    266: }
                    267:
                    268: void _mulmdm_dup(mod,p,m0,pr)
                    269: int mod;
                    270: DP p;
                    271: MP m0;
                    272: DP *pr;
                    273: {
                    274:        MP m,mr,mr0;
                    275:        DL d,dt,dm;
                    276:        int c,n,r,i;
                    277:        int *pt,*p1,*p2;
                    278:
                    279:        if ( !p )
                    280:                *pr = 0;
                    281:        else {
                    282:                for ( mr0 = 0, m = BDY(p), c = ITOS(C(m0)), d = m0->dl, n = NV(p);
                    283:                        m; m = NEXT(m) ) {
                    284:                        _NEXTMP(mr0,mr);
                    285:                        C(mr) = STOI(dmar(ITOS(C(m)),c,0,mod));
                    286:                        _NEWDL_NOINIT(dt,n); mr->dl = dt;
                    287:                        dm = m->dl;
                    288:                        dt->td = d->td + dm->td;
                    289:                        for ( i = 0, pt = dt->d, p1=d->d, p2 = dm->d; i < n; i++ )
                    290:                                *pt++ = *p1++ + *p2++;
                    291:                }
                    292:                NEXT(mr) = 0; _MKDP(NV(p),mr0,*pr);
                    293:                if ( *pr )
                    294:                        (*pr)->sugar = p->sugar + m0->dl->td;
                    295:        }
                    296: }
                    297:
                    298: void _weyl_mulmdm_dup(mod,p,m0,pr)
                    299: int mod;
                    300: DP p;
                    301: MP m0;
                    302: DP *pr;
                    303: {
                    304:        DP r,t,t1;
                    305:        MP m;
                    306:        int n,l,i;
                    307:        static MP *w;
                    308:        static int wlen;
                    309:
                    310:        if ( !p )
                    311:                *pr = 0;
                    312:        else {
                    313:                for ( m = BDY(p), l = 0; m; m = NEXT(m), l++ );
                    314:                if ( l > wlen ) {
                    315:                        if ( w ) GC_free(w);
                    316:                        w = (MP *)MALLOC(l*sizeof(MP));
                    317:                        wlen = l;
                    318:                }
                    319:                for ( m = BDY(p), i = 0; i < l; m = NEXT(m), i++ )
                    320:                        w[i] = m;
                    321:                for ( r = 0, i = l-1, n = NV(p); i >= 0; i-- ) {
                    322:                        _weyl_mulmmm_dup(mod,w[i],m0,n,&t);
                    323:                        _addmd_destructive(mod,r,t,&t1); r = t1;
                    324:                }
                    325:                bzero(w,l*sizeof(MP));
                    326:                if ( r )
                    327:                        r->sugar = p->sugar + m0->dl->td;
                    328:                *pr = r;
                    329:        }
                    330: }
                    331: /* m0 = x0^d0*x1^d1*... * dx0^d(n/2)*dx1^d(n/2+1)*... */
                    332:
                    333: void _weyl_mulmmm_dup(mod,m0,m1,n,pr)
                    334: int mod;
                    335: MP m0,m1;
                    336: int n;
                    337: DP *pr;
                    338: {
                    339:        MP m,mr,mr0;
                    340:        DP r,t,t1;
                    341:        int c,c0,c1,cc;
                    342:        DL d,d0,d1;
                    343:        int i,j,a,b,k,l,n2,s,min,h;
                    344:        static int *tab;
                    345:        static int tablen;
                    346:
                    347:        if ( !m0 || !m1 )
                    348:                *pr = 0;
                    349:        else {
                    350:                c0 = ITOS(C(m0)); c1 = ITOS(C(m1));
                    351:                c = dmar(c0,c1,0,mod);
                    352:                d0 = m0->dl; d1 = m1->dl;
                    353:                n2 = n>>1;
                    354:
                    355:                _NEWDL(d,n);
                    356:                if ( n & 1 )
                    357:                        /* offset of h-degree */
                    358:                        d->td = d->d[n-1] = d0->d[n-1]+d1->d[n-1];
                    359:                else
                    360:                        d->td = 0;
                    361:                _NEWMP(mr); mr->c = STOI(c); mr->dl = d; NEXT(mr) = 0;
                    362:                _MKDP(n,mr,r); r->sugar = d->td;
                    363:
                    364:                /* homogenized computation; dx-xd=h^2 */
                    365:                for ( i = 0; i < n2; i++ ) {
                    366:                        a = d0->d[i]; b = d1->d[n2+i];
                    367:                        k = d0->d[n2+i]; l = d1->d[i];
                    368:                        /* degree of xi^a*(Di^k*xi^l)*Di^b */
                    369:                        s = a+k+l+b;
                    370:                        /* compute xi^a*(Di^k*xi^l)*Di^b */
                    371:                        min = MIN(k,l);
                    372:
                    373:                        if ( min+1 > tablen ) {
                    374:                                if ( tab ) GC_free(tab);
                    375:                                tab = (int *)MALLOC((min+1)*sizeof(int));
                    376:                                tablen = min+1;
                    377:                        }
                    378:                        mkwcm(k,l,mod,tab);
                    379:                        if ( n & 1 )
                    380:                                for ( mr0 = 0, j = 0; j <= min; j++ ) {
                    381:                                        _NEXTMP(mr0,mr); _NEWDL(d,n);
                    382:                                        d->d[i] = l-j+a; d->d[n2+i] = k-j+b;
                    383:                                        d->td = s;
                    384:                                        d->d[n-1] = s-(d->d[i]+d->d[n2+i]);
                    385:                                        mr->c = STOI(tab[j]); mr->dl = d;
                    386:                                }
                    387:                        else
                    388:                                for ( mr0 = 0, s = 0, j = 0; j <= min; j++ ) {
                    389:                                        _NEXTMP(mr0,mr); _NEWDL(d,n);
                    390:                                        d->d[i] = l-j+a; d->d[n2+i] = k-j+b;
                    391:                                        d->td = d->d[i]+d->d[n2+i]; /* XXX */
                    392:                                        s = MAX(s,d->td); /* XXX */
                    393:                                        mr->c = STOI(tab[j]); mr->dl = d;
                    394:                                }
                    395:                        bzero(tab,(min+1)*sizeof(int));
                    396:                        if ( mr0 )
                    397:                                NEXT(mr) = 0;
                    398:                        _MKDP(n,mr0,t);
                    399:                        if ( t )
                    400:                                t->sugar = s;
                    401:                        _comm_mulmd_dup(mod,r,t,&t1);
                    402:                        _free_dp(r); _free_dp(t);
                    403:                        r = t1;
                    404:                }
                    405:                *pr = r;
                    406:        }
                    407: }
                    408:
                    409: void _dp_red_mod_destructive(p1,p2,mod,rp)
                    410: DP p1,p2;
                    411: int mod;
                    412: DP *rp;
                    413: {
                    414:        int i,n;
                    415:        DL d1,d2,d;
                    416:        MP m;
                    417:        DP t,s;
                    418:        int c,c1;
                    419:
                    420:        n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
                    421:        _NEWDL(d,n); d->td = d1->td - d2->td;
                    422:        for ( i = 0; i < n; i++ )
                    423:                d->d[i] = d1->d[i]-d2->d[i];
                    424:        c = invm(ITOS(BDY(p2)->c),mod); c1 = dmar(c,ITOS(BDY(p1)->c),0,mod);
                    425:        _NEWMP(m); m->dl = d; m->c = STOI(mod-c1); NEXT(m) = 0;
                    426:        _MKDP(n,m,s); s->sugar = d->td;
                    427:        _mulmd_dup(mod,s,p2,&t); _free_dp(s);
                    428:        _addmd_destructive(mod,p1,t,rp);
                    429: }
                    430:
                    431: void _dp_sp_mod_dup(p1,p2,mod,rp)
                    432: DP p1,p2;
                    433: int mod;
                    434: DP *rp;
                    435: {
                    436:        int i,n,td;
                    437:        int *w;
                    438:        DL d1,d2,d;
                    439:        MP m;
                    440:        DP t,s,u;
                    441:
                    442:        n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
                    443:        w = (int *)ALLOCA(n*sizeof(int));
                    444:        for ( i = 0, td = 0; i < n; i++ ) {
                    445:                w[i] = MAX(d1->d[i],d2->d[i]); td += w[i];
                    446:        }
                    447:        _NEWDL(d,n); d->td = td - d1->td;
                    448:        for ( i = 0; i < n; i++ )
                    449:                d->d[i] = w[i] - d1->d[i];
                    450:        _NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
                    451:        _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p1,&t); _free_dp(s);
                    452:        _NEWDL(d,n); d->td = td - d2->td;
                    453:        for ( i = 0; i < n; i++ )
                    454:                d->d[i] = w[i] - d2->d[i];
                    455:        _NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
                    456:        _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p2,&u); _free_dp(s);
                    457:        _addmd_destructive(mod,t,u,rp);
                    458: }
                    459:
                    460: void _dltodl(d,dr)
                    461: DL d;
                    462: DL *dr;
                    463: {
                    464:        int i,n;
                    465:        DL t;
                    466:
                    467:        n = current_dl_length;
                    468:        NEWDL(t,n); *dr = t;
                    469:        t->td = d->td;
                    470:        for ( i = 0; i < n; i++ )
                    471:                t->d[i] = d->d[i];
                    472: }
                    473:
                    474: void adddl_dup(n,d1,d2,dr)
                    475: int n;
                    476: DL d1,d2;
                    477: DL *dr;
                    478: {
                    479:        DL dt;
                    480:        int i;
                    481:
                    482:        _NEWDL(dt,n); *dr = dt;
                    483:        dt->td = d1->td + d2->td;
                    484:        for ( i = 0; i < n; i++ )
                    485:                dt->d[i] = d1->d[i]+d2->d[i];
                    486: }
                    487:
                    488: void _free_dp(f)
                    489: DP f;
                    490: {
                    491:        MP m,s;
                    492:
                    493:        if ( !f )
                    494:                return;
                    495:        m = f->body;
                    496:        while ( m ) {
                    497:                s = m; m = NEXT(m); FREEDL(s->dl); FREEMP(s);
                    498:        }
                    499:        FREEDP(f);
                    500: }
                    501:
                    502: void _dptodp(p,r)
                    503: DP p;
                    504: DP *r;
                    505: {
                    506:        MP m,mr0,mr;
                    507:
                    508:        if ( !p )
                    509:                *r = 0;
                    510:        else {
                    511:                for ( m = BDY(p), mr0 = 0; m; m = NEXT(m) ) {
                    512:                        NEXTMP(mr0,mr);
                    513:                        _dltodl(m->dl,&mr->dl);
                    514:                        mr->c = m->c;
                    515:                }
                    516:                NEXT(mr) = 0;
                    517:                MKDP(p->nv,mr0,*r);
                    518:                (*r)->sugar = p->sugar;
                    519:        }
                    520: }
                    521:
                    522: void _dp_nf_mod_destructive(b,g,ps,mod,full,rp)
                    523: NODE b;
                    524: DP g;
                    525: DP *ps;
                    526: int mod,full;
                    527: DP *rp;
                    528: {
                    529:        DP u,p,d,s,t;
                    530:        NODE l;
                    531:        MP m,mr,mrd;
                    532:        int sugar,psugar,n,h_reducible,i;
                    533:
                    534:        if ( !g ) {
                    535:                *rp = 0; return;
                    536:        }
                    537:        sugar = g->sugar;
                    538:        n = g->nv;
                    539:        for ( d = 0; g; ) {
                    540:                for ( h_reducible = 0, l = b; l; l = NEXT(l) ) {
                    541:                        if ( dp_redble(g,p = ps[(int)BDY(l)]) ) {
                    542:                                h_reducible = 1;
                    543:                                psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
                    544:                                _dp_red_mod_destructive(g,p,mod,&u); g = u;
                    545:                                sugar = MAX(sugar,psugar);
                    546:                                if ( !g ) {
                    547:                                        if ( d )
                    548:                                                d->sugar = sugar;
                    549:                                        _dptodp(d,rp); _free_dp(d); return;
                    550:                                }
                    551:                                break;
                    552:                        }
                    553:                }
                    554:                if ( !h_reducible ) {
                    555:                        /* head term is not reducible */
                    556:                        if ( !full ) {
                    557:                                if ( g )
                    558:                                        g->sugar = sugar;
                    559:                                _dptodp(g,rp); _free_dp(g); return;
                    560:                        } else {
                    561:                                m = BDY(g);
                    562:                                if ( NEXT(m) ) {
                    563:                                        BDY(g) = NEXT(m); NEXT(m) = 0;
                    564:                                } else {
                    565:                                        FREEDP(g); g = 0;
                    566:                                }
                    567:                                if ( d ) {
                    568:                                        for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) );
                    569:                                        NEXT(mrd) = m;
                    570:                                } else {
                    571:                                        _MKDP(n,m,d);
                    572:                                }
                    573:                        }
                    574:                }
                    575:        }
                    576:        if ( d )
                    577:                d->sugar = sugar;
                    578:        _dptodp(d,rp); _free_dp(d);
                    579: }

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