Annotation of OpenXM_contrib2/asir2018/builtin/dp-supp.c, Revision 1.1
1.1 ! 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@sec.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$
! 49: */
! 50: #include "ca.h"
! 51: #include "base.h"
! 52: #include "inline.h"
! 53: #include "parse.h"
! 54: #include "ox.h"
! 55:
! 56: #define HMAG(p) (p_mag((P)BDY(p)->c))
! 57:
! 58: extern int (*cmpdl)();
! 59: extern double pz_t_e,pz_t_d,pz_t_d1,pz_t_c;
! 60: extern int dp_nelim,dp_fcoeffs;
! 61: extern int NoGCD;
! 62: extern int GenTrace;
! 63: extern NODE TraceList;
! 64:
! 65: int show_orderspec;
! 66:
! 67: void print_composite_order_spec(struct order_spec *spec);
! 68: void dpm_rest(DPM,DPM *);
! 69:
! 70: /*
! 71: * content reduction
! 72: *
! 73: */
! 74:
! 75: static NODE RatDenomList;
! 76:
! 77: void init_denomlist()
! 78: {
! 79: RatDenomList = 0;
! 80: }
! 81:
! 82: void add_denomlist(P f)
! 83: {
! 84: NODE n;
! 85:
! 86: if ( OID(f)==O_P ) {
! 87: MKNODE(n,f,RatDenomList); RatDenomList = n;
! 88: }
! 89: }
! 90:
! 91: LIST get_denomlist()
! 92: {
! 93: LIST l;
! 94:
! 95: MKLIST(l,RatDenomList); RatDenomList = 0;
! 96: return l;
! 97: }
! 98:
! 99: void dp_ptozp(DP p,DP *rp)
! 100: {
! 101: MP m,mr,mr0;
! 102: int i,n;
! 103: Q *w;
! 104: Z dvr;
! 105: P t;
! 106:
! 107: if ( !p )
! 108: *rp = 0;
! 109: else {
! 110: for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ );
! 111: w = (Q *)ALLOCA(n*sizeof(Q));
! 112: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
! 113: if ( NUM(m->c) )
! 114: w[i] = (Q)m->c;
! 115: else
! 116: ptozp((P)m->c,1,&w[i],&t);
! 117: sortbynm(w,n);
! 118: qltozl(w,n,&dvr);
! 119: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 120: NEXTMP(mr0,mr); divsp(CO,(P)m->c,(P)dvr,(P *)&mr->c); mr->dl = m->dl;
! 121: }
! 122: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 123: }
! 124: }
! 125:
! 126: void dp_ptozp2(DP p0,DP p1,DP *hp,DP *rp)
! 127: {
! 128: DP t,s,h,r;
! 129: MP m,mr,mr0,m0;
! 130:
! 131: addd(CO,p0,p1,&t); dp_ptozp(t,&s);
! 132: if ( !p0 ) {
! 133: h = 0; r = s;
! 134: } else if ( !p1 ) {
! 135: h = s; r = 0;
! 136: } else {
! 137: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
! 138: m = NEXT(m), m0 = NEXT(m0) ) {
! 139: NEXTMP(mr0,mr); mr->c = m->c; mr->dl = m->dl;
! 140: }
! 141: NEXT(mr) = 0; MKDP(p0->nv,mr0,h); MKDP(p0->nv,m,r);
! 142: }
! 143: if ( h )
! 144: h->sugar = p0->sugar;
! 145: if ( r )
! 146: r->sugar = p1->sugar;
! 147: *hp = h; *rp = r;
! 148: }
! 149:
! 150: void dpm_ptozp(DPM p,DPM *rp)
! 151: {
! 152: DMM m,mr,mr0;
! 153: int i,n;
! 154: Q *w;
! 155: Z dvr;
! 156: P t;
! 157:
! 158: if ( !p )
! 159: *rp = 0;
! 160: else {
! 161: for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ );
! 162: w = (Q *)ALLOCA(n*sizeof(Q));
! 163: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
! 164: if ( NUM(m->c) )
! 165: w[i] = (Q)m->c;
! 166: else
! 167: ptozp((P)m->c,1,&w[i],&t);
! 168: sortbynm(w,n);
! 169: qltozl(w,n,&dvr);
! 170: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 171: NEXTDMM(mr0,mr); divsp(CO,(P)m->c,(P)dvr,(P *)&mr->c); mr->dl = m->dl; mr->pos = m->pos;
! 172: }
! 173: NEXT(mr) = 0; MKDPM(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 174: }
! 175: }
! 176:
! 177: void dpm_ptozp2(DPM p0,DPM p1,DPM *hp,DPM *rp)
! 178: {
! 179: DPM t,s,h,r;
! 180: DMM m,mr,mr0,m0;
! 181:
! 182: adddpm(CO,p0,p1,&t); dpm_ptozp(t,&s);
! 183: if ( !p0 ) {
! 184: h = 0; r = s;
! 185: } else if ( !p1 ) {
! 186: h = s; r = 0;
! 187: } else {
! 188: for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0;
! 189: m = NEXT(m), m0 = NEXT(m0) ) {
! 190: NEXTDMM(mr0,mr); mr->c = m->c; mr->dl = m->dl; mr->pos = m->pos;
! 191: }
! 192: NEXT(mr) = 0; MKDPM(p0->nv,mr0,h); MKDPM(p0->nv,m,r);
! 193: }
! 194: if ( h )
! 195: h->sugar = p0->sugar;
! 196: if ( r )
! 197: r->sugar = p1->sugar;
! 198: *hp = h; *rp = r;
! 199: }
! 200:
! 201:
! 202: void dp_ptozp3(DP p,Z *dvr,DP *rp)
! 203: {
! 204: MP m,mr,mr0;
! 205: int i,n;
! 206: Q *w;
! 207: P t;
! 208:
! 209: if ( !p ) {
! 210: *rp = 0; *dvr = 0;
! 211: }else {
! 212: for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ );
! 213: w = (Q *)ALLOCA(n*sizeof(Q));
! 214: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
! 215: if ( NUM(m->c) )
! 216: w[i] = (Q)m->c;
! 217: else
! 218: ptozp((P)m->c,1,&w[i],&t);
! 219: sortbynm(w,n);
! 220: qltozl(w,n,dvr);
! 221: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 222: NEXTMP(mr0,mr); divsp(CO,(P)m->c,(P)(*dvr),(P *)&mr->c); mr->dl = m->dl;
! 223: }
! 224: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 225: }
! 226: }
! 227:
! 228: void dp_idiv(DP p,Z c,DP *rp)
! 229: {
! 230: MP mr0,m,mr;
! 231:
! 232: if ( !p )
! 233: *rp = 0;
! 234: else if ( MUNIQ((Q)c) )
! 235: *rp = p;
! 236: else if ( MUNIQ((Q)c) )
! 237: chsgnd(p,rp);
! 238: else {
! 239: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 240: NEXTMP(mr0,mr);
! 241: divz((Z)(m->c),c,(Z *)&mr->c);
! 242: mr->dl = m->dl;
! 243: }
! 244: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
! 245: if ( *rp )
! 246: (*rp)->sugar = p->sugar;
! 247: }
! 248: }
! 249:
! 250: void dp_mbase(NODE hlist,NODE *mbase)
! 251: {
! 252: DL *dl;
! 253: DL d;
! 254: int *t;
! 255: int i,j,k,n,nvar,td;
! 256:
! 257: n = length(hlist); nvar = ((DP)BDY(hlist))->nv;
! 258: dl = (DL *)MALLOC(n*sizeof(DL));
! 259: NEWDL(d,nvar); *mbase = 0;
! 260: for ( i = 0; i < n; i++, hlist = NEXT(hlist) ) {
! 261: dl[i] = BDY((DP)BDY(hlist))->dl;
! 262: /* trivial ideal check */
! 263: if ( (*cmpdl)(nvar,d,dl[i]) == 0 ) {
! 264: return;
! 265: }
! 266: }
! 267: /* zero-dim. ideal check */
! 268: for ( i = 0; i < nvar; i++ ) {
! 269: for ( j = 0; j < n; j++ ) {
! 270: for ( k = 0, t = dl[j]->d; k < nvar; k++ )
! 271: if ( k != i && t[k] != 0 ) break;
! 272: if ( k == nvar ) break;
! 273: }
! 274: if ( j == n )
! 275: error("dp_mbase : input ideal is not zero-dimensional");
! 276: }
! 277: while ( 1 ) {
! 278: insert_to_node(d,mbase,nvar);
! 279: for ( i = nvar-1; i >= 0; ) {
! 280: d->d[i]++;
! 281: d->td += MUL_WEIGHT(1,i);
! 282: for ( j = 0; j < n; j++ ) {
! 283: if ( _dl_redble(dl[j],d,nvar) )
! 284: break;
! 285: }
! 286: if ( j < n ) {
! 287: for ( j = nvar-1; j >= i; j-- )
! 288: d->d[j] = 0;
! 289: for ( j = 0, td = 0; j < i; j++ )
! 290: td += MUL_WEIGHT(d->d[j],j);
! 291: d->td = td;
! 292: i--;
! 293: } else
! 294: break;
! 295: }
! 296: if ( i < 0 )
! 297: break;
! 298: }
! 299: }
! 300:
! 301: int _dl_redble(DL d1,DL d2,int nvar)
! 302: {
! 303: int i;
! 304:
! 305: if ( d1->td > d2->td )
! 306: return 0;
! 307: for ( i = 0; i < nvar; i++ )
! 308: if ( d1->d[i] > d2->d[i] )
! 309: break;
! 310: if ( i < nvar )
! 311: return 0;
! 312: else
! 313: return 1;
! 314: }
! 315:
! 316: void insert_to_node(DL d,NODE *n,int nvar)
! 317: {
! 318: DL d1;
! 319: MP m;
! 320: DP dp;
! 321: NODE n0,n1,n2;
! 322:
! 323: NEWDL(d1,nvar); d1->td = d->td;
! 324: bcopy((char *)d->d,(char *)d1->d,nvar*sizeof(int));
! 325: NEWMP(m); m->dl = d1; m->c = (Obj)ONE; NEXT(m) = 0;
! 326: MKDP(nvar,m,dp); dp->sugar = d->td;
! 327: if ( !(*n) ) {
! 328: MKNODE(n1,dp,0); *n = n1;
! 329: } else {
! 330: for ( n1 = *n, n0 = 0; n1; n0 = n1, n1 = NEXT(n1) )
! 331: if ( (*cmpdl)(nvar,d,BDY((DP)BDY(n1))->dl) > 0 ) {
! 332: MKNODE(n2,dp,n1);
! 333: if ( !n0 )
! 334: *n = n2;
! 335: else
! 336: NEXT(n0) = n2;
! 337: break;
! 338: }
! 339: if ( !n1 ) {
! 340: MKNODE(n2,dp,0); NEXT(n0) = n2;
! 341: }
! 342: }
! 343: }
! 344:
! 345: void dp_vtod(Q *c,DP p,DP *rp)
! 346: {
! 347: MP mr0,m,mr;
! 348: int i;
! 349:
! 350: if ( !p )
! 351: *rp = 0;
! 352: else {
! 353: for ( mr0 = 0, m = BDY(p), i = 0; m; m = NEXT(m), i++ ) {
! 354: NEXTMP(mr0,mr); mr->c = (Obj)c[i]; mr->dl = m->dl;
! 355: }
! 356: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp);
! 357: (*rp)->sugar = p->sugar;
! 358: }
! 359: }
! 360:
! 361: int have_sf_coef(P p)
! 362: {
! 363: DCP dc;
! 364:
! 365: if ( !p )
! 366: return 0;
! 367: else if ( NUM(p) )
! 368: return NID((Num)p) == N_GFS ? 1 : 0;
! 369: else {
! 370: for ( dc = DC(p); dc; dc = NEXT(dc) )
! 371: if ( have_sf_coef(COEF(dc)) )
! 372: return 1;
! 373: return 0;
! 374: }
! 375: }
! 376:
! 377: void head_coef(P p,Num *c)
! 378: {
! 379: if ( !p )
! 380: *c = 0;
! 381: else if ( NUM(p) )
! 382: *c = (Num)p;
! 383: else
! 384: head_coef(COEF(DC(p)),c);
! 385: }
! 386:
! 387: void dp_monic_sf(DP p,DP *rp)
! 388: {
! 389: Num c;
! 390:
! 391: if ( !p )
! 392: *rp = 0;
! 393: else {
! 394: head_coef((P)BDY(p)->c,&c);
! 395: divsdc(CO,p,(P)c,rp);
! 396: }
! 397: }
! 398:
! 399: void dp_prim(DP p,DP *rp)
! 400: {
! 401: P t,g;
! 402: DP p1;
! 403: MP m,mr,mr0;
! 404: int i,n;
! 405: P *w;
! 406: Q *c;
! 407: Z dvr;
! 408: NODE tn;
! 409:
! 410: if ( !p )
! 411: *rp = 0;
! 412: else if ( dp_fcoeffs == N_GFS ) {
! 413: for ( m = BDY(p); m; m = NEXT(m) )
! 414: if ( OID(m->c) == O_N ) {
! 415: /* GCD of coeffs = 1 */
! 416: dp_monic_sf(p,rp);
! 417: return;
! 418: } else break;
! 419: /* compute GCD over the finite fieid */
! 420: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
! 421: w = (P *)ALLOCA(n*sizeof(P));
! 422: for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ )
! 423: w[i] = (P)m->c;
! 424: gcdsf(CO,w,n,&g);
! 425: if ( NUM(g) )
! 426: dp_monic_sf(p,rp);
! 427: else {
! 428: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 429: NEXTMP(mr0,mr); divsp(CO,(P)m->c,g,(P *)&mr->c); mr->dl = m->dl;
! 430: }
! 431: NEXT(mr) = 0; MKDP(p->nv,mr0,p1); p1->sugar = p->sugar;
! 432: dp_monic_sf(p1,rp);
! 433: }
! 434: return;
! 435: } else if ( dp_fcoeffs )
! 436: *rp = p;
! 437: else if ( NoGCD )
! 438: dp_ptozp(p,rp);
! 439: else {
! 440: dp_ptozp(p,&p1); p = p1;
! 441: for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ );
! 442: if ( n == 1 ) {
! 443: m = BDY(p);
! 444: NEWMP(mr); mr->dl = m->dl; mr->c = (Obj)ONE; NEXT(mr) = 0;
! 445: MKDP(p->nv,mr,*rp); (*rp)->sugar = p->sugar;
! 446: return;
! 447: }
! 448: w = (P *)ALLOCA(n*sizeof(P));
! 449: c = (Q *)ALLOCA(n*sizeof(Q));
! 450: for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ )
! 451: if ( NUM(m->c) ) {
! 452: c[i] = (Q)m->c; w[i] = (P)ONE;
! 453: } else
! 454: ptozp((P)m->c,1,&c[i],&w[i]);
! 455: qltozl(c,n,&dvr); heu_nezgcdnpz(CO,w,n,&t); mulp(CO,t,(P)dvr,&g);
! 456: if ( NUM(g) )
! 457: *rp = p;
! 458: else {
! 459: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 460: NEXTMP(mr0,mr); divsp(CO,(P)m->c,g,(P *)&mr->c); mr->dl = m->dl;
! 461: }
! 462: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 463: add_denomlist(g);
! 464: }
! 465: }
! 466: }
! 467:
! 468: void heu_nezgcdnpz(VL vl,P *pl,int m,P *pr)
! 469: {
! 470: int i,r;
! 471: P gcd,t,s1,s2,u;
! 472: Z rq;
! 473: DCP dc;
! 474: extern int DP_Print;
! 475:
! 476: while ( 1 ) {
! 477: for ( i = 0, s1 = 0; i < m; i++ ) {
! 478: r = random(); UTOQ(r,rq);
! 479: mulp(vl,pl[i],(P)rq,&t); addp(vl,s1,t,&u); s1 = u;
! 480: }
! 481: for ( i = 0, s2 = 0; i < m; i++ ) {
! 482: r = random(); UTOQ(r,rq);
! 483: mulp(vl,pl[i],(P)rq,&t); addp(vl,s2,t,&u); s2 = u;
! 484: }
! 485: ezgcdp(vl,s1,s2,&gcd);
! 486: if ( DP_Print > 2 )
! 487: { fprintf(asir_out,"(%d)",nmonop(gcd)); fflush(asir_out); }
! 488: for ( i = 0; i < m; i++ ) {
! 489: if ( !divtpz(vl,pl[i],gcd,&t) )
! 490: break;
! 491: }
! 492: if ( i == m )
! 493: break;
! 494: }
! 495: *pr = gcd;
! 496: }
! 497:
! 498: void dp_prim_mod(DP p,int mod,DP *rp)
! 499: {
! 500: P t,g;
! 501: MP m,mr,mr0;
! 502:
! 503: if ( !p )
! 504: *rp = 0;
! 505: else if ( NoGCD )
! 506: *rp = p;
! 507: else {
! 508: for ( m = BDY(p), g = (P)m->c, m = NEXT(m); m; m = NEXT(m) ) {
! 509: gcdprsmp(CO,mod,g,(P)m->c,&t); g = t;
! 510: }
! 511: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 512: NEXTMP(mr0,mr); divsmp(CO,mod,(P)m->c,g,(P *)&mr->c); mr->dl = m->dl;
! 513: }
! 514: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 515: }
! 516: }
! 517:
! 518: void dp_cont(DP p,Z *rp)
! 519: {
! 520: VECT v;
! 521:
! 522: dp_dtov(p,&v); gcdvz(v,rp);
! 523: }
! 524:
! 525: void dp_dtov(DP dp,VECT *rp)
! 526: {
! 527: MP m,t;
! 528: int i,n;
! 529: VECT v;
! 530: pointer *p;
! 531:
! 532: m = BDY(dp);
! 533: for ( t = m, n = 0; t; t = NEXT(t), n++ );
! 534: MKVECT(v,n);
! 535: for ( i = 0, p = BDY(v), t = m; i < n; t = NEXT(t), i++ )
! 536: p[i] = (pointer)(t->c);
! 537: *rp = v;
! 538: }
! 539:
! 540: /*
! 541: * s-poly computation
! 542: *
! 543: */
! 544:
! 545: void dp_sp(DP p1,DP p2,DP *rp)
! 546: {
! 547: int i,n,td;
! 548: int *w;
! 549: DL d1,d2,d;
! 550: MP m;
! 551: DP t,s1,s2,u;
! 552: Z c,c1,c2;
! 553: Z gn;
! 554:
! 555: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 556: w = (int *)ALLOCA(n*sizeof(int));
! 557: for ( i = 0, td = 0; i < n; i++ ) {
! 558: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
! 559: }
! 560:
! 561: NEWDL(d,n); d->td = td - d1->td;
! 562: for ( i = 0; i < n; i++ )
! 563: d->d[i] = w[i] - d1->d[i];
! 564: c1 = (Z)BDY(p1)->c; c2 = (Z)BDY(p2)->c;
! 565: if ( INT(c1) && INT(c2) ) {
! 566: gcdz(c1,c2,&gn);
! 567: if ( !UNIQ(gn) ) {
! 568: divz(c1,gn,&c); c1 = c;
! 569: divz(c2,gn,&c);c2 = c;
! 570: }
! 571: }
! 572:
! 573: NEWMP(m); m->dl = d; m->c = (Obj)c2; NEXT(m) = 0;
! 574: MKDP(n,m,s1); s1->sugar = d->td; muld(CO,s1,p1,&t);
! 575:
! 576: NEWDL(d,n); d->td = td - d2->td;
! 577: for ( i = 0; i < n; i++ )
! 578: d->d[i] = w[i] - d2->d[i];
! 579: NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0;
! 580: MKDP(n,m,s2); s2->sugar = d->td; muld(CO,s2,p2,&u);
! 581:
! 582: subd(CO,t,u,rp);
! 583: if ( GenTrace ) {
! 584: LIST hist;
! 585: NODE node;
! 586:
! 587: node = mknode(4,ONE,NULLP,s1,ONE);
! 588: MKLIST(hist,node);
! 589: MKNODE(TraceList,hist,0);
! 590:
! 591: node = mknode(4,ONE,NULLP,NULLP,ONE);
! 592: chsgnd(s2,(DP *)&ARG2(node));
! 593: MKLIST(hist,node);
! 594: MKNODE(node,hist,TraceList); TraceList = node;
! 595: }
! 596: }
! 597:
! 598: void dpm_sp(DPM p1,DPM p2,DPM *rp)
! 599: {
! 600: int i,n,td;
! 601: int *w;
! 602: DL d1,d2,d;
! 603: MP m;
! 604: DP s1,s2;
! 605: DPM t,u;
! 606: Z c,c1,c2;
! 607: Z gn;
! 608:
! 609: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 610: if ( BDY(p1)->pos != BDY(p2)->pos ) {
! 611: *rp = 0;
! 612: return;
! 613: }
! 614: w = (int *)ALLOCA(n*sizeof(int));
! 615: for ( i = 0, td = 0; i < n; i++ ) {
! 616: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
! 617: }
! 618:
! 619: NEWDL(d,n); d->td = td - d1->td;
! 620: for ( i = 0; i < n; i++ )
! 621: d->d[i] = w[i] - d1->d[i];
! 622: c1 = (Z)BDY(p1)->c; c2 = (Z)BDY(p2)->c;
! 623: if ( INT(c1) && INT(c2) ) {
! 624: gcdz(c1,c2,&gn);
! 625: if ( !UNIQ(gn) ) {
! 626: divz(c1,gn,&c); c1 = c;
! 627: divz(c2,gn,&c);c2 = c;
! 628: }
! 629: }
! 630:
! 631: NEWMP(m); m->dl = d; m->c = (Obj)c2; NEXT(m) = 0;
! 632: MKDP(n,m,s1); s1->sugar = d->td; mulobjdpm(CO,(Obj)s1,p1,&t);
! 633:
! 634: NEWDL(d,n); d->td = td - d2->td;
! 635: for ( i = 0; i < n; i++ )
! 636: d->d[i] = w[i] - d2->d[i];
! 637: NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0;
! 638: MKDP(n,m,s2); s2->sugar = d->td; mulobjdpm(CO,(Obj)s2,p2,&u);
! 639:
! 640: subdpm(CO,t,u,rp);
! 641: if ( GenTrace ) {
! 642: LIST hist;
! 643: NODE node;
! 644:
! 645: node = mknode(4,ONE,NULLP,s1,ONE);
! 646: MKLIST(hist,node);
! 647: MKNODE(TraceList,hist,0);
! 648:
! 649: node = mknode(4,ONE,NULLP,NULLP,ONE);
! 650: chsgnd(s2,(DP *)&ARG2(node));
! 651: MKLIST(hist,node);
! 652: MKNODE(node,hist,TraceList); TraceList = node;
! 653: }
! 654: }
! 655:
! 656: void _dp_sp_dup(DP p1,DP p2,DP *rp)
! 657: {
! 658: int i,n,td;
! 659: int *w;
! 660: DL d1,d2,d;
! 661: MP m;
! 662: DP t,s1,s2,u;
! 663: Z c,c1,c2;
! 664: Z gn;
! 665:
! 666: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 667: w = (int *)ALLOCA(n*sizeof(int));
! 668: for ( i = 0, td = 0; i < n; i++ ) {
! 669: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
! 670: }
! 671:
! 672: _NEWDL(d,n); d->td = td - d1->td;
! 673: for ( i = 0; i < n; i++ )
! 674: d->d[i] = w[i] - d1->d[i];
! 675: c1 = (Z)BDY(p1)->c; c2 = (Z)BDY(p2)->c;
! 676: if ( INT(c1) && INT(c2) ) {
! 677: gcdz(c1,c2,&gn);
! 678: if ( !UNIQ(gn) ) {
! 679: divz(c1,gn,&c); c1 = c;
! 680: divz(c2,gn,&c);c2 = c;
! 681: }
! 682: }
! 683:
! 684: _NEWMP(m); m->dl = d; m->c = (Obj)c2; NEXT(m) = 0;
! 685: _MKDP(n,m,s1); s1->sugar = d->td; _muld_dup(CO,s1,p1,&t); _free_dp(s1);
! 686:
! 687: _NEWDL(d,n); d->td = td - d2->td;
! 688: for ( i = 0; i < n; i++ )
! 689: d->d[i] = w[i] - d2->d[i];
! 690: _NEWMP(m); m->dl = d; chsgnp((P)c1,(P *)&m->c); NEXT(m) = 0;
! 691: _MKDP(n,m,s2); s2->sugar = d->td; _muld_dup(CO,s2,p2,&u); _free_dp(s2);
! 692:
! 693: _addd_destructive(CO,t,u,rp);
! 694: if ( GenTrace ) {
! 695: LIST hist;
! 696: NODE node;
! 697:
! 698: node = mknode(4,ONE,NULLP,s1,ONE);
! 699: MKLIST(hist,node);
! 700: MKNODE(TraceList,hist,0);
! 701:
! 702: node = mknode(4,ONE,NULLP,NULLP,ONE);
! 703: chsgnd(s2,(DP *)&ARG2(node));
! 704: MKLIST(hist,node);
! 705: MKNODE(node,hist,TraceList); TraceList = node;
! 706: }
! 707: }
! 708:
! 709: void dp_sp_mod(DP p1,DP p2,int mod,DP *rp)
! 710: {
! 711: int i,n,td;
! 712: int *w;
! 713: DL d1,d2,d;
! 714: MP m;
! 715: DP t,s,u;
! 716:
! 717: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 718: w = (int *)ALLOCA(n*sizeof(int));
! 719: for ( i = 0, td = 0; i < n; i++ ) {
! 720: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
! 721: }
! 722: NEWDL_NOINIT(d,n); d->td = td - d1->td;
! 723: for ( i = 0; i < n; i++ )
! 724: d->d[i] = w[i] - d1->d[i];
! 725: NEWMP(m); m->dl = d; m->c = (Obj)BDY(p2)->c; NEXT(m) = 0;
! 726: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p1,s,&t);
! 727: NEWDL_NOINIT(d,n); d->td = td - d2->td;
! 728: for ( i = 0; i < n; i++ )
! 729: d->d[i] = w[i] - d2->d[i];
! 730: NEWMP(m); m->dl = d; m->c = (Obj)BDY(p1)->c; NEXT(m) = 0;
! 731: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p2,s,&u);
! 732: submd(CO,mod,t,u,rp);
! 733: }
! 734:
! 735: void _dp_sp_mod_dup(DP p1,DP p2,int mod,DP *rp)
! 736: {
! 737: int i,n,td;
! 738: int *w;
! 739: DL d1,d2,d;
! 740: MP m;
! 741: DP t,s,u;
! 742:
! 743: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 744: w = (int *)ALLOCA(n*sizeof(int));
! 745: for ( i = 0, td = 0; i < n; i++ ) {
! 746: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
! 747: }
! 748: _NEWDL(d,n); d->td = td - d1->td;
! 749: for ( i = 0; i < n; i++ )
! 750: d->d[i] = w[i] - d1->d[i];
! 751: _NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
! 752: _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p1,&t); _free_dp(s);
! 753: _NEWDL(d,n); d->td = td - d2->td;
! 754: for ( i = 0; i < n; i++ )
! 755: d->d[i] = w[i] - d2->d[i];
! 756: _NEWMP(m); m->dl = d; m->c = (Obj)STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
! 757: _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p2,&u); _free_dp(s);
! 758: _addmd_destructive(mod,t,u,rp);
! 759: }
! 760:
! 761: void _dp_sp_mod(DP p1,DP p2,int mod,DP *rp)
! 762: {
! 763: int i,n,td;
! 764: int *w;
! 765: DL d1,d2,d;
! 766: MP m;
! 767: DP t,s,u;
! 768:
! 769: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 770: w = (int *)ALLOCA(n*sizeof(int));
! 771: for ( i = 0, td = 0; i < n; i++ ) {
! 772: w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i);
! 773: }
! 774: NEWDL(d,n); d->td = td - d1->td;
! 775: for ( i = 0; i < n; i++ )
! 776: d->d[i] = w[i] - d1->d[i];
! 777: NEWMP(m); m->dl = d; m->c = BDY(p2)->c; NEXT(m) = 0;
! 778: MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p1,&t);
! 779: NEWDL(d,n); d->td = td - d2->td;
! 780: for ( i = 0; i < n; i++ )
! 781: d->d[i] = w[i] - d2->d[i];
! 782: NEWMP(m); m->dl = d; m->c = (Obj)STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0;
! 783: MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p2,&u);
! 784: addmd_destructive(mod,t,u,rp);
! 785: }
! 786:
! 787: /*
! 788: * m-reduction
! 789: * do content reduction over Z or Q(x,...)
! 790: * do nothing over finite fields
! 791: *
! 792: */
! 793:
! 794: void dp_red(DP p0,DP p1,DP p2,DP *head,DP *rest,P *dnp,DP *multp)
! 795: {
! 796: int i,n;
! 797: DL d1,d2,d;
! 798: MP m;
! 799: DP t,s,r,h;
! 800: Z c,c1,c2,gn;
! 801: P g,a;
! 802: P p[2];
! 803:
! 804: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 805: NEWDL(d,n); d->td = d1->td - d2->td;
! 806: for ( i = 0; i < n; i++ )
! 807: d->d[i] = d1->d[i]-d2->d[i];
! 808: c1 = (Z)BDY(p1)->c; c2 = (Z)BDY(p2)->c;
! 809: if ( dp_fcoeffs == N_GFS ) {
! 810: p[0] = (P)c1; p[1] = (P)c2;
! 811: gcdsf(CO,p,2,&g);
! 812: divsp(CO,(P)c1,g,&a); c1 = (Z)a; divsp(CO,(P)c2,g,&a); c2 = (Z)a;
! 813: } else if ( dp_fcoeffs ) {
! 814: /* do nothing */
! 815: } else if ( INT(c1) && INT(c2) ) {
! 816: gcdz(c1,c2,&gn);
! 817: if ( !UNIQ(gn) ) {
! 818: divz(c1,gn,&c); c1 = c;
! 819: divz(c2,gn,&c); c2 = c;
! 820: }
! 821: } else {
! 822: ezgcdpz(CO,(P)c1,(P)c2,&g);
! 823: divsp(CO,(P)c1,g,&a); c1 = (Z)a; divsp(CO,(P)c2,g,&a); c2 = (Z)a;
! 824: add_denomlist(g);
! 825: }
! 826: NEWMP(m); m->dl = d; chsgnp((P)c1,(P *)&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
! 827: *multp = s;
! 828: muld(CO,s,p2,&t); muldc(CO,p1,(Obj)c2,&s); addd(CO,s,t,&r);
! 829: muldc(CO,p0,(Obj)c2,&h);
! 830: *head = h; *rest = r; *dnp = (P)c2;
! 831: }
! 832:
! 833: void dpm_red(DPM p0,DPM p1,DPM p2,DPM *head,DPM *rest,P *dnp,DP *multp)
! 834: {
! 835: int i,n,pos;
! 836: DL d1,d2,d;
! 837: MP m;
! 838: DP s;
! 839: DPM t,r,h,u,w;
! 840: Z c,c1,c2,gn;
! 841: P g,a;
! 842: P p[2];
! 843:
! 844: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl; pos = BDY(p1)->pos;
! 845: if ( pos != BDY(p2)->pos )
! 846: error("dpm_red : cannot happen");
! 847: NEWDL(d,n); d->td = d1->td - d2->td;
! 848: for ( i = 0; i < n; i++ )
! 849: d->d[i] = d1->d[i]-d2->d[i];
! 850: c1 = (Z)BDY(p1)->c; c2 = (Z)BDY(p2)->c;
! 851: if ( dp_fcoeffs == N_GFS ) {
! 852: p[0] = (P)c1; p[1] = (P)c2;
! 853: gcdsf(CO,p,2,&g);
! 854: divsp(CO,(P)c1,g,&a); c1 = (Z)a; divsp(CO,(P)c2,g,&a); c2 = (Z)a;
! 855: } else if ( dp_fcoeffs ) {
! 856: /* do nothing */
! 857: } else if ( INT(c1) && INT(c2) ) {
! 858: gcdz(c1,c2,&gn);
! 859: if ( !UNIQ(gn) ) {
! 860: divz(c1,gn,&c); c1 = c;
! 861: divz(c2,gn,&c); c2 = c;
! 862: }
! 863: } else {
! 864: ezgcdpz(CO,(P)c1,(P)c2,&g);
! 865: divsp(CO,(P)c1,g,&a); c1 = (Z)a; divsp(CO,(P)c2,g,&a); c2 = (Z)a;
! 866: add_denomlist(g);
! 867: }
! 868: NEWMP(m); m->dl = d; chsgnp((P)c1,(P *)&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
! 869: *multp = s;
! 870: mulobjdpm(CO,(Obj)s,p2,&u); mulobjdpm(CO,(Obj)c2,p1,&w); adddpm(CO,u,w,&r);
! 871: mulobjdpm(CO,(Obj)c2,p0,&h);
! 872: *head = h; *rest = r; *dnp = (P)c2;
! 873: }
! 874:
! 875:
! 876: /*
! 877: * m-reduction by a marked poly
! 878: * do content reduction over Z or Q(x,...)
! 879: * do nothing over finite fields
! 880: *
! 881: */
! 882:
! 883:
! 884: void dp_red_marked(DP p0,DP p1,DP p2,DP hp2,DP *head,DP *rest,P *dnp,DP *multp)
! 885: {
! 886: int i,n;
! 887: DL d1,d2,d;
! 888: MP m;
! 889: DP t,s,r,h;
! 890: Z c,c1,c2,gn;
! 891: P g,a;
! 892: P p[2];
! 893:
! 894: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(hp2)->dl;
! 895: NEWDL(d,n); d->td = d1->td - d2->td;
! 896: for ( i = 0; i < n; i++ )
! 897: d->d[i] = d1->d[i]-d2->d[i];
! 898: c1 = (Z)BDY(p1)->c; c2 = (Z)BDY(hp2)->c;
! 899: if ( dp_fcoeffs == N_GFS ) {
! 900: p[0] = (P)c1; p[1] = (P)c2;
! 901: gcdsf(CO,p,2,&g);
! 902: divsp(CO,(P)c1,g,&a); c1 = (Z)a; divsp(CO,(P)c2,g,&a); c2 = (Z)a;
! 903: } else if ( dp_fcoeffs ) {
! 904: /* do nothing */
! 905: } else if ( INT(c1) && INT(c2) ) {
! 906: gcdz(c1,c2,&gn);
! 907: if ( !UNIQ(gn) ) {
! 908: divz(c1,gn,&c); c1 = c;
! 909: divz(c2,gn,&c); c2 = c;
! 910: }
! 911: } else {
! 912: ezgcdpz(CO,(P)c1,(P)c2,&g);
! 913: divsp(CO,(P)c1,g,&a); c1 = (Z)a; divsp(CO,(P)c2,g,&a); c2 = (Z)a;
! 914: }
! 915: NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
! 916: *multp = s;
! 917: muld(CO,s,p2,&t); muldc(CO,p1,(Obj)c2,&s); subd(CO,s,t,&r);
! 918: muldc(CO,p0,(Obj)c2,&h);
! 919: *head = h; *rest = r; *dnp = (P)c2;
! 920: }
! 921:
! 922: void dp_red_marked_mod(DP p0,DP p1,DP p2,DP hp2,int mod,DP *head,DP *rest,P *dnp,DP *multp)
! 923: {
! 924: int i,n;
! 925: DL d1,d2,d;
! 926: MP m;
! 927: DP t,s,r,h;
! 928: P c1,c2,g,u;
! 929:
! 930: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(hp2)->dl;
! 931: NEWDL(d,n); d->td = d1->td - d2->td;
! 932: for ( i = 0; i < n; i++ )
! 933: d->d[i] = d1->d[i]-d2->d[i];
! 934: c1 = (P)BDY(p1)->c; c2 = (P)BDY(hp2)->c;
! 935: gcdprsmp(CO,mod,c1,c2,&g);
! 936: divsmp(CO,mod,c1,g,&u); c1 = u; divsmp(CO,mod,c2,g,&u); c2 = u;
! 937: if ( NUM(c2) ) {
! 938: divsmp(CO,mod,c1,c2,&u); c1 = u; c2 = (P)ONEM;
! 939: }
! 940: NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0;
! 941: MKDP(n,m,s); s->sugar = d->td;
! 942: *multp = s;
! 943: mulmd(CO,mod,s,p2,&t);
! 944: if ( NUM(c2) ) {
! 945: submd(CO,mod,p1,t,&r); h = p0;
! 946: } else {
! 947: mulmdc(CO,mod,p1,c2,&s); submd(CO,mod,s,t,&r); mulmdc(CO,mod,p0,c2,&h);
! 948: }
! 949: *head = h; *rest = r; *dnp = c2;
! 950: }
! 951:
! 952: /* m-reduction over a field */
! 953:
! 954: void dp_red_f(DP p1,DP p2,DP *rest)
! 955: {
! 956: int i,n;
! 957: DL d1,d2,d;
! 958: MP m;
! 959: DP t,s;
! 960: Obj a,b;
! 961:
! 962: n = p1->nv;
! 963: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 964:
! 965: NEWDL(d,n); d->td = d1->td - d2->td;
! 966: for ( i = 0; i < n; i++ )
! 967: d->d[i] = d1->d[i]-d2->d[i];
! 968:
! 969: NEWMP(m); m->dl = d;
! 970: divr(CO,(Obj)BDY(p1)->c,(Obj)BDY(p2)->c,&a); chsgnr(a,&b);
! 971: C(m) = (Obj)b;
! 972: NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
! 973:
! 974: muld(CO,s,p2,&t); addd(CO,p1,t,rest);
! 975: }
! 976:
! 977: void dpm_red_f(DPM p1,DPM p2,DPM *rest)
! 978: {
! 979: int i,n;
! 980: DL d1,d2,d;
! 981: MP m;
! 982: DPM t;
! 983: DP s;
! 984: Obj a,b;
! 985:
! 986: n = p1->nv;
! 987: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 988:
! 989: NEWDL(d,n); d->td = d1->td - d2->td;
! 990: for ( i = 0; i < n; i++ )
! 991: d->d[i] = d1->d[i]-d2->d[i];
! 992:
! 993: NEWMP(m); m->dl = d;
! 994: arf_div(CO,(Obj)BDY(p1)->c,(Obj)BDY(p2)->c,&a); arf_chsgn(a,&b);
! 995: C(m) = b;
! 996: NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td;
! 997:
! 998: mulobjdpm(CO,(Obj)s,p2,&t); adddpm(CO,p1,t,rest);
! 999: }
! 1000:
! 1001:
! 1002: void dp_red_mod(DP p0,DP p1,DP p2,int mod,DP *head,DP *rest,P *dnp)
! 1003: {
! 1004: int i,n;
! 1005: DL d1,d2,d;
! 1006: MP m;
! 1007: DP t,s,r,h;
! 1008: P c1,c2,g,u;
! 1009:
! 1010: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 1011: NEWDL(d,n); d->td = d1->td - d2->td;
! 1012: for ( i = 0; i < n; i++ )
! 1013: d->d[i] = d1->d[i]-d2->d[i];
! 1014: c1 = (P)BDY(p1)->c; c2 = (P)BDY(p2)->c;
! 1015: gcdprsmp(CO,mod,c1,c2,&g);
! 1016: divsmp(CO,mod,c1,g,&u); c1 = u; divsmp(CO,mod,c2,g,&u); c2 = u;
! 1017: if ( NUM(c2) ) {
! 1018: divsmp(CO,mod,c1,c2,&u); c1 = u; c2 = (P)ONEM;
! 1019: }
! 1020: NEWMP(m); m->dl = d; chsgnmp(mod,(P)c1,(P *)&m->c); NEXT(m) = 0;
! 1021: MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,s,p2,&t);
! 1022: if ( NUM(c2) ) {
! 1023: addmd(CO,mod,p1,t,&r); h = p0;
! 1024: } else {
! 1025: mulmdc(CO,mod,p1,c2,&s); addmd(CO,mod,s,t,&r); mulmdc(CO,mod,p0,c2,&h);
! 1026: }
! 1027: *head = h; *rest = r; *dnp = c2;
! 1028: }
! 1029:
! 1030: struct oEGT eg_red_mod;
! 1031:
! 1032: void _dp_red_mod_destructive(DP p1,DP p2,int mod,DP *rp)
! 1033: {
! 1034: int i,n;
! 1035: DL d1,d2,d;
! 1036: MP m;
! 1037: DP t,s;
! 1038: int c,c1,c2;
! 1039: extern int do_weyl;
! 1040:
! 1041: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 1042: _NEWDL(d,n); d->td = d1->td - d2->td;
! 1043: for ( i = 0; i < n; i++ )
! 1044: d->d[i] = d1->d[i]-d2->d[i];
! 1045: c = invm(ITOS(BDY(p2)->c),mod);
! 1046: c2 = ITOS(BDY(p1)->c);
! 1047: DMAR(c,c2,0,mod,c1);
! 1048: _NEWMP(m); m->dl = d; m->c = (Obj)STOI(mod-c1); NEXT(m) = 0;
! 1049: #if 0
! 1050: _MKDP(n,m,s); s->sugar = d->td;
! 1051: _mulmd_dup(mod,s,p2,&t); _free_dp(s);
! 1052: #else
! 1053: if ( do_weyl ) {
! 1054: _MKDP(n,m,s); s->sugar = d->td;
! 1055: _mulmd_dup(mod,s,p2,&t); _free_dp(s);
! 1056: } else {
! 1057: _mulmdm_dup(mod,p2,m,&t); _FREEMP(m);
! 1058: }
! 1059: #endif
! 1060: /* get_eg(&t0); */
! 1061: _addmd_destructive(mod,p1,t,rp);
! 1062: /* get_eg(&t1); add_eg(&eg_red_mod,&t0,&t1); */
! 1063: }
! 1064:
! 1065: /*
! 1066: * normal form computation
! 1067: *
! 1068: */
! 1069:
! 1070: void dp_true_nf(NODE b,DP g,DP *ps,int full,DP *rp,P *dnp)
! 1071: {
! 1072: DP u,p,d,s,t,dmy;
! 1073: NODE l;
! 1074: MP m,mr;
! 1075: int i,n;
! 1076: int *wb;
! 1077: int sugar,psugar;
! 1078: P dn,tdn,tdn1;
! 1079:
! 1080: dn = (P)ONE;
! 1081: if ( !g ) {
! 1082: *rp = 0; *dnp = dn; return;
! 1083: }
! 1084: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1085: wb = (int *)ALLOCA(n*sizeof(int));
! 1086: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1087: wb[i] = QTOS((Q)BDY(l));
! 1088: sugar = g->sugar;
! 1089: for ( d = 0; g; ) {
! 1090: for ( u = 0, i = 0; i < n; i++ ) {
! 1091: if ( dp_redble(g,p = ps[wb[i]]) ) {
! 1092: dp_red(d,g,p,&t,&u,&tdn,&dmy);
! 1093: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1094: sugar = MAX(sugar,psugar);
! 1095: if ( !u ) {
! 1096: if ( d )
! 1097: d->sugar = sugar;
! 1098: *rp = d; *dnp = dn; return;
! 1099: } else {
! 1100: d = t;
! 1101: mulp(CO,dn,tdn,&tdn1); dn = tdn1;
! 1102: }
! 1103: break;
! 1104: }
! 1105: }
! 1106: if ( u )
! 1107: g = u;
! 1108: else if ( !full ) {
! 1109: if ( g ) {
! 1110: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1111: }
! 1112: *rp = g; *dnp = dn; return;
! 1113: } else {
! 1114: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1115: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1116: addd(CO,d,t,&s); d = s;
! 1117: dp_rest(g,&t); g = t;
! 1118: }
! 1119: }
! 1120: if ( d )
! 1121: d->sugar = sugar;
! 1122: *rp = d; *dnp = dn;
! 1123: }
! 1124:
! 1125: void dp_removecont2(DP p1,DP p2,DP *r1p,DP *r2p,Z *contp)
! 1126: {
! 1127: struct oVECT v;
! 1128: int i,n1,n2,n;
! 1129: MP m,m0,t;
! 1130: Z *w;
! 1131: Z h;
! 1132:
! 1133: if ( p1 ) {
! 1134: for ( i = 0, m = BDY(p1); m; m = NEXT(m), i++ );
! 1135: n1 = i;
! 1136: } else
! 1137: n1 = 0;
! 1138: if ( p2 ) {
! 1139: for ( i = 0, m = BDY(p2); m; m = NEXT(m), i++ );
! 1140: n2 = i;
! 1141: } else
! 1142: n2 = 0;
! 1143: n = n1+n2;
! 1144: if ( !n ) {
! 1145: *r1p = 0; *r2p = 0; *contp = ONE; return;
! 1146: }
! 1147: w = (Z *)ALLOCA(n*sizeof(Q));
! 1148: v.len = n;
! 1149: v.body = (pointer *)w;
! 1150: i = 0;
! 1151: if ( p1 )
! 1152: for ( m = BDY(p1); i < n1; m = NEXT(m), i++ ) w[i] = (Z)m->c;
! 1153: if ( p2 )
! 1154: for ( m = BDY(p2); i < n; m = NEXT(m), i++ ) w[i] = (Z)m->c;
! 1155: h = w[0]; removecont_array((P *)w,n,1); divz(h,w[0],contp);
! 1156: i = 0;
! 1157: if ( p1 ) {
! 1158: for ( m0 = 0, t = BDY(p1); i < n1; i++, t = NEXT(t) ) {
! 1159: NEXTMP(m0,m); m->c = (Obj)w[i]; m->dl = t->dl;
! 1160: }
! 1161: NEXT(m) = 0;
! 1162: MKDP(p1->nv,m0,*r1p); (*r1p)->sugar = p1->sugar;
! 1163: } else
! 1164: *r1p = 0;
! 1165: if ( p2 ) {
! 1166: for ( m0 = 0, t = BDY(p2); i < n; i++, t = NEXT(t) ) {
! 1167: NEXTMP(m0,m); m->c = (Obj)w[i]; m->dl = t->dl;
! 1168: }
! 1169: NEXT(m) = 0;
! 1170: MKDP(p2->nv,m0,*r2p); (*r2p)->sugar = p2->sugar;
! 1171: } else
! 1172: *r2p = 0;
! 1173: }
! 1174:
! 1175: /* true nf by a marked GB */
! 1176:
! 1177: void dp_true_nf_marked(NODE b,DP g,DP *ps,DP *hps,DP *rp,P *nmp,P *dnp)
! 1178: {
! 1179: DP u,p,d,s,t,dmy,hp;
! 1180: NODE l;
! 1181: MP m,mr;
! 1182: int i,n,hmag;
! 1183: int *wb;
! 1184: int sugar,psugar,multiple;
! 1185: P nm,tnm1,dn,tdn,tdn1;
! 1186: Z cont;
! 1187:
! 1188: multiple = 0;
! 1189: hmag = multiple*HMAG(g);
! 1190: nm = (P)ONE;
! 1191: dn = (P)ONE;
! 1192: if ( !g ) {
! 1193: *rp = 0; *dnp = dn; return;
! 1194: }
! 1195: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1196: wb = (int *)ALLOCA(n*sizeof(int));
! 1197: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1198: wb[i] = QTOS((Z)BDY(l));
! 1199: sugar = g->sugar;
! 1200: for ( d = 0; g; ) {
! 1201: for ( u = 0, i = 0; i < n; i++ ) {
! 1202: if ( dp_redble(g,hp = hps[wb[i]]) ) {
! 1203: p = ps[wb[i]];
! 1204: dp_red_marked(d,g,p,hp,&t,&u,&tdn,&dmy);
! 1205: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1206: sugar = MAX(sugar,psugar);
! 1207: if ( !u ) {
! 1208: goto last;
! 1209: } else {
! 1210: d = t;
! 1211: mulp(CO,dn,tdn,&tdn1); dn = tdn1;
! 1212: }
! 1213: break;
! 1214: }
! 1215: }
! 1216: if ( u ) {
! 1217: g = u;
! 1218: if ( multiple && ((d && HMAG(d)>hmag) || (HMAG(g)>hmag)) ) {
! 1219: dp_removecont2(d,g,&t,&u,&cont); d = t; g = u;
! 1220: mulp(CO,nm,(P)cont,&tnm1); nm = tnm1;
! 1221: if ( d )
! 1222: hmag = multiple*HMAG(d);
! 1223: else
! 1224: hmag = multiple*HMAG(g);
! 1225: }
! 1226: } else {
! 1227: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1228: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1229: addd(CO,d,t,&s); d = s;
! 1230: dp_rest(g,&t); g = t;
! 1231: }
! 1232: }
! 1233: last:
! 1234: if ( d ) {
! 1235: dp_removecont2(d,0,&t,&u,&cont); d = t;
! 1236: mulp(CO,nm,(P)cont,&tnm1); nm = tnm1;
! 1237: d->sugar = sugar;
! 1238: }
! 1239: *rp = d; *nmp = nm; *dnp = dn;
! 1240: }
! 1241:
! 1242: void dp_true_nf_marked_mod(NODE b,DP g,DP *ps,DP *hps,int mod,DP *rp,P *dnp)
! 1243: {
! 1244: DP hp,u,p,d,s,t,dmy;
! 1245: NODE l;
! 1246: MP m,mr;
! 1247: int i,n;
! 1248: int *wb;
! 1249: int sugar,psugar;
! 1250: P dn,tdn,tdn1;
! 1251:
! 1252: dn = (P)ONEM;
! 1253: if ( !g ) {
! 1254: *rp = 0; *dnp = dn; return;
! 1255: }
! 1256: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1257: wb = (int *)ALLOCA(n*sizeof(int));
! 1258: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1259: wb[i] = QTOS((Q)BDY(l));
! 1260: sugar = g->sugar;
! 1261: for ( d = 0; g; ) {
! 1262: for ( u = 0, i = 0; i < n; i++ ) {
! 1263: if ( dp_redble(g,hp = hps[wb[i]]) ) {
! 1264: p = ps[wb[i]];
! 1265: dp_red_marked_mod(d,g,p,hp,mod,&t,&u,&tdn,&dmy);
! 1266: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1267: sugar = MAX(sugar,psugar);
! 1268: if ( !u ) {
! 1269: if ( d )
! 1270: d->sugar = sugar;
! 1271: *rp = d; *dnp = dn; return;
! 1272: } else {
! 1273: d = t;
! 1274: mulmp(CO,mod,dn,tdn,&tdn1); dn = tdn1;
! 1275: }
! 1276: break;
! 1277: }
! 1278: }
! 1279: if ( u )
! 1280: g = u;
! 1281: else {
! 1282: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1283: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1284: addmd(CO,mod,d,t,&s); d = s;
! 1285: dp_rest(g,&t); g = t;
! 1286: }
! 1287: }
! 1288: if ( d )
! 1289: d->sugar = sugar;
! 1290: *rp = d; *dnp = dn;
! 1291: }
! 1292:
! 1293: /* true nf by a marked GB and collect quotients */
! 1294:
! 1295: DP *dp_true_nf_and_quotient_marked (NODE b,DP g,DP *ps,DP *hps,DP *rp,P *dnp)
! 1296: {
! 1297: DP u,p,d,s,t,dmy,hp,mult;
! 1298: DP *q;
! 1299: NODE l;
! 1300: MP m,mr;
! 1301: int i,n,j;
! 1302: int *wb;
! 1303: int sugar,psugar,multiple;
! 1304: P nm,tnm1,dn,tdn,tdn1;
! 1305: Q cont;
! 1306:
! 1307: dn = (P)ONE;
! 1308: if ( !g ) {
! 1309: *rp = 0; *dnp = dn; return 0;
! 1310: }
! 1311: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1312: wb = (int *)ALLOCA(n*sizeof(int));
! 1313: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1314: wb[i] = QTOS((Q)BDY(l));
! 1315: q = (DP *)MALLOC(n*sizeof(DP));
! 1316: for ( i = 0; i < n; i++ ) q[i] = 0;
! 1317: sugar = g->sugar;
! 1318: for ( d = 0; g; ) {
! 1319: for ( u = 0, i = 0; i < n; i++ ) {
! 1320: if ( dp_redble(g,hp = hps[wb[i]]) ) {
! 1321: p = ps[wb[i]];
! 1322: dp_red_marked(d,g,p,hp,&t,&u,&tdn,&mult);
! 1323: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1324: sugar = MAX(sugar,psugar);
! 1325: for ( j = 0; j < n; j++ ) {
! 1326: muldc(CO,q[j],(Obj)tdn,&dmy); q[j] = dmy;
! 1327: }
! 1328: addd(CO,q[wb[i]],mult,&dmy); q[wb[i]] = dmy;
! 1329: mulp(CO,dn,tdn,&tdn1); dn = tdn1;
! 1330: d = t;
! 1331: if ( !u ) goto last;
! 1332: break;
! 1333: }
! 1334: }
! 1335: if ( u ) {
! 1336: g = u;
! 1337: } else {
! 1338: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1339: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1340: addd(CO,d,t,&s); d = s;
! 1341: dp_rest(g,&t); g = t;
! 1342: }
! 1343: }
! 1344: last:
! 1345: if ( d ) d->sugar = sugar;
! 1346: *rp = d; *dnp = dn;
! 1347: return q;
! 1348: }
! 1349:
! 1350: DP *dp_true_nf_and_quotient_marked_mod(NODE b,DP g,DP *ps,DP *hps,int mod,DP *rp,P *dnp)
! 1351: {
! 1352: DP u,p,d,s,t,dmy,hp,mult;
! 1353: DP *q;
! 1354: NODE l;
! 1355: MP m,mr;
! 1356: int i,n,j;
! 1357: int *wb;
! 1358: int sugar,psugar;
! 1359: P dn,tdn,tdn1;
! 1360:
! 1361: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1362: q = (DP *)MALLOC(n*sizeof(DP));
! 1363: for ( i = 0; i < n; i++ ) q[i] = 0;
! 1364: dn = (P)ONEM;
! 1365: if ( !g ) {
! 1366: *rp = 0; *dnp = dn; return 0;
! 1367: }
! 1368: wb = (int *)ALLOCA(n*sizeof(int));
! 1369: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1370: wb[i] = QTOS((Q)BDY(l));
! 1371: sugar = g->sugar;
! 1372: for ( d = 0; g; ) {
! 1373: for ( u = 0, i = 0; i < n; i++ ) {
! 1374: if ( dp_redble(g,hp = hps[wb[i]]) ) {
! 1375: p = ps[wb[i]];
! 1376: dp_red_marked_mod(d,g,p,hp,mod,&t,&u,&tdn,&mult);
! 1377: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1378: sugar = MAX(sugar,psugar);
! 1379: for ( j = 0; j < n; j++ ) {
! 1380: mulmdc(CO,mod,q[j],(P)tdn,&dmy); q[j] = dmy;
! 1381: }
! 1382: addmd(CO,mod,q[wb[i]],mult,&dmy); q[wb[i]] = dmy;
! 1383: mulmp(CO,mod,dn,tdn,&tdn1); dn = tdn1;
! 1384: d = t;
! 1385: if ( !u ) goto last;
! 1386: break;
! 1387: }
! 1388: }
! 1389: if ( u )
! 1390: g = u;
! 1391: else {
! 1392: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1393: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1394: addmd(CO,mod,d,t,&s); d = s;
! 1395: dp_rest(g,&t); g = t;
! 1396: }
! 1397: }
! 1398: last:
! 1399: if ( d )
! 1400: d->sugar = sugar;
! 1401: *rp = d; *dnp = dn;
! 1402: return q;
! 1403: }
! 1404:
! 1405: /* nf computation over Z */
! 1406:
! 1407: void dp_nf_z(NODE b,DP g,DP *ps,int full,int multiple,DP *rp)
! 1408: {
! 1409: DP u,p,d,s,t,dmy1;
! 1410: P dmy;
! 1411: NODE l;
! 1412: MP m,mr;
! 1413: int i,n;
! 1414: int *wb;
! 1415: int hmag;
! 1416: int sugar,psugar;
! 1417:
! 1418: if ( !g ) {
! 1419: *rp = 0; return;
! 1420: }
! 1421: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1422: wb = (int *)ALLOCA(n*sizeof(int));
! 1423: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1424: wb[i] = QTOS((Q)BDY(l));
! 1425:
! 1426: hmag = multiple*HMAG(g);
! 1427: sugar = g->sugar;
! 1428:
! 1429: for ( d = 0; g; ) {
! 1430: for ( u = 0, i = 0; i < n; i++ ) {
! 1431: if ( dp_redble(g,p = ps[wb[i]]) ) {
! 1432: dp_red(d,g,p,&t,&u,&dmy,&dmy1);
! 1433: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1434: sugar = MAX(sugar,psugar);
! 1435: if ( !u ) {
! 1436: if ( d )
! 1437: d->sugar = sugar;
! 1438: *rp = d; return;
! 1439: }
! 1440: d = t;
! 1441: break;
! 1442: }
! 1443: }
! 1444: if ( u ) {
! 1445: g = u;
! 1446: if ( d ) {
! 1447: if ( multiple && HMAG(d) > hmag ) {
! 1448: dp_ptozp2(d,g,&t,&u); d = t; g = u;
! 1449: hmag = multiple*HMAG(d);
! 1450: }
! 1451: } else {
! 1452: if ( multiple && HMAG(g) > hmag ) {
! 1453: dp_ptozp(g,&t); g = t;
! 1454: hmag = multiple*HMAG(g);
! 1455: }
! 1456: }
! 1457: }
! 1458: else if ( !full ) {
! 1459: if ( g ) {
! 1460: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1461: }
! 1462: *rp = g; return;
! 1463: } else {
! 1464: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1465: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1466: addd(CO,d,t,&s); d = s;
! 1467: dp_rest(g,&t); g = t;
! 1468:
! 1469: }
! 1470: }
! 1471: if ( d )
! 1472: d->sugar = sugar;
! 1473: *rp = d;
! 1474: }
! 1475:
! 1476: void dpm_nf_z(NODE b,DPM g,DPM *ps,int full,int multiple,DPM *rp)
! 1477: {
! 1478: DPM u,p,d,s,t;
! 1479: DP dmy1;
! 1480: P dmy;
! 1481: NODE l;
! 1482: DMM m,mr;
! 1483: int i,n;
! 1484: int *wb;
! 1485: int hmag;
! 1486: int sugar,psugar;
! 1487:
! 1488: if ( !g ) {
! 1489: *rp = 0; return;
! 1490: }
! 1491: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1492: wb = (int *)ALLOCA(n*sizeof(int));
! 1493: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1494: wb[i] = QTOS((Q)BDY(l));
! 1495:
! 1496: hmag = multiple*HMAG(g);
! 1497: sugar = g->sugar;
! 1498:
! 1499: for ( d = 0; g; ) {
! 1500: for ( u = 0, i = 0; i < n; i++ ) {
! 1501: if ( dpm_redble(g,p = ps[wb[i]]) ) {
! 1502: dpm_red(d,g,p,&t,&u,&dmy,&dmy1);
! 1503: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1504: sugar = MAX(sugar,psugar);
! 1505: if ( !u ) {
! 1506: if ( d )
! 1507: d->sugar = sugar;
! 1508: *rp = d; return;
! 1509: }
! 1510: d = t;
! 1511: break;
! 1512: }
! 1513: }
! 1514: if ( u ) {
! 1515: g = u;
! 1516: if ( d ) {
! 1517: if ( multiple && HMAG(d) > hmag ) {
! 1518: dpm_ptozp2(d,g,&t,&u); d = t; g = u;
! 1519: hmag = multiple*HMAG(d);
! 1520: }
! 1521: } else {
! 1522: if ( multiple && HMAG(g) > hmag ) {
! 1523: dpm_ptozp(g,&t); g = t;
! 1524: hmag = multiple*HMAG(g);
! 1525: }
! 1526: }
! 1527: }
! 1528: else if ( !full ) {
! 1529: if ( g ) {
! 1530: MKDPM(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1531: }
! 1532: *rp = g; return;
! 1533: } else {
! 1534: m = BDY(g); NEWDMM(mr); mr->dl = m->dl; mr->c = m->c; mr->pos = m->pos;
! 1535: NEXT(mr) = 0; MKDPM(g->nv,mr,t); t->sugar = mr->dl->td;
! 1536: adddpm(CO,d,t,&s); d = s;
! 1537: dpm_rest(g,&t); g = t;
! 1538: }
! 1539: }
! 1540: if ( d )
! 1541: d->sugar = sugar;
! 1542: *rp = d;
! 1543: }
! 1544:
! 1545: /* nf computation over a field */
! 1546:
! 1547: void dp_nf_f(NODE b,DP g,DP *ps,int full,DP *rp)
! 1548: {
! 1549: DP u,p,d,s,t;
! 1550: NODE l;
! 1551: MP m,mr;
! 1552: int i,n;
! 1553: int *wb;
! 1554: int sugar,psugar;
! 1555:
! 1556: if ( !g ) {
! 1557: *rp = 0; return;
! 1558: }
! 1559: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1560: wb = (int *)ALLOCA(n*sizeof(int));
! 1561: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1562: wb[i] = QTOS((Q)BDY(l));
! 1563:
! 1564: sugar = g->sugar;
! 1565: for ( d = 0; g; ) {
! 1566: for ( u = 0, i = 0; i < n; i++ ) {
! 1567: if ( dp_redble(g,p = ps[wb[i]]) ) {
! 1568: dp_red_f(g,p,&u);
! 1569: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1570: sugar = MAX(sugar,psugar);
! 1571: if ( !u ) {
! 1572: if ( d )
! 1573: d->sugar = sugar;
! 1574: *rp = d; return;
! 1575: }
! 1576: break;
! 1577: }
! 1578: }
! 1579: if ( u )
! 1580: g = u;
! 1581: else if ( !full ) {
! 1582: if ( g ) {
! 1583: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1584: }
! 1585: *rp = g; return;
! 1586: } else {
! 1587: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1588: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1589: addd(CO,d,t,&s); d = s;
! 1590: dp_rest(g,&t); g = t;
! 1591: }
! 1592: }
! 1593: if ( d )
! 1594: d->sugar = sugar;
! 1595: *rp = d;
! 1596: }
! 1597:
! 1598: void dpm_nf_f(NODE b,DPM g,DPM *ps,int full,DPM *rp)
! 1599: {
! 1600: DPM u,p,d,s,t;
! 1601: NODE l;
! 1602: DMM m,mr;
! 1603: int i,n;
! 1604: int *wb;
! 1605: int sugar,psugar;
! 1606:
! 1607: if ( !g ) {
! 1608: *rp = 0; return;
! 1609: }
! 1610: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1611: wb = (int *)ALLOCA(n*sizeof(int));
! 1612: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1613: wb[i] = QTOS((Q)BDY(l));
! 1614:
! 1615: sugar = g->sugar;
! 1616: for ( d = 0; g; ) {
! 1617: for ( u = 0, i = 0; i < n; i++ ) {
! 1618: if ( dpm_redble(g,p = ps[wb[i]]) ) {
! 1619: dpm_red_f(g,p,&u);
! 1620: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1621: sugar = MAX(sugar,psugar);
! 1622: if ( !u ) {
! 1623: if ( d )
! 1624: d->sugar = sugar;
! 1625: *rp = d; return;
! 1626: }
! 1627: break;
! 1628: }
! 1629: }
! 1630: if ( u )
! 1631: g = u;
! 1632: else if ( !full ) {
! 1633: if ( g ) {
! 1634: MKDPM(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1635: }
! 1636: *rp = g; return;
! 1637: } else {
! 1638: m = BDY(g); NEWDMM(mr); mr->dl = m->dl; mr->c = m->c; mr->pos = m->pos;
! 1639: NEXT(mr) = 0; MKDPM(g->nv,mr,t); t->sugar = mr->dl->td;
! 1640: adddpm(CO,d,t,&s); d = s;
! 1641: dpm_rest(g,&t); g = t;
! 1642: }
! 1643: }
! 1644: if ( d )
! 1645: d->sugar = sugar;
! 1646: *rp = d;
! 1647: }
! 1648:
! 1649: /* nf computation over GF(mod) (only for internal use) */
! 1650:
! 1651: void dp_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp)
! 1652: {
! 1653: DP u,p,d,s,t;
! 1654: P dmy;
! 1655: NODE l;
! 1656: MP m,mr;
! 1657: int sugar,psugar;
! 1658:
! 1659: if ( !g ) {
! 1660: *rp = 0; return;
! 1661: }
! 1662: sugar = g->sugar;
! 1663: for ( d = 0; g; ) {
! 1664: for ( u = 0, l = b; l; l = NEXT(l) ) {
! 1665: if ( dp_redble(g,p = ps[(long)BDY(l)]) ) {
! 1666: dp_red_mod(d,g,p,mod,&t,&u,&dmy);
! 1667: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1668: sugar = MAX(sugar,psugar);
! 1669: if ( !u ) {
! 1670: if ( d )
! 1671: d->sugar = sugar;
! 1672: *rp = d; return;
! 1673: }
! 1674: d = t;
! 1675: break;
! 1676: }
! 1677: }
! 1678: if ( u )
! 1679: g = u;
! 1680: else if ( !full ) {
! 1681: if ( g ) {
! 1682: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1683: }
! 1684: *rp = g; return;
! 1685: } else {
! 1686: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1687: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1688: addmd(CO,mod,d,t,&s); d = s;
! 1689: dp_rest(g,&t); g = t;
! 1690: }
! 1691: }
! 1692: if ( d )
! 1693: d->sugar = sugar;
! 1694: *rp = d;
! 1695: }
! 1696:
! 1697: void dp_true_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp,P *dnp)
! 1698: {
! 1699: DP u,p,d,s,t;
! 1700: NODE l;
! 1701: MP m,mr;
! 1702: int i,n;
! 1703: int *wb;
! 1704: int sugar,psugar;
! 1705: P dn,tdn,tdn1;
! 1706:
! 1707: dn = (P)ONEM;
! 1708: if ( !g ) {
! 1709: *rp = 0; *dnp = dn; return;
! 1710: }
! 1711: for ( n = 0, l = b; l; l = NEXT(l), n++ );
! 1712: wb = (int *)ALLOCA(n*sizeof(int));
! 1713: for ( i = 0, l = b; i < n; l = NEXT(l), i++ )
! 1714: wb[i] = QTOS((Q)BDY(l));
! 1715: sugar = g->sugar;
! 1716: for ( d = 0; g; ) {
! 1717: for ( u = 0, i = 0; i < n; i++ ) {
! 1718: if ( dp_redble(g,p = ps[wb[i]]) ) {
! 1719: dp_red_mod(d,g,p,mod,&t,&u,&tdn);
! 1720: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1721: sugar = MAX(sugar,psugar);
! 1722: if ( !u ) {
! 1723: if ( d )
! 1724: d->sugar = sugar;
! 1725: *rp = d; *dnp = dn; return;
! 1726: } else {
! 1727: d = t;
! 1728: mulmp(CO,mod,dn,tdn,&tdn1); dn = tdn1;
! 1729: }
! 1730: break;
! 1731: }
! 1732: }
! 1733: if ( u )
! 1734: g = u;
! 1735: else if ( !full ) {
! 1736: if ( g ) {
! 1737: MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
! 1738: }
! 1739: *rp = g; *dnp = dn; return;
! 1740: } else {
! 1741: m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c;
! 1742: NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td;
! 1743: addmd(CO,mod,d,t,&s); d = s;
! 1744: dp_rest(g,&t); g = t;
! 1745: }
! 1746: }
! 1747: if ( d )
! 1748: d->sugar = sugar;
! 1749: *rp = d; *dnp = dn;
! 1750: }
! 1751:
! 1752: void _dp_nf_mod_destructive(NODE b,DP g,DP *ps,int mod,int full,DP *rp)
! 1753: {
! 1754: DP u,p,d;
! 1755: NODE l;
! 1756: MP m,mrd;
! 1757: int sugar,psugar,n,h_reducible;
! 1758:
! 1759: if ( !g ) {
! 1760: *rp = 0; return;
! 1761: }
! 1762: sugar = g->sugar;
! 1763: n = g->nv;
! 1764: for ( d = 0; g; ) {
! 1765: for ( h_reducible = 0, l = b; l; l = NEXT(l) ) {
! 1766: if ( dp_redble(g,p = ps[(long)BDY(l)]) ) {
! 1767: h_reducible = 1;
! 1768: psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
! 1769: _dp_red_mod_destructive(g,p,mod,&u); g = u;
! 1770: sugar = MAX(sugar,psugar);
! 1771: if ( !g ) {
! 1772: if ( d )
! 1773: d->sugar = sugar;
! 1774: _dptodp(d,rp); _free_dp(d); return;
! 1775: }
! 1776: break;
! 1777: }
! 1778: }
! 1779: if ( !h_reducible ) {
! 1780: /* head term is not reducible */
! 1781: if ( !full ) {
! 1782: if ( g )
! 1783: g->sugar = sugar;
! 1784: _dptodp(g,rp); _free_dp(g); return;
! 1785: } else {
! 1786: m = BDY(g);
! 1787: if ( NEXT(m) ) {
! 1788: BDY(g) = NEXT(m); NEXT(m) = 0;
! 1789: } else {
! 1790: _FREEDP(g); g = 0;
! 1791: }
! 1792: if ( d ) {
! 1793: for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) );
! 1794: NEXT(mrd) = m;
! 1795: } else {
! 1796: _MKDP(n,m,d);
! 1797: }
! 1798: }
! 1799: }
! 1800: }
! 1801: if ( d )
! 1802: d->sugar = sugar;
! 1803: _dptodp(d,rp); _free_dp(d);
! 1804: }
! 1805:
! 1806: /* reduction by linear base over a field */
! 1807:
! 1808: void dp_lnf_f(DP p1,DP p2,NODE g,DP *r1p,DP *r2p)
! 1809: {
! 1810: DP r1,r2,b1,b2,t,s;
! 1811: Obj c,c1,c2;
! 1812: NODE l,b;
! 1813: int n;
! 1814:
! 1815: if ( !p1 ) {
! 1816: *r1p = p1; *r2p = p2; return;
! 1817: }
! 1818: n = p1->nv;
! 1819: for ( l = g, r1 = p1, r2 = p2; l; l = NEXT(l) ) {
! 1820: if ( !r1 ) {
! 1821: *r1p = r1; *r2p = r2; return;
! 1822: }
! 1823: b = BDY((LIST)BDY(l)); b1 = (DP)BDY(b);
! 1824: if ( dl_equal(n,BDY(r1)->dl,BDY(b1)->dl) ) {
! 1825: b2 = (DP)BDY(NEXT(b));
! 1826: divr(CO,(Obj)ONE,(Obj)BDY(b1)->c,&c1);
! 1827: mulr(CO,c1,(Obj)BDY(r1)->c,&c2); chsgnr(c2,&c);
! 1828: muldc(CO,b1,(Obj)c,&t); addd(CO,r1,t,&s); r1 = s;
! 1829: muldc(CO,b2,(Obj)c,&t); addd(CO,r2,t,&s); r2 = s;
! 1830: }
! 1831: }
! 1832: *r1p = r1; *r2p = r2;
! 1833: }
! 1834:
! 1835: /* reduction by linear base over GF(mod) */
! 1836:
! 1837: void dp_lnf_mod(DP p1,DP p2,NODE g,int mod,DP *r1p,DP *r2p)
! 1838: {
! 1839: DP r1,r2,b1,b2,t,s;
! 1840: P c;
! 1841: MQ c1,c2;
! 1842: NODE l,b;
! 1843: int n;
! 1844:
! 1845: if ( !p1 ) {
! 1846: *r1p = p1; *r2p = p2; return;
! 1847: }
! 1848: n = p1->nv;
! 1849: for ( l = g, r1 = p1, r2 = p2; l; l = NEXT(l) ) {
! 1850: if ( !r1 ) {
! 1851: *r1p = r1; *r2p = r2; return;
! 1852: }
! 1853: b = BDY((LIST)BDY(l)); b1 = (DP)BDY(b);
! 1854: if ( dl_equal(n,BDY(r1)->dl,BDY(b1)->dl) ) {
! 1855: b2 = (DP)BDY(NEXT(b));
! 1856: invmq(mod,(MQ)BDY(b1)->c,&c1);
! 1857: mulmq(mod,c1,(MQ)BDY(r1)->c,&c2); chsgnmp(mod,(P)c2,&c);
! 1858: mulmdc(CO,mod,b1,c,&t); addmd(CO,mod,r1,t,&s); r1 = s;
! 1859: mulmdc(CO,mod,b2,c,&t); addmd(CO,mod,r2,t,&s); r2 = s;
! 1860: }
! 1861: }
! 1862: *r1p = r1; *r2p = r2;
! 1863: }
! 1864:
! 1865: void dp_nf_tab_mod(DP p,LIST *tab,int mod,DP *rp)
! 1866: {
! 1867: DP s,t,u;
! 1868: MP m;
! 1869: DL h;
! 1870: int i,n;
! 1871:
! 1872: if ( !p ) {
! 1873: *rp = p; return;
! 1874: }
! 1875: n = p->nv;
! 1876: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
! 1877: h = m->dl;
! 1878: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
! 1879: i++;
! 1880: mulmdc(CO,mod,(DP)BDY(NEXT(BDY(tab[i]))),(P)m->c,&t);
! 1881: addmd(CO,mod,s,t,&u); s = u;
! 1882: }
! 1883: *rp = s;
! 1884: }
! 1885:
! 1886: void dp_nf_tab_f(DP p,LIST *tab,DP *rp)
! 1887: {
! 1888: DP s,t,u;
! 1889: MP m;
! 1890: DL h;
! 1891: int i,n;
! 1892:
! 1893: if ( !p ) {
! 1894: *rp = p; return;
! 1895: }
! 1896: n = p->nv;
! 1897: for ( s = 0, i = 0, m = BDY(p); m; m = NEXT(m) ) {
! 1898: h = m->dl;
! 1899: while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) )
! 1900: i++;
! 1901: muldc(CO,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t);
! 1902: addd(CO,s,t,&u); s = u;
! 1903: }
! 1904: *rp = s;
! 1905: }
! 1906:
! 1907: /*
! 1908: * setting flags
! 1909: * call create_order_spec with vl=0 to set old type order.
! 1910: *
! 1911: */
! 1912:
! 1913: int create_order_spec(VL vl,Obj obj,struct order_spec **specp)
! 1914: {
! 1915: int i,j,n,s,row,col,ret,wlen;
! 1916: struct order_spec *spec;
! 1917: struct order_pair *l;
! 1918: Obj wp,wm;
! 1919: NODE node,t,tn,wpair;
! 1920: MAT m;
! 1921: VECT v;
! 1922: pointer **b,*bv;
! 1923: int **w;
! 1924:
! 1925: if ( vl && obj && OID(obj) == O_LIST ) {
! 1926: ret = create_composite_order_spec(vl,(LIST)obj,specp);
! 1927: if ( show_orderspec )
! 1928: print_composite_order_spec(*specp);
! 1929: return ret;
! 1930: }
! 1931:
! 1932: *specp = spec = (struct order_spec *)MALLOC(sizeof(struct order_spec));
! 1933: if ( !obj || NUM(obj) ) {
! 1934: spec->id = 0; spec->obj = obj;
! 1935: spec->ord.simple = QTOS((Q)obj);
! 1936: return 1;
! 1937: } else if ( OID(obj) == O_LIST ) {
! 1938: /* module order; obj = [0|1,w,ord] or [0|1,ord] */
! 1939: node = BDY((LIST)obj);
! 1940: if ( !BDY(node) || NUM(BDY(node)) ) {
! 1941: switch ( length(node) ) {
! 1942: case 2:
! 1943: create_order_spec(0,(Obj)BDY(NEXT(node)),&spec);
! 1944: spec->id += 256; spec->obj = obj;
! 1945: spec->top_weight = 0;
! 1946: spec->module_rank = 0;
! 1947: spec->module_top_weight = 0;
! 1948: spec->ispot = (BDY(node)!=0);
! 1949: if ( spec->ispot ) {
! 1950: n = QTOS((Q)BDY(node));
! 1951: if ( n < 0 )
! 1952: spec->pot_nelim = -n;
! 1953: else
! 1954: spec->pot_nelim = 0;
! 1955: }
! 1956: break;
! 1957:
! 1958: case 3:
! 1959: create_order_spec(0,(Obj)BDY(NEXT(NEXT(node))),&spec);
! 1960: spec->id += 256; spec->obj = obj;
! 1961: spec->ispot = (BDY(node)!=0);
! 1962: node = NEXT(node);
! 1963: if ( !BDY(node) || OID(BDY(node)) != O_LIST )
! 1964: error("create_order_spec : [weight_for_poly,weight_for_modlue] must be specified as a module topweight");
! 1965: wpair = BDY((LIST)BDY(node));
! 1966: if ( length(wpair) != 2 )
! 1967: error("create_order_spec : [weight_for_poly,weight_for_modlue] must be specified as a module topweight");
! 1968:
! 1969: wp = BDY(wpair);
! 1970: wm = BDY(NEXT(wpair));
! 1971: if ( !wp || OID(wp) != O_LIST || !wm || OID(wm) != O_LIST )
! 1972: error("create_order_spec : [weight_for_poly,weight_for_modlue] must be specified as a module topweight");
! 1973: spec->nv = length(BDY((LIST)wp));
! 1974: spec->top_weight = (int *)MALLOC_ATOMIC(spec->nv*sizeof(int));
! 1975: for ( i = 0, t = BDY((LIST)wp); i < spec->nv; t = NEXT(t), i++ )
! 1976: spec->top_weight[i] = QTOS((Q)BDY(t));
! 1977:
! 1978: spec->module_rank = length(BDY((LIST)wm));
! 1979: spec->module_top_weight = (int *)MALLOC_ATOMIC(spec->module_rank*sizeof(int));
! 1980: for ( i = 0, t = BDY((LIST)wm); i < spec->module_rank; t = NEXT(t), i++ )
! 1981: spec->module_top_weight[i] = QTOS((Q)BDY(t));
! 1982: break;
! 1983: default:
! 1984: error("create_order_spec : invalid arguments for module order");
! 1985: }
! 1986:
! 1987: *specp = spec;
! 1988: return 1;
! 1989: } else {
! 1990: /* block order in polynomial ring */
! 1991: for ( n = 0, t = node; t; t = NEXT(t), n++ );
! 1992: l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair));
! 1993: for ( i = 0, t = node, s = 0; i < n; t = NEXT(t), i++ ) {
! 1994: tn = BDY((LIST)BDY(t)); l[i].order = QTOS((Q)BDY(tn));
! 1995: tn = NEXT(tn); l[i].length = QTOS((Q)BDY(tn));
! 1996: s += l[i].length;
! 1997: }
! 1998: spec->id = 1; spec->obj = obj;
! 1999: spec->ord.block.order_pair = l;
! 2000: spec->ord.block.length = n; spec->nv = s;
! 2001: return 1;
! 2002: }
! 2003: } else if ( OID(obj) == O_MAT ) {
! 2004: m = (MAT)obj; row = m->row; col = m->col; b = BDY(m);
! 2005: w = almat(row,col);
! 2006: for ( i = 0; i < row; i++ )
! 2007: for ( j = 0; j < col; j++ )
! 2008: w[i][j] = QTOS((Q)b[i][j]);
! 2009: spec->id = 2; spec->obj = obj;
! 2010: spec->nv = col; spec->ord.matrix.row = row;
! 2011: spec->ord.matrix.matrix = w;
! 2012: return 1;
! 2013: } else
! 2014: return 0;
! 2015: }
! 2016:
! 2017: void print_composite_order_spec(struct order_spec *spec)
! 2018: {
! 2019: int nv,n,len,i,j,k,start;
! 2020: struct weight_or_block *worb;
! 2021:
! 2022: nv = spec->nv;
! 2023: n = spec->ord.composite.length;
! 2024: worb = spec->ord.composite.w_or_b;
! 2025: for ( i = 0; i < n; i++, worb++ ) {
! 2026: len = worb->length;
! 2027: printf("[ ");
! 2028: switch ( worb->type ) {
! 2029: case IS_DENSE_WEIGHT:
! 2030: for ( j = 0; j < len; j++ )
! 2031: printf("%d ",worb->body.dense_weight[j]);
! 2032: for ( ; j < nv; j++ )
! 2033: printf("0 ");
! 2034: break;
! 2035: case IS_SPARSE_WEIGHT:
! 2036: for ( j = 0, k = 0; j < nv; j++ )
! 2037: if ( j == worb->body.sparse_weight[k].pos )
! 2038: printf("%d ",worb->body.sparse_weight[k++].value);
! 2039: else
! 2040: printf("0 ");
! 2041: break;
! 2042: case IS_BLOCK:
! 2043: start = worb->body.block.start;
! 2044: for ( j = 0; j < start; j++ ) printf("0 ");
! 2045: switch ( worb->body.block.order ) {
! 2046: case 0:
! 2047: for ( k = 0; k < len; k++, j++ ) printf("R ");
! 2048: break;
! 2049: case 1:
! 2050: for ( k = 0; k < len; k++, j++ ) printf("G ");
! 2051: break;
! 2052: case 2:
! 2053: for ( k = 0; k < len; k++, j++ ) printf("L ");
! 2054: break;
! 2055: }
! 2056: for ( ; j < nv; j++ ) printf("0 ");
! 2057: break;
! 2058: }
! 2059: printf("]\n");
! 2060: }
! 2061: }
! 2062:
! 2063: struct order_spec *append_block(struct order_spec *spec,
! 2064: int nv,int nalg,int ord)
! 2065: {
! 2066: MAT m,mat;
! 2067: int i,j,row,col,n;
! 2068: Z **b,**wp;
! 2069: int **w;
! 2070: NODE t,s,s0;
! 2071: struct order_pair *l,*l0;
! 2072: int n0,nv0;
! 2073: LIST list0,list1,list;
! 2074: Z oq,nq;
! 2075: struct order_spec *r;
! 2076:
! 2077: r = (struct order_spec *)MALLOC(sizeof(struct order_spec));
! 2078: switch ( spec->id ) {
! 2079: case 0:
! 2080: STOQ(spec->ord.simple,oq); STOQ(nv,nq);
! 2081: t = mknode(2,oq,nq); MKLIST(list0,t);
! 2082: STOQ(ord,oq); STOQ(nalg,nq);
! 2083: t = mknode(2,oq,nq); MKLIST(list1,t);
! 2084: t = mknode(2,list0,list1); MKLIST(list,t);
! 2085: l = (struct order_pair *)MALLOC_ATOMIC(2*sizeof(struct order_pair));
! 2086: l[0].order = spec->ord.simple; l[0].length = nv;
! 2087: l[1].order = ord; l[1].length = nalg;
! 2088: r->id = 1; r->obj = (Obj)list;
! 2089: r->ord.block.order_pair = l;
! 2090: r->ord.block.length = 2;
! 2091: r->nv = nv+nalg;
! 2092: break;
! 2093: case 1:
! 2094: if ( spec->nv != nv )
! 2095: error("append_block : number of variables mismatch");
! 2096: l0 = spec->ord.block.order_pair;
! 2097: n0 = spec->ord.block.length;
! 2098: nv0 = spec->nv;
! 2099: list0 = (LIST)spec->obj;
! 2100: n = n0+1;
! 2101: l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair));
! 2102: for ( i = 0; i < n0; i++ )
! 2103: l[i] = l0[i];
! 2104: l[i].order = ord; l[i].length = nalg;
! 2105: for ( t = BDY(list0), s0 = 0; t; t = NEXT(t) ) {
! 2106: NEXTNODE(s0,s); BDY(s) = BDY(t);
! 2107: }
! 2108: STOQ(ord,oq); STOQ(nalg,nq);
! 2109: t = mknode(2,oq,nq); MKLIST(list,t);
! 2110: NEXTNODE(s0,s); BDY(s) = (pointer)list; NEXT(s) = 0;
! 2111: MKLIST(list,s0);
! 2112: r->id = 1; r->obj = (Obj)list;
! 2113: r->ord.block.order_pair = l;
! 2114: r->ord.block.length = n;
! 2115: r->nv = nv+nalg;
! 2116: break;
! 2117: case 2:
! 2118: if ( spec->nv != nv )
! 2119: error("append_block : number of variables mismatch");
! 2120: m = (MAT)spec->obj;
! 2121: row = m->row; col = m->col; b = (Z **)BDY(m);
! 2122: w = almat(row+nalg,col+nalg);
! 2123: MKMAT(mat,row+nalg,col+nalg); wp = (Z **)BDY(mat);
! 2124: for ( i = 0; i < row; i++ )
! 2125: for ( j = 0; j < col; j++ ) {
! 2126: w[i][j] = QTOS(b[i][j]);
! 2127: wp[i][j] = b[i][j];
! 2128: }
! 2129: for ( i = 0; i < nalg; i++ ) {
! 2130: w[i+row][i+col] = 1;
! 2131: wp[i+row][i+col] = ONE;
! 2132: }
! 2133: r->id = 2; r->obj = (Obj)mat;
! 2134: r->nv = col+nalg; r->ord.matrix.row = row+nalg;
! 2135: r->ord.matrix.matrix = w;
! 2136: break;
! 2137: case 3:
! 2138: default:
! 2139: /* XXX */
! 2140: error("append_block : not implemented yet");
! 2141: }
! 2142: return r;
! 2143: }
! 2144:
! 2145: int comp_sw(struct sparse_weight *a, struct sparse_weight *b)
! 2146: {
! 2147: if ( a->pos > b->pos ) return 1;
! 2148: else if ( a->pos < b->pos ) return -1;
! 2149: else return 0;
! 2150: }
! 2151:
! 2152: /* order = [w_or_b, w_or_b, ... ] */
! 2153: /* w_or_b = w or b */
! 2154: /* w = [1,2,...] or [x,1,y,2,...] */
! 2155: /* b = [@lex,x,y,...,z] etc */
! 2156:
! 2157: int create_composite_order_spec(VL vl,LIST order,struct order_spec **specp)
! 2158: {
! 2159: NODE wb,t,p;
! 2160: struct order_spec *spec;
! 2161: VL tvl;
! 2162: int n,i,j,k,l,start,end,len,w;
! 2163: int *dw;
! 2164: struct sparse_weight *sw;
! 2165: struct weight_or_block *w_or_b;
! 2166: Obj a0;
! 2167: NODE a;
! 2168: V v,sv,ev;
! 2169: SYMBOL sym;
! 2170: int *top;
! 2171:
! 2172: /* l = number of vars in vl */
! 2173: for ( l = 0, tvl = vl; tvl; tvl = NEXT(tvl), l++ );
! 2174: /* n = number of primitives in order */
! 2175: wb = BDY(order);
! 2176: n = length(wb);
! 2177: *specp = spec = (struct order_spec *)MALLOC(sizeof(struct order_spec));
! 2178: spec->id = 3;
! 2179: spec->obj = (Obj)order;
! 2180: spec->nv = l;
! 2181: spec->ord.composite.length = n;
! 2182: w_or_b = spec->ord.composite.w_or_b = (struct weight_or_block *)
! 2183: MALLOC(sizeof(struct weight_or_block)*(n+1));
! 2184:
! 2185: /* top : register the top variable in each w_or_b specification */
! 2186: top = (int *)ALLOCA(l*sizeof(int));
! 2187: for ( i = 0; i < l; i++ ) top[i] = 0;
! 2188:
! 2189: for ( t = wb, i = 0; t; t = NEXT(t), i++ ) {
! 2190: if ( !BDY(t) || OID((Obj)BDY(t)) != O_LIST )
! 2191: error("a list of lists must be specified for the key \"order\"");
! 2192: a = BDY((LIST)BDY(t));
! 2193: len = length(a);
! 2194: a0 = (Obj)BDY(a);
! 2195: if ( !a0 || OID(a0) == O_N ) {
! 2196: /* a is a dense weight vector */
! 2197: dw = (int *)MALLOC(sizeof(int)*len);
! 2198: for ( j = 0, p = a; j < len; p = NEXT(p), j++ ) {
! 2199: if ( !INT((Q)BDY(p)) )
! 2200: error("a dense weight vector must be specified as a list of integers");
! 2201: dw[j] = QTOS((Q)BDY(p));
! 2202: }
! 2203: w_or_b[i].type = IS_DENSE_WEIGHT;
! 2204: w_or_b[i].length = len;
! 2205: w_or_b[i].body.dense_weight = dw;
! 2206:
! 2207: /* find the top */
! 2208: for ( k = 0; k < len && !dw[k]; k++ );
! 2209: if ( k < len ) top[k] = 1;
! 2210:
! 2211: } else if ( OID(a0) == O_P ) {
! 2212: /* a is a sparse weight vector */
! 2213: len >>= 1;
! 2214: sw = (struct sparse_weight *)
! 2215: MALLOC(sizeof(struct sparse_weight)*len);
! 2216: for ( j = 0, p = a; j < len; j++ ) {
! 2217: if ( !BDY(p) || OID((P)BDY(p)) != O_P )
! 2218: error("a sparse weight vector must be specified as [var1,weight1,...]");
! 2219: v = VR((P)BDY(p)); p = NEXT(p);
! 2220: for ( tvl = vl, k = 0; tvl && tvl->v != v;
! 2221: k++, tvl = NEXT(tvl) );
! 2222: if ( !tvl )
! 2223: error("invalid variable name in a sparse weight vector");
! 2224: sw[j].pos = k;
! 2225: if ( !INT((Q)BDY(p)) )
! 2226: error("a sparse weight vector must be specified as [var1,weight1,...]");
! 2227: sw[j].value = QTOS((Q)BDY(p)); p = NEXT(p);
! 2228: }
! 2229: qsort(sw,len,sizeof(struct sparse_weight),
! 2230: (int (*)(const void *,const void *))comp_sw);
! 2231: w_or_b[i].type = IS_SPARSE_WEIGHT;
! 2232: w_or_b[i].length = len;
! 2233: w_or_b[i].body.sparse_weight = sw;
! 2234:
! 2235: /* find the top */
! 2236: for ( k = 0; k < len && !sw[k].value; k++ );
! 2237: if ( k < len ) top[sw[k].pos] = 1;
! 2238: } else if ( OID(a0) == O_RANGE ) {
! 2239: /* [range(v1,v2),w] */
! 2240: sv = VR((P)(((RANGE)a0)->start));
! 2241: ev = VR((P)(((RANGE)a0)->end));
! 2242: for ( tvl = vl, start = 0; tvl && tvl->v != sv; start++, tvl = NEXT(tvl) );
! 2243: if ( !tvl )
! 2244: error("invalid range");
! 2245: for ( end = start; tvl && tvl->v != ev; end++, tvl = NEXT(tvl) );
! 2246: if ( !tvl )
! 2247: error("invalid range");
! 2248: len = end-start+1;
! 2249: sw = (struct sparse_weight *)
! 2250: MALLOC(sizeof(struct sparse_weight)*len);
! 2251: w = QTOS((Q)BDY(NEXT(a)));
! 2252: for ( tvl = vl, k = 0; k < start; k++, tvl = NEXT(tvl) );
! 2253: for ( j = 0 ; k <= end; k++, tvl = NEXT(tvl), j++ ) {
! 2254: sw[j].pos = k;
! 2255: sw[j].value = w;
! 2256: }
! 2257: w_or_b[i].type = IS_SPARSE_WEIGHT;
! 2258: w_or_b[i].length = len;
! 2259: w_or_b[i].body.sparse_weight = sw;
! 2260:
! 2261: /* register the top */
! 2262: if ( w ) top[start] = 1;
! 2263: } else if ( OID(a0) == O_SYMBOL ) {
! 2264: /* a is a block */
! 2265: sym = (SYMBOL)a0; a = NEXT(a); len--;
! 2266: if ( OID((Obj)BDY(a)) == O_RANGE ) {
! 2267: sv = VR((P)(((RANGE)BDY(a))->start));
! 2268: ev = VR((P)(((RANGE)BDY(a))->end));
! 2269: for ( tvl = vl, start = 0; tvl && tvl->v != sv; start++, tvl = NEXT(tvl) );
! 2270: if ( !tvl )
! 2271: error("invalid range");
! 2272: for ( end = start; tvl && tvl->v != ev; end++, tvl = NEXT(tvl) );
! 2273: if ( !tvl )
! 2274: error("invalid range");
! 2275: len = end-start+1;
! 2276: } else {
! 2277: for ( start = 0, tvl = vl; tvl->v != VR((P)BDY(a));
! 2278: tvl = NEXT(tvl), start++ );
! 2279: for ( p = NEXT(a), tvl = NEXT(tvl); p;
! 2280: p = NEXT(p), tvl = NEXT(tvl) ) {
! 2281: if ( !BDY(p) || OID((P)BDY(p)) != O_P )
! 2282: error("a block must be specified as [ordsymbol,var1,var2,...]");
! 2283: if ( tvl->v != VR((P)BDY(p)) ) break;
! 2284: }
! 2285: if ( p )
! 2286: error("a block must be contiguous in the variable list");
! 2287: }
! 2288: w_or_b[i].type = IS_BLOCK;
! 2289: w_or_b[i].length = len;
! 2290: w_or_b[i].body.block.start = start;
! 2291: if ( !strcmp(sym->name,"@grlex") )
! 2292: w_or_b[i].body.block.order = 0;
! 2293: else if ( !strcmp(sym->name,"@glex") )
! 2294: w_or_b[i].body.block.order = 1;
! 2295: else if ( !strcmp(sym->name,"@lex") )
! 2296: w_or_b[i].body.block.order = 2;
! 2297: else
! 2298: error("invalid ordername");
! 2299: /* register the tops */
! 2300: for ( j = 0, k = start; j < len; j++, k++ )
! 2301: top[k] = 1;
! 2302: }
! 2303: }
! 2304: for ( k = 0; k < l && top[k]; k++ );
! 2305: if ( k < l ) {
! 2306: /* incomplete order specification; add @grlex */
! 2307: w_or_b[n].type = IS_BLOCK;
! 2308: w_or_b[n].length = l;
! 2309: w_or_b[n].body.block.start = 0;
! 2310: w_or_b[n].body.block.order = 0;
! 2311: spec->ord.composite.length = n+1;
! 2312: }
! 2313: }
! 2314:
! 2315: /* module order spec */
! 2316:
! 2317: void create_modorder_spec(int id,LIST shift,struct modorder_spec **s)
! 2318: {
! 2319: struct modorder_spec *spec;
! 2320: NODE n,t;
! 2321: LIST list;
! 2322: int *ds;
! 2323: int i,l;
! 2324: Z q;
! 2325:
! 2326: *s = spec = (struct modorder_spec *)MALLOC(sizeof(struct modorder_spec));
! 2327: spec->id = id;
! 2328: if ( shift ) {
! 2329: n = BDY(shift);
! 2330: spec->len = l = length(n);
! 2331: spec->degree_shift = ds = (int *)MALLOC_ATOMIC(l*sizeof(int));
! 2332: for ( t = n, i = 0; t; t = NEXT(t), i++ )
! 2333: ds[i] = QTOS((Q)BDY(t));
! 2334: } else {
! 2335: spec->len = 0;
! 2336: spec->degree_shift = 0;
! 2337: }
! 2338: STOQ(id,q);
! 2339: n = mknode(2,q,shift);
! 2340: MKLIST(list,n);
! 2341: spec->obj = (Obj)list;
! 2342: }
! 2343:
! 2344: /*
! 2345: * converters
! 2346: *
! 2347: */
! 2348:
! 2349: void dp_homo(DP p,DP *rp)
! 2350: {
! 2351: MP m,mr,mr0;
! 2352: int i,n,nv,td;
! 2353: DL dl,dlh;
! 2354:
! 2355: if ( !p )
! 2356: *rp = 0;
! 2357: else {
! 2358: n = p->nv; nv = n + 1;
! 2359: m = BDY(p); td = sugard(m);
! 2360: for ( mr0 = 0; m; m = NEXT(m) ) {
! 2361: NEXTMP(mr0,mr); mr->c = m->c;
! 2362: dl = m->dl;
! 2363: mr->dl = dlh = (DL)MALLOC_ATOMIC((nv+1)*sizeof(int));
! 2364: dlh->td = td;
! 2365: for ( i = 0; i < n; i++ )
! 2366: dlh->d[i] = dl->d[i];
! 2367: dlh->d[n] = td - dl->td;
! 2368: }
! 2369: NEXT(mr) = 0; MKDP(nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 2370: }
! 2371: }
! 2372:
! 2373: void dp_dehomo(DP p,DP *rp)
! 2374: {
! 2375: MP m,mr,mr0;
! 2376: int i,n,nv;
! 2377: DL dl,dlh;
! 2378:
! 2379: if ( !p )
! 2380: *rp = 0;
! 2381: else {
! 2382: n = p->nv; nv = n - 1;
! 2383: m = BDY(p);
! 2384: for ( mr0 = 0; m; m = NEXT(m) ) {
! 2385: NEXTMP(mr0,mr); mr->c = m->c;
! 2386: dlh = m->dl;
! 2387: mr->dl = dl = (DL)MALLOC_ATOMIC((nv+1)*sizeof(int));
! 2388: dl->td = dlh->td - dlh->d[nv];
! 2389: for ( i = 0; i < nv; i++ )
! 2390: dl->d[i] = dlh->d[i];
! 2391: }
! 2392: NEXT(mr) = 0; MKDP(nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 2393: }
! 2394: }
! 2395:
! 2396: void dp_mod(DP p,int mod,NODE subst,DP *rp)
! 2397: {
! 2398: MP m,mr,mr0;
! 2399: P t,s,s1;
! 2400: V v;
! 2401: NODE tn;
! 2402:
! 2403: if ( !p )
! 2404: *rp = 0;
! 2405: else {
! 2406: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 2407: for ( tn = subst, s = (P)m->c; tn; tn = NEXT(tn) ) {
! 2408: v = VR((P)BDY(tn)); tn = NEXT(tn);
! 2409: substp(CO,s,v,(P)BDY(tn),&s1); s = s1;
! 2410: }
! 2411: ptomp(mod,s,&t);
! 2412: if ( t ) {
! 2413: NEXTMP(mr0,mr); mr->c = (Obj)t; mr->dl = m->dl;
! 2414: }
! 2415: }
! 2416: if ( mr0 ) {
! 2417: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 2418: } else
! 2419: *rp = 0;
! 2420: }
! 2421: }
! 2422:
! 2423: void dp_rat(DP p,DP *rp)
! 2424: {
! 2425: MP m,mr,mr0;
! 2426:
! 2427: if ( !p )
! 2428: *rp = 0;
! 2429: else {
! 2430: for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) {
! 2431: NEXTMP(mr0,mr); mptop((P)m->c,(P *)&mr->c); mr->dl = m->dl;
! 2432: }
! 2433: if ( mr0 ) {
! 2434: NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar;
! 2435: } else
! 2436: *rp = 0;
! 2437: }
! 2438: }
! 2439:
! 2440:
! 2441: void homogenize_order(struct order_spec *old,int n,struct order_spec **newp)
! 2442: {
! 2443: struct order_pair *l;
! 2444: int length,nv,row,i,j;
! 2445: int **newm,**oldm;
! 2446: struct order_spec *new;
! 2447: int onv,nnv,nlen,olen,owlen;
! 2448: struct weight_or_block *owb,*nwb;
! 2449:
! 2450: *newp = new = (struct order_spec *)MALLOC(sizeof(struct order_spec));
! 2451: switch ( old->id ) {
! 2452: case 0:
! 2453: switch ( old->ord.simple ) {
! 2454: case 0:
! 2455: new->id = 0; new->ord.simple = 0; break;
! 2456: case 1:
! 2457: l = (struct order_pair *)
! 2458: MALLOC_ATOMIC(2*sizeof(struct order_pair));
! 2459: l[0].length = n; l[0].order = 1;
! 2460: l[1].length = 1; l[1].order = 2;
! 2461: new->id = 1;
! 2462: new->ord.block.order_pair = l;
! 2463: new->ord.block.length = 2; new->nv = n+1;
! 2464: break;
! 2465: case 2:
! 2466: new->id = 0; new->ord.simple = 1; break;
! 2467: case 3: case 4: case 5:
! 2468: new->id = 0; new->ord.simple = old->ord.simple+3;
! 2469: dp_nelim = n-1; break;
! 2470: case 6: case 7: case 8: case 9:
! 2471: new->id = 0; new->ord.simple = old->ord.simple; break;
! 2472: default:
! 2473: error("homogenize_order : invalid input");
! 2474: }
! 2475: break;
! 2476: case 1: case 257:
! 2477: length = old->ord.block.length;
! 2478: l = (struct order_pair *)
! 2479: MALLOC_ATOMIC((length+1)*sizeof(struct order_pair));
! 2480: bcopy((char *)old->ord.block.order_pair,(char *)l,length*sizeof(struct order_pair));
! 2481: l[length].order = 2; l[length].length = 1;
! 2482: new->id = old->id; new->nv = n+1;
! 2483: new->ord.block.order_pair = l;
! 2484: new->ord.block.length = length+1;
! 2485: new->ispot = old->ispot;
! 2486: break;
! 2487: case 2: case 258:
! 2488: nv = old->nv; row = old->ord.matrix.row;
! 2489: oldm = old->ord.matrix.matrix; newm = almat(row+1,nv+1);
! 2490: for ( i = 0; i <= nv; i++ )
! 2491: newm[0][i] = 1;
! 2492: for ( i = 0; i < row; i++ ) {
! 2493: for ( j = 0; j < nv; j++ )
! 2494: newm[i+1][j] = oldm[i][j];
! 2495: newm[i+1][j] = 0;
! 2496: }
! 2497: new->id = old->id; new->nv = nv+1;
! 2498: new->ord.matrix.row = row+1; new->ord.matrix.matrix = newm;
! 2499: new->ispot = old->ispot;
! 2500: break;
! 2501: case 3: case 259:
! 2502: onv = old->nv;
! 2503: nnv = onv+1;
! 2504: olen = old->ord.composite.length;
! 2505: nlen = olen+1;
! 2506: owb = old->ord.composite.w_or_b;
! 2507: nwb = (struct weight_or_block *)
! 2508: MALLOC(nlen*sizeof(struct weight_or_block));
! 2509: for ( i = 0; i < olen; i++ ) {
! 2510: nwb[i].type = owb[i].type;
! 2511: switch ( owb[i].type ) {
! 2512: case IS_DENSE_WEIGHT:
! 2513: owlen = owb[i].length;
! 2514: nwb[i].length = owlen+1;
! 2515: nwb[i].body.dense_weight = (int *)MALLOC((owlen+1)*sizeof(int));
! 2516: for ( j = 0; j < owlen; j++ )
! 2517: nwb[i].body.dense_weight[j] = owb[i].body.dense_weight[j];
! 2518: nwb[i].body.dense_weight[owlen] = 0;
! 2519: break;
! 2520: case IS_SPARSE_WEIGHT:
! 2521: nwb[i].length = owb[i].length;
! 2522: nwb[i].body.sparse_weight = owb[i].body.sparse_weight;
! 2523: break;
! 2524: case IS_BLOCK:
! 2525: nwb[i].length = owb[i].length;
! 2526: nwb[i].body.block = owb[i].body.block;
! 2527: break;
! 2528: }
! 2529: }
! 2530: nwb[i].type = IS_SPARSE_WEIGHT;
! 2531: nwb[i].body.sparse_weight =
! 2532: (struct sparse_weight *)MALLOC(sizeof(struct sparse_weight));
! 2533: nwb[i].body.sparse_weight[0].pos = onv;
! 2534: nwb[i].body.sparse_weight[0].value = 1;
! 2535: new->id = old->id;
! 2536: new->nv = nnv;
! 2537: new->ord.composite.length = nlen;
! 2538: new->ord.composite.w_or_b = nwb;
! 2539: new->ispot = old->ispot;
! 2540: print_composite_order_spec(new);
! 2541: break;
! 2542: case 256: /* simple module order */
! 2543: switch ( old->ord.simple ) {
! 2544: case 0:
! 2545: new->id = 256; new->ord.simple = 0; break;
! 2546: case 1:
! 2547: l = (struct order_pair *)
! 2548: MALLOC_ATOMIC(2*sizeof(struct order_pair));
! 2549: l[0].length = n; l[0].order = old->ord.simple;
! 2550: l[1].length = 1; l[1].order = 2;
! 2551: new->id = 257;
! 2552: new->ord.block.order_pair = l;
! 2553: new->ord.block.length = 2; new->nv = n+1;
! 2554: break;
! 2555: case 2:
! 2556: new->id = 256; new->ord.simple = 1; break;
! 2557: default:
! 2558: error("homogenize_order : invalid input");
! 2559: }
! 2560: new->ispot = old->ispot;
! 2561: break;
! 2562: default:
! 2563: error("homogenize_order : invalid input");
! 2564: }
! 2565: }
! 2566:
! 2567: int comp_nm(Q *a,Q *b)
! 2568: {
! 2569: Z z,nma,nmb;
! 2570:
! 2571: nmq(*a,&z); absz(z,&nma);
! 2572: nmq(*b,&z); absz(z,&nmb);
! 2573: return cmpz(nma,nmb);
! 2574: }
! 2575:
! 2576: void sortbynm(Q *w,int n)
! 2577: {
! 2578: qsort(w,n,sizeof(Q),(int (*)(const void *,const void *))comp_nm);
! 2579: }
! 2580:
! 2581:
! 2582: /*
! 2583: * simple operations
! 2584: *
! 2585: */
! 2586:
! 2587: int dp_redble(DP p1,DP p2)
! 2588: {
! 2589: int i,n;
! 2590: DL d1,d2;
! 2591:
! 2592: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 2593: if ( d1->td < d2->td )
! 2594: return 0;
! 2595: else {
! 2596: for ( i = 0, n = p1->nv; i < n; i++ )
! 2597: if ( d1->d[i] < d2->d[i] )
! 2598: return 0;
! 2599: return 1;
! 2600: }
! 2601: }
! 2602:
! 2603: int dpm_redble(DPM p1,DPM p2)
! 2604: {
! 2605: int i,n;
! 2606: DL d1,d2;
! 2607:
! 2608: if ( BDY(p1)->pos != BDY(p2)->pos ) return 0;
! 2609: d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 2610: if ( d1->td < d2->td )
! 2611: return 0;
! 2612: else {
! 2613: for ( i = 0, n = p1->nv; i < n; i++ )
! 2614: if ( d1->d[i] < d2->d[i] )
! 2615: return 0;
! 2616: return 1;
! 2617: }
! 2618: }
! 2619:
! 2620:
! 2621: void dp_subd(DP p1,DP p2,DP *rp)
! 2622: {
! 2623: int i,n;
! 2624: DL d1,d2,d;
! 2625: MP m;
! 2626: DP s;
! 2627:
! 2628: n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl;
! 2629: NEWDL(d,n); d->td = d1->td - d2->td;
! 2630: for ( i = 0; i < n; i++ )
! 2631: d->d[i] = d1->d[i]-d2->d[i];
! 2632: NEWMP(m); m->dl = d; m->c = (Obj)ONE; NEXT(m) = 0;
! 2633: MKDP(n,m,s); s->sugar = d->td;
! 2634: *rp = s;
! 2635: }
! 2636:
! 2637: void dltod(DL d,int n,DP *rp)
! 2638: {
! 2639: MP m;
! 2640: DP s;
! 2641:
! 2642: NEWMP(m); m->dl = d; m->c = (Obj)ONE; NEXT(m) = 0;
! 2643: MKDP(n,m,s); s->sugar = d->td;
! 2644: *rp = s;
! 2645: }
! 2646:
! 2647: void dp_hm(DP p,DP *rp)
! 2648: {
! 2649: MP m,mr;
! 2650:
! 2651: if ( !p )
! 2652: *rp = 0;
! 2653: else {
! 2654: m = BDY(p);
! 2655: NEWMP(mr); mr->dl = m->dl; mr->c = m->c; NEXT(mr) = 0;
! 2656: MKDP(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
! 2657: }
! 2658: }
! 2659:
! 2660: void dp_ht(DP p,DP *rp)
! 2661: {
! 2662: MP m,mr;
! 2663:
! 2664: if ( !p )
! 2665: *rp = 0;
! 2666: else {
! 2667: m = BDY(p);
! 2668: NEWMP(mr); mr->dl = m->dl; mr->c = (Obj)ONE; NEXT(mr) = 0;
! 2669: MKDP(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
! 2670: }
! 2671: }
! 2672:
! 2673: void dpm_hm(DPM p,DPM *rp)
! 2674: {
! 2675: DMM m,mr;
! 2676:
! 2677: if ( !p )
! 2678: *rp = 0;
! 2679: else {
! 2680: m = BDY(p);
! 2681: NEWDMM(mr); mr->dl = m->dl; mr->c = m->c; mr->pos = m->pos; NEXT(mr) = 0;
! 2682: MKDPM(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
! 2683: }
! 2684: }
! 2685:
! 2686: void dpm_ht(DPM p,DPM *rp)
! 2687: {
! 2688: DMM m,mr;
! 2689:
! 2690: if ( !p )
! 2691: *rp = 0;
! 2692: else {
! 2693: m = BDY(p);
! 2694: NEWDMM(mr); mr->dl = m->dl; mr->pos = m->pos; mr->c = (Obj)ONE; NEXT(mr) = 0;
! 2695: MKDPM(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */
! 2696: }
! 2697: }
! 2698:
! 2699:
! 2700: void dp_rest(DP p,DP *rp)
! 2701: {
! 2702: MP m;
! 2703:
! 2704: m = BDY(p);
! 2705: if ( !NEXT(m) )
! 2706: *rp = 0;
! 2707: else {
! 2708: MKDP(p->nv,NEXT(m),*rp);
! 2709: if ( *rp )
! 2710: (*rp)->sugar = p->sugar;
! 2711: }
! 2712: }
! 2713:
! 2714: void dpm_rest(DPM p,DPM *rp)
! 2715: {
! 2716: DMM m;
! 2717:
! 2718: m = BDY(p);
! 2719: if ( !NEXT(m) )
! 2720: *rp = 0;
! 2721: else {
! 2722: MKDPM(p->nv,NEXT(m),*rp);
! 2723: if ( *rp )
! 2724: (*rp)->sugar = p->sugar;
! 2725: }
! 2726: }
! 2727:
! 2728: DL lcm_of_DL(int nv,DL dl1,DL dl2,DL dl)
! 2729: {
! 2730: register int i, *d1, *d2, *d, td;
! 2731:
! 2732: if ( !dl ) NEWDL(dl,nv);
! 2733: d = dl->d, d1 = dl1->d, d2 = dl2->d;
! 2734: for ( td = 0, i = 0; i < nv; d1++, d2++, d++, i++ ) {
! 2735: *d = *d1 > *d2 ? *d1 : *d2;
! 2736: td += MUL_WEIGHT(*d,i);
! 2737: }
! 2738: dl->td = td;
! 2739: return dl;
! 2740: }
! 2741:
! 2742: int dl_equal(int nv,DL dl1,DL dl2)
! 2743: {
! 2744: register int *d1, *d2, n;
! 2745:
! 2746: if ( dl1->td != dl2->td ) return 0;
! 2747: for ( d1 = dl1->d, d2 = dl2->d, n = nv; --n >= 0; d1++, d2++ )
! 2748: if ( *d1 != *d2 ) return 0;
! 2749: return 1;
! 2750: }
! 2751:
! 2752: int dp_nt(DP p)
! 2753: {
! 2754: int i;
! 2755: MP m;
! 2756:
! 2757: if ( !p )
! 2758: return 0;
! 2759: else {
! 2760: for ( i = 0, m = BDY(p); m; m = NEXT(m), i++ );
! 2761: return i;
! 2762: }
! 2763: }
! 2764:
! 2765: int dp_homogeneous(DP p)
! 2766: {
! 2767: MP m;
! 2768: int d;
! 2769:
! 2770: if ( !p )
! 2771: return 1;
! 2772: else {
! 2773: m = BDY(p);
! 2774: d = m->dl->td;
! 2775: m = NEXT(m);
! 2776: for ( ; m; m = NEXT(m) ) {
! 2777: if ( m->dl->td != d )
! 2778: return 0;
! 2779: }
! 2780: return 1;
! 2781: }
! 2782: }
! 2783:
! 2784: void _print_mp(int nv,MP m)
! 2785: {
! 2786: int i;
! 2787:
! 2788: if ( !m )
! 2789: return;
! 2790: for ( ; m; m = NEXT(m) ) {
! 2791: fprintf(stderr,"%d<",ITOS(C(m)));
! 2792: for ( i = 0; i < nv; i++ ) {
! 2793: fprintf(stderr,"%d",m->dl->d[i]);
! 2794: if ( i != nv-1 )
! 2795: fprintf(stderr," ");
! 2796: }
! 2797: fprintf(stderr,">",C(m));
! 2798: }
! 2799: fprintf(stderr,"\n");
! 2800: }
! 2801:
! 2802: static int cmp_mp_nvar;
! 2803:
! 2804: int comp_mp(MP *a,MP *b)
! 2805: {
! 2806: return -(*cmpdl)(cmp_mp_nvar,(*a)->dl,(*b)->dl);
! 2807: }
! 2808:
! 2809: void dp_sort(DP p,DP *rp)
! 2810: {
! 2811: MP t,mp,mp0;
! 2812: int i,n;
! 2813: DP r;
! 2814: MP *w;
! 2815:
! 2816: if ( !p ) {
! 2817: *rp = 0;
! 2818: return;
! 2819: }
! 2820: for ( t = BDY(p), n = 0; t; t = NEXT(t), n++ );
! 2821: w = (MP *)ALLOCA(n*sizeof(MP));
! 2822: for ( t = BDY(p), i = 0; i < n; t = NEXT(t), i++ )
! 2823: w[i] = t;
! 2824: cmp_mp_nvar = NV(p);
! 2825: qsort(w,n,sizeof(MP),(int (*)(const void *,const void *))comp_mp);
! 2826: mp0 = 0;
! 2827: for ( i = n-1; i >= 0; i-- ) {
! 2828: NEWMP(mp); mp->dl = w[i]->dl; C(mp) = C(w[i]);
! 2829: NEXT(mp) = mp0; mp0 = mp;
! 2830: }
! 2831: MKDP(p->nv,mp0,r);
! 2832: r->sugar = p->sugar;
! 2833: *rp = r;
! 2834: }
! 2835:
! 2836: DP extract_initial_term_from_dp(DP p,int *weight,int n);
! 2837: LIST extract_initial_term(LIST f,int *weight,int n);
! 2838:
! 2839: DP extract_initial_term_from_dp(DP p,int *weight,int n)
! 2840: {
! 2841: int w,t,i,top;
! 2842: MP m,r0,r;
! 2843: DP dp;
! 2844:
! 2845: if ( !p ) return 0;
! 2846: top = 1;
! 2847: for ( m = BDY(p); m; m = NEXT(m) ) {
! 2848: for ( i = 0, t = 0; i < n; i++ )
! 2849: t += weight[i]*m->dl->d[i];
! 2850: if ( top || t > w ) {
! 2851: r0 = 0;
! 2852: w = t;
! 2853: top = 0;
! 2854: }
! 2855: if ( t == w ) {
! 2856: NEXTMP(r0,r);
! 2857: r->dl = m->dl;
! 2858: r->c = m->c;
! 2859: }
! 2860: }
! 2861: NEXT(r) = 0;
! 2862: MKDP(p->nv,r0,dp);
! 2863: return dp;
! 2864: }
! 2865:
! 2866: LIST extract_initial_term(LIST f,int *weight,int n)
! 2867: {
! 2868: NODE nd,r0,r;
! 2869: Obj p;
! 2870: LIST l;
! 2871:
! 2872: nd = BDY(f);
! 2873: for ( r0 = 0; nd; nd = NEXT(nd) ) {
! 2874: NEXTNODE(r0,r);
! 2875: p = (Obj)BDY(nd);
! 2876: BDY(r) = (pointer)extract_initial_term_from_dp((DP)p,weight,n);
! 2877: }
! 2878: if ( r0 ) NEXT(r) = 0;
! 2879: MKLIST(l,r0);
! 2880: return l;
! 2881: }
! 2882:
! 2883: LIST dp_initial_term(LIST f,struct order_spec *ord)
! 2884: {
! 2885: int n,l,i;
! 2886: struct weight_or_block *worb;
! 2887: int *weight;
! 2888:
! 2889: switch ( ord->id ) {
! 2890: case 2: /* matrix order */
! 2891: /* extract the first row */
! 2892: n = ord->nv;
! 2893: weight = ord->ord.matrix.matrix[0];
! 2894: return extract_initial_term(f,weight,n);
! 2895: case 3: /* composite order */
! 2896: /* the first w_or_b */
! 2897: worb = ord->ord.composite.w_or_b;
! 2898: switch ( worb->type ) {
! 2899: case IS_DENSE_WEIGHT:
! 2900: n = worb->length;
! 2901: weight = worb->body.dense_weight;
! 2902: return extract_initial_term(f,weight,n);
! 2903: case IS_SPARSE_WEIGHT:
! 2904: n = ord->nv;
! 2905: weight = (int *)ALLOCA(n*sizeof(int));
! 2906: for ( i = 0; i < n; i++ ) weight[i] = 0;
! 2907: l = worb->length;
! 2908: for ( i = 0; i < l; i++ )
! 2909: weight[worb->body.sparse_weight[i].pos]
! 2910: = worb->body.sparse_weight[i].value;
! 2911: return extract_initial_term(f,weight,n);
! 2912: default:
! 2913: error("dp_initial_term : unsupported order");
! 2914: }
! 2915: default:
! 2916: error("dp_initial_term : unsupported order");
! 2917: }
! 2918: }
! 2919:
! 2920: int highest_order_dp(DP p,int *weight,int n);
! 2921: LIST highest_order(LIST f,int *weight,int n);
! 2922:
! 2923: int highest_order_dp(DP p,int *weight,int n)
! 2924: {
! 2925: int w,t,i,top;
! 2926: MP m;
! 2927:
! 2928: if ( !p ) return -1;
! 2929: top = 1;
! 2930: for ( m = BDY(p); m; m = NEXT(m) ) {
! 2931: for ( i = 0, t = 0; i < n; i++ )
! 2932: t += weight[i]*m->dl->d[i];
! 2933: if ( top || t > w ) {
! 2934: w = t;
! 2935: top = 0;
! 2936: }
! 2937: }
! 2938: return w;
! 2939: }
! 2940:
! 2941: LIST highest_order(LIST f,int *weight,int n)
! 2942: {
! 2943: int h;
! 2944: NODE nd,r0,r;
! 2945: Obj p;
! 2946: LIST l;
! 2947: Z q;
! 2948:
! 2949: nd = BDY(f);
! 2950: for ( r0 = 0; nd; nd = NEXT(nd) ) {
! 2951: NEXTNODE(r0,r);
! 2952: p = (Obj)BDY(nd);
! 2953: h = highest_order_dp((DP)p,weight,n);
! 2954: STOQ(h,q);
! 2955: BDY(r) = (pointer)q;
! 2956: }
! 2957: if ( r0 ) NEXT(r) = 0;
! 2958: MKLIST(l,r0);
! 2959: return l;
! 2960: }
! 2961:
! 2962: LIST dp_order(LIST f,struct order_spec *ord)
! 2963: {
! 2964: int n,l,i;
! 2965: struct weight_or_block *worb;
! 2966: int *weight;
! 2967:
! 2968: switch ( ord->id ) {
! 2969: case 2: /* matrix order */
! 2970: /* extract the first row */
! 2971: n = ord->nv;
! 2972: weight = ord->ord.matrix.matrix[0];
! 2973: return highest_order(f,weight,n);
! 2974: case 3: /* composite order */
! 2975: /* the first w_or_b */
! 2976: worb = ord->ord.composite.w_or_b;
! 2977: switch ( worb->type ) {
! 2978: case IS_DENSE_WEIGHT:
! 2979: n = worb->length;
! 2980: weight = worb->body.dense_weight;
! 2981: return highest_order(f,weight,n);
! 2982: case IS_SPARSE_WEIGHT:
! 2983: n = ord->nv;
! 2984: weight = (int *)ALLOCA(n*sizeof(int));
! 2985: for ( i = 0; i < n; i++ ) weight[i] = 0;
! 2986: l = worb->length;
! 2987: for ( i = 0; i < l; i++ )
! 2988: weight[worb->body.sparse_weight[i].pos]
! 2989: = worb->body.sparse_weight[i].value;
! 2990: return highest_order(f,weight,n);
! 2991: default:
! 2992: error("dp_initial_term : unsupported order");
! 2993: }
! 2994: default:
! 2995: error("dp_initial_term : unsupported order");
! 2996: }
! 2997: }
! 2998:
! 2999: int dpv_ht(DPV p,DP *h)
! 3000: {
! 3001: int len,max,maxi,i,t;
! 3002: DP *e;
! 3003: MP m,mr;
! 3004:
! 3005: len = p->len;
! 3006: e = p->body;
! 3007: max = -1;
! 3008: maxi = -1;
! 3009: for ( i = 0; i < len; i++ )
! 3010: if ( e[i] && (t = BDY(e[i])->dl->td) > max ) {
! 3011: max = t;
! 3012: maxi = i;
! 3013: }
! 3014: if ( max < 0 ) {
! 3015: *h = 0;
! 3016: return -1;
! 3017: } else {
! 3018: m = BDY(e[maxi]);
! 3019: NEWMP(mr); mr->dl = m->dl; mr->c = (Obj)ONE; NEXT(mr) = 0;
! 3020: MKDP(e[maxi]->nv,mr,*h); (*h)->sugar = mr->dl->td; /* XXX */
! 3021: return maxi;
! 3022: }
! 3023: }
! 3024:
! 3025: /* return 1 if 0 <_w1 v && v <_w2 0 */
! 3026:
! 3027: int in_c12(int n,int *v,int row1,int **w1,int row2, int **w2)
! 3028: {
! 3029: int t1,t2;
! 3030:
! 3031: t1 = compare_zero(n,v,row1,w1);
! 3032: t2 = compare_zero(n,v,row2,w2);
! 3033: if ( t1 > 0 && t2 < 0 ) return 1;
! 3034: else return 0;
! 3035: }
! 3036:
! 3037: /* 0 < u => 1, 0 > u => -1 */
! 3038:
! 3039: int compare_zero(int n,int *u,int row,int **w)
! 3040: {
! 3041: int i,j,t;
! 3042: int *wi;
! 3043:
! 3044: for ( i = 0; i < row; i++ ) {
! 3045: wi = w[i];
! 3046: for ( j = 0, t = 0; j < n; j++ ) t += u[j]*wi[j];
! 3047: if ( t > 0 ) return 1;
! 3048: else if ( t < 0 ) return -1;
! 3049: }
! 3050: return 0;
! 3051: }
! 3052:
! 3053: /* functions for generic groebner walk */
! 3054: /* u=0 means u=-infty */
! 3055:
! 3056: int compare_facet_preorder(int n,int *u,int *v,
! 3057: int row1,int **w1,int row2,int **w2)
! 3058: {
! 3059: int i,j,s,t,tu,tv;
! 3060: int *w2i,*uv;
! 3061:
! 3062: if ( !u ) return 1;
! 3063: uv = W_ALLOC(n);
! 3064: for ( i = 0; i < row2; i++ ) {
! 3065: w2i = w2[i];
! 3066: for ( j = 0, tu = tv = 0; j < n; j++ )
! 3067: if ( s = w2i[j] ) {
! 3068: tu += s*u[j]; tv += s*v[j];
! 3069: }
! 3070: for ( j = 0; j < n; j++ ) uv[j] = u[j]*tv-v[j]*tu;
! 3071: t = compare_zero(n,uv,row1,w1);
! 3072: if ( t > 0 ) return 1;
! 3073: else if ( t < 0 ) return 0;
! 3074: }
! 3075: return 1;
! 3076: }
! 3077:
! 3078: Q inner_product_with_small_vector(VECT w,int *v)
! 3079: {
! 3080: int n,i;
! 3081: Z q;
! 3082: Q s,t,u;
! 3083:
! 3084: n = w->len;
! 3085: s = 0;
! 3086: for ( i = 0; i < n; i++ ) {
! 3087: STOQ(v[i],q); mulq((Q)w->body[i],(Q)q,&t); addq(t,s,&u); s = u;
! 3088: }
! 3089: return s;
! 3090: }
! 3091:
! 3092: Q compute_last_t(NODE g,NODE gh,Q t,VECT w1,VECT w2,NODE *homo,VECT *wp)
! 3093: {
! 3094: int n,i;
! 3095: int *wt;
! 3096: Q last,d1,d2,dn,nm,s,t1;
! 3097: VECT wd,wt1,wt2,w;
! 3098: NODE tg,tgh;
! 3099: MP f;
! 3100: int *h;
! 3101: NODE r0,r;
! 3102: MP m0,m;
! 3103: DP d;
! 3104:
! 3105: n = w1->len;
! 3106: wt = W_ALLOC(n);
! 3107: last = (Q)ONE;
! 3108: /* t1 = 1-t */
! 3109: for ( tg = g, tgh = gh; tg; tg = NEXT(tg), tgh = NEXT(tgh ) ) {
! 3110: f = BDY((DP)BDY(tg));
! 3111: h = BDY((DP)BDY(tgh))->dl->d;
! 3112: for ( ; f; f = NEXT(f) ) {
! 3113: for ( i = 0; i < n; i++ ) wt[i] = h[i]-f->dl->d[i];
! 3114: for ( i = 0; i < n && !wt[i]; i++ );
! 3115: if ( i == n ) continue;
! 3116: d1 = inner_product_with_small_vector(w1,wt);
! 3117: d2 = inner_product_with_small_vector(w2,wt);
! 3118: nm = d1; subq(d1,d2,&dn);
! 3119: /* if d1=d2 then nothing happens */
! 3120: if ( !dn ) continue;
! 3121: /* s satisfies ds = 0*/
! 3122: divq(nm,dn,&s);
! 3123:
! 3124: if ( cmpq(s,t) > 0 && cmpq(s,last) < 0 )
! 3125: last = s;
! 3126: else if ( !cmpq(s,t) ) {
! 3127: if ( cmpq(d2,0) < 0 ) {
! 3128: last = t;
! 3129: break;
! 3130: }
! 3131: }
! 3132: }
! 3133: }
! 3134: nmq(last,(Z *)&nm);
! 3135: dnq(last,(Z *)&dn);
! 3136: /* (1-n/d)*w1+n/d*w2 -> w=(d-n)*w1+n*w2 */
! 3137: subq(dn,nm,&t1); mulvect(CO,(Obj)w1,(Obj)t1,(Obj *)&wt1);
! 3138: mulvect(CO,(Obj)w2,(Obj)nm,(Obj *)&wt2); addvect(CO,wt1,wt2,&w);
! 3139:
! 3140: r0 = 0;
! 3141: for ( tg = g, tgh = gh; tg; tg = NEXT(tg), tgh = NEXT(tgh ) ) {
! 3142: f = BDY((DP)BDY(tg));
! 3143: h = BDY((DP)BDY(tgh))->dl->d;
! 3144: for ( m0 = 0; f; f = NEXT(f) ) {
! 3145: for ( i = 0; i < n; i++ ) wt[i] = h[i]-f->dl->d[i];
! 3146: for ( i = 0; i < n && !wt[i]; i++ );
! 3147: if ( !inner_product_with_small_vector(w,wt) ) {
! 3148: NEXTMP(m0,m); m->c = f->c; m->dl = f->dl;
! 3149: }
! 3150: }
! 3151: NEXT(m) = 0;
! 3152: MKDP(((DP)BDY(tg))->nv,m0,d); d->sugar = ((DP)BDY(tg))->sugar;
! 3153: NEXTNODE(r0,r); BDY(r) = (pointer)d;
! 3154: }
! 3155: NEXT(r) = 0;
! 3156: *homo = r0;
! 3157: *wp = w;
! 3158: return last;
! 3159: }
! 3160:
! 3161: /* return 0 if last_w = infty */
! 3162:
! 3163: NODE compute_last_w(NODE g,NODE gh,int n,int **w,
! 3164: int row1,int **w1,int row2,int **w2)
! 3165: {
! 3166: DP d;
! 3167: MP f,m0,m;
! 3168: int *wt,*v,*h;
! 3169: NODE t,s,n0,tn,n1,r0,r;
! 3170: int i;
! 3171:
! 3172: wt = W_ALLOC(n);
! 3173: n0 = 0;
! 3174: for ( t = g, s = gh; t; t = NEXT(t), s = NEXT(s) ) {
! 3175: f = BDY((DP)BDY(t));
! 3176: h = BDY((DP)BDY(s))->dl->d;
! 3177: for ( ; f; f = NEXT(f) ) {
! 3178: for ( i = 0; i < n; i++ ) wt[i] = h[i]-f->dl->d[i];
! 3179: for ( i = 0; i < n && !wt[i]; i++ );
! 3180: if ( i == n ) continue;
! 3181:
! 3182: if ( in_c12(n,wt,row1,w1,row2,w2) &&
! 3183: compare_facet_preorder(n,*w,wt,row1,w1,row2,w2) ) {
! 3184: v = (int *)MALLOC_ATOMIC(n*sizeof(int));
! 3185: for ( i = 0; i < n; i++ ) v[i] = wt[i];
! 3186: MKNODE(n1,v,n0); n0 = n1;
! 3187: }
! 3188: }
! 3189: }
! 3190: if ( !n0 ) return 0;
! 3191: for ( t = n0; t; t = NEXT(t) ) {
! 3192: v = (int *)BDY(t);
! 3193: for ( s = n0; s; s = NEXT(s) )
! 3194: if ( !compare_facet_preorder(n,v,(int *)BDY(s),row1,w1,row2,w2) )
! 3195: break;
! 3196: if ( !s ) {
! 3197: *w = v;
! 3198: break;
! 3199: }
! 3200: }
! 3201: if ( !t )
! 3202: error("compute_last_w : cannot happen");
! 3203: r0 = 0;
! 3204: for ( t = g, s = gh; t; t = NEXT(t), s = NEXT(s) ) {
! 3205: f = BDY((DP)BDY(t));
! 3206: h = BDY((DP)BDY(s))->dl->d;
! 3207: for ( m0 = 0; f; f = NEXT(f) ) {
! 3208: for ( i = 0; i < n; i++ ) wt[i] = h[i]-f->dl->d[i];
! 3209: for ( i = 0; i < n && !wt[i]; i++ );
! 3210: if ( i == n ||
! 3211: (compare_facet_preorder(n,wt,*w,row1,w1,row2,w2)
! 3212: && compare_facet_preorder(n,*w,wt,row1,w1,row2,w2)) ) {
! 3213: NEXTMP(m0,m); m->c = f->c; m->dl = f->dl;
! 3214: }
! 3215: }
! 3216: NEXT(m) = 0;
! 3217: MKDP(((DP)BDY(t))->nv,m0,d); d->sugar = ((DP)BDY(t))->sugar;
! 3218: NEXTNODE(r0,r); BDY(r) = (pointer)d;
! 3219: }
! 3220: NEXT(r) = 0;
! 3221: return r0;
! 3222: }
! 3223:
! 3224: /* compute a sufficient set of d(f)=u-v */
! 3225:
! 3226: NODE compute_essential_df(DP *g,DP *gh,int ng)
! 3227: {
! 3228: int nv,i,j,k,t,lj;
! 3229: NODE r,r1,ri,rt,r0;
! 3230: MP m;
! 3231: MP *mj;
! 3232: DL di,hj,dl,dlt;
! 3233: int *d,*dt;
! 3234: LIST l;
! 3235: Z q;
! 3236:
! 3237: nv = g[0]->nv;
! 3238: r = 0;
! 3239: for ( j = 0; j < ng; j++ ) {
! 3240: for ( m = BDY(g[j]), lj = 0; m; m = NEXT(m), lj++ );
! 3241: mj = (MP *)ALLOCA(lj*sizeof(MP));
! 3242: for ( m = BDY(g[j]), k = 0; m; m = NEXT(m), k++ )
! 3243: mj[k] = m;
! 3244: for ( i = 0; i < lj; i++ ) {
! 3245: for ( di = mj[i]->dl, k = i+1; k < lj; k++ )
! 3246: if ( _dl_redble(di,mj[k]->dl,nv) ) break;
! 3247: if ( k < lj ) mj[i] = 0;
! 3248: }
! 3249: hj = BDY(gh[j])->dl;
! 3250: _NEWDL(dl,nv); d = dl->d;
! 3251: r0 = r;
! 3252: for ( i = 0; i < lj; i++ ) {
! 3253: if ( mj[i] && !dl_equal(nv,di=mj[i]->dl,hj) ) {
! 3254: for ( k = 0, t = 0; k < nv; k++ ) {
! 3255: d[k] = hj->d[k]-di->d[k];
! 3256: t += d[k];
! 3257: }
! 3258: dl->td = t;
! 3259: #if 1
! 3260: for ( rt = r0; rt; rt = NEXT(rt) ) {
! 3261: dlt = (DL)BDY(rt);
! 3262: if ( dlt->td != dl->td ) continue;
! 3263: for ( dt = dlt->d, k = 0; k < nv; k++ )
! 3264: if ( d[k] != dt[k] ) break;
! 3265: if ( k == nv ) break;
! 3266: }
! 3267: #else
! 3268: rt = 0;
! 3269: #endif
! 3270: if ( !rt ) {
! 3271: MKNODE(r1,dl,r); r = r1;
! 3272: _NEWDL(dl,nv); d = dl->d;
! 3273: }
! 3274: }
! 3275: }
! 3276: }
! 3277: for ( rt = r; rt; rt = NEXT(rt) ) {
! 3278: dl = (DL)BDY(rt); d = dl->d;
! 3279: ri = 0;
! 3280: for ( k = nv-1; k >= 0; k-- ) {
! 3281: STOQ(d[k],q);
! 3282: MKNODE(r1,q,ri); ri = r1;
! 3283: }
! 3284: MKNODE(r1,0,ri); MKLIST(l,r1);
! 3285: BDY(rt) = (pointer)l;
! 3286: }
! 3287: return r;
! 3288: }
! 3289:
! 3290: int comp_bits_divisible(int *a,int *b,int n)
! 3291: {
! 3292: int bpi,i,wi,bi;
! 3293:
! 3294: bpi = (sizeof(int)/sizeof(char))*8;
! 3295: for ( i = 0; i < n; i++ ) {
! 3296: wi = i/bpi; bi = i%bpi;
! 3297: if ( !(a[wi]&(1<<bi)) && (b[wi]&(1<<bi)) ) return 0;
! 3298: }
! 3299: return 1;
! 3300: }
! 3301:
! 3302: int comp_bits_lex(int *a,int *b,int n)
! 3303: {
! 3304: int bpi,i,wi,ba,bb,bi;
! 3305:
! 3306: bpi = (sizeof(int)/sizeof(char))*8;
! 3307: for ( i = 0; i < n; i++ ) {
! 3308: wi = i/bpi; bi = i%bpi;
! 3309: ba = (a[wi]&(1<<bi))?1:0;
! 3310: bb = (b[wi]&(1<<bi))?1:0;
! 3311: if ( ba > bb ) return 1;
! 3312: else if ( ba < bb ) return -1;
! 3313: }
! 3314: return 0;
! 3315: }
! 3316:
! 3317: NODE mono_raddec(NODE ideal)
! 3318: {
! 3319: DP p;
! 3320: int nv,w,i,bpi,di,c,len;
! 3321: int *d,*s,*u,*new;
! 3322: NODE t,t1,v,r,rem,prev;
! 3323:
! 3324: if( !ideal ) return 0;
! 3325: p = (DP)BDY(ideal);
! 3326: nv = NV(p);
! 3327: bpi = (sizeof(int)/sizeof(char))*8;
! 3328: w = (nv+(bpi-1))/bpi;
! 3329: d = p->body->dl->d;
! 3330: if ( !NEXT(ideal) ) {
! 3331: for ( t = 0, i = nv-1; i >= 0; i-- ) {
! 3332: if ( d[i] ) {
! 3333: s = (int *)CALLOC(w,sizeof(int));
! 3334: s[i/bpi] |= 1<<(i%bpi);
! 3335: MKNODE(t1,s,t);
! 3336: t = t1;
! 3337: }
! 3338: }
! 3339: return t;
! 3340: }
! 3341: rem = mono_raddec(NEXT(ideal));
! 3342: r = 0;
! 3343: len = w*sizeof(int);
! 3344: u = (int *)CALLOC(w,sizeof(int));
! 3345: for ( i = nv-1; i >= 0; i-- ) {
! 3346: if ( d[i] ) {
! 3347: for ( t = rem; t; t = NEXT(t) ) {
! 3348: bcopy((char *)BDY(t),(char *)u,len);
! 3349: u[i/bpi] |= 1<<(i%bpi);
! 3350: for ( v = r; v; v = NEXT(v) ) {
! 3351: if ( comp_bits_divisible(u,(int *)BDY(v),nv) ) break;
! 3352: }
! 3353: if ( v ) continue;
! 3354: for ( v = r, prev = 0; v; v = NEXT(v) ) {
! 3355: if ( comp_bits_divisible((int *)BDY(v),u,nv) ) {
! 3356: if ( prev ) NEXT(prev) = NEXT(v);
! 3357: else r = NEXT(r);
! 3358: } else prev =v;
! 3359: }
! 3360: for ( v = r, prev = 0; v; prev = v, v = NEXT(v) ) {
! 3361: if ( comp_bits_lex(u,(int *)BDY(v),nv) < 0 ) break;
! 3362: }
! 3363: new = (int *)CALLOC(w,sizeof(int));
! 3364: bcopy((char *)u,(char *)new,len);
! 3365: MKNODE(t1,new,v);
! 3366: if ( prev ) NEXT(prev) = t1;
! 3367: else r = t1;
! 3368: }
! 3369: }
! 3370: }
! 3371: return r;
! 3372: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to dp-supp.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2018 / builtin