=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/dp-supp.c,v retrieving revision 1.49 retrieving revision 1.66 diff -u -p -r1.49 -r1.66 --- OpenXM_contrib2/asir2000/builtin/dp-supp.c 2009/01/04 05:44:51 1.49 +++ OpenXM_contrib2/asir2000/builtin/dp-supp.c 2017/08/31 02:36:20 1.66 @@ -45,7 +45,7 @@ * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. * - * $OpenXM: OpenXM_contrib2/asir2000/builtin/dp-supp.c,v 1.48 2007/10/21 07:47:59 noro Exp $ + * $OpenXM: OpenXM_contrib2/asir2000/builtin/dp-supp.c,v 1.65 2017/03/27 09:05:46 noro Exp $ */ #include "ca.h" #include "base.h" @@ -53,7 +53,7 @@ #include "parse.h" #include "ox.h" -#define HMAG(p) (p_mag(BDY(p)->c)) +#define HMAG(p) (p_mag((P)BDY(p)->c)) extern int (*cmpdl)(); extern double pz_t_e,pz_t_d,pz_t_d1,pz_t_c; @@ -65,6 +65,7 @@ extern NODE TraceList; int show_orderspec; void print_composite_order_spec(struct order_spec *spec); +void dpm_rest(DPM,DPM *); /* * content reduction @@ -112,11 +113,11 @@ void dp_ptozp(DP p,DP *rp) if ( NUM(m->c) ) w[i] = (Q)m->c; else - ptozp(m->c,1,&w[i],&t); + ptozp((P)m->c,1,&w[i],&t); sortbynm(w,n); qltozl(w,n,&dvr); for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - NEXTMP(mr0,mr); divsp(CO,m->c,(P)dvr,&mr->c); mr->dl = m->dl; + NEXTMP(mr0,mr); divsp(CO,(P)m->c,(P)dvr,(P *)&mr->c); mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; } @@ -146,6 +147,58 @@ void dp_ptozp2(DP p0,DP p1,DP *hp,DP *rp) *hp = h; *rp = r; } +void dpm_ptozp(DPM p,DPM *rp) +{ + DMM m,mr,mr0; + int i,n; + Q *w; + Q dvr; + P t; + + if ( !p ) + *rp = 0; + else { + for ( m =BDY(p), n = 0; m; m = NEXT(m), n++ ); + w = (Q *)ALLOCA(n*sizeof(Q)); + for ( m =BDY(p), i = 0; i < n; m = NEXT(m), i++ ) + if ( NUM(m->c) ) + w[i] = (Q)m->c; + else + ptozp((P)m->c,1,&w[i],&t); + sortbynm(w,n); + qltozl(w,n,&dvr); + for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { + NEXTDMM(mr0,mr); divsp(CO,(P)m->c,(P)dvr,(P *)&mr->c); mr->dl = m->dl; mr->pos = m->pos; + } + NEXT(mr) = 0; MKDPM(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; + } +} + +void dpm_ptozp2(DPM p0,DPM p1,DPM *hp,DPM *rp) +{ + DPM t,s,h,r; + DMM m,mr,mr0,m0; + + adddpm(CO,p0,p1,&t); dpm_ptozp(t,&s); + if ( !p0 ) { + h = 0; r = s; + } else if ( !p1 ) { + h = s; r = 0; + } else { + for ( mr0 = 0, m = BDY(s), m0 = BDY(p0); m0; + m = NEXT(m), m0 = NEXT(m0) ) { + NEXTDMM(mr0,mr); mr->c = m->c; mr->dl = m->dl; mr->pos = m->pos; + } + NEXT(mr) = 0; MKDPM(p0->nv,mr0,h); MKDPM(p0->nv,m,r); + } + if ( h ) + h->sugar = p0->sugar; + if ( r ) + r->sugar = p1->sugar; + *hp = h; *rp = r; +} + + void dp_ptozp3(DP p,Q *dvr,DP *rp) { MP m,mr,mr0; @@ -162,11 +215,11 @@ void dp_ptozp3(DP p,Q *dvr,DP *rp) if ( NUM(m->c) ) w[i] = (Q)m->c; else - ptozp(m->c,1,&w[i],&t); + ptozp((P)m->c,1,&w[i],&t); sortbynm(w,n); qltozl(w,n,dvr); for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - NEXTMP(mr0,mr); divsp(CO,m->c,(P)(*dvr),&mr->c); mr->dl = m->dl; + NEXTMP(mr0,mr); divsp(CO,(P)m->c,(P)(*dvr),(P *)&mr->c); mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; } @@ -193,7 +246,7 @@ void dp_idiv(DP p,Q c,DP *rp) divsn(NM((Q)(m->c)),nm,&q); s = sgn*SGN((Q)(m->c)); NTOQ(q,s,t); - mr->c = (P)t; + mr->c = (Obj)t; mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); @@ -206,13 +259,29 @@ void dp_mbase(NODE hlist,NODE *mbase) { DL *dl; DL d; - int i,j,n,nvar,td; + int *t; + int i,j,k,n,nvar,td; n = length(hlist); nvar = ((DP)BDY(hlist))->nv; dl = (DL *)MALLOC(n*sizeof(DL)); - for ( i = 0; i < n; i++, hlist = NEXT(hlist) ) - dl[i] = BDY((DP)BDY(hlist))->dl; NEWDL(d,nvar); *mbase = 0; + for ( i = 0; i < n; i++, hlist = NEXT(hlist) ) { + dl[i] = BDY((DP)BDY(hlist))->dl; + /* trivial ideal check */ + if ( (*cmpdl)(nvar,d,dl[i]) == 0 ) { + return; + } + } + /* zero-dim. ideal check */ + for ( i = 0; i < nvar; i++ ) { + for ( j = 0; j < n; j++ ) { + for ( k = 0, t = dl[j]->d; k < nvar; k++ ) + if ( k != i && t[k] != 0 ) break; + if ( k == nvar ) break; + } + if ( j == n ) + error("dp_mbase : input ideal is not zero-dimensional"); + } while ( 1 ) { insert_to_node(d,mbase,nvar); for ( i = nvar-1; i >= 0; ) { @@ -261,7 +330,7 @@ void insert_to_node(DL d,NODE *n,int nvar) NEWDL(d1,nvar); d1->td = d->td; bcopy((char *)d->d,(char *)d1->d,nvar*sizeof(int)); - NEWMP(m); m->dl = d1; m->c = (P)ONE; NEXT(m) = 0; + NEWMP(m); m->dl = d1; m->c = (Obj)ONE; NEXT(m) = 0; MKDP(nvar,m,dp); dp->sugar = d->td; if ( !(*n) ) { MKNODE(n1,dp,0); *n = n1; @@ -290,7 +359,7 @@ void dp_vtod(Q *c,DP p,DP *rp) *rp = 0; else { for ( mr0 = 0, m = BDY(p), i = 0; m; m = NEXT(m), i++ ) { - NEXTMP(mr0,mr); mr->c = (P)c[i]; mr->dl = m->dl; + NEXTMP(mr0,mr); mr->c = (Obj)c[i]; mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; @@ -332,7 +401,7 @@ void dp_ptozp_d(DP p,DP *rp) if ( PCoeffs ) { dp_ptozp(p,rp); return; } - if ( !Dist || p_mag(BDY(p)->c) <= mpi_mag ) { + if ( !Dist || p_mag((P)BDY(p)->c) <= mpi_mag ) { dist = 0; ndist = 0; if ( DP_NFStat ) fprintf(asir_out,"L"); } else { @@ -458,7 +527,7 @@ void dp_monic_sf(DP p,DP *rp) if ( !p ) *rp = 0; else { - head_coef(BDY(p)->c,&c); + head_coef((P)BDY(p)->c,&c); divsdc(CO,p,(P)c,rp); } } @@ -487,13 +556,13 @@ void dp_prim(DP p,DP *rp) for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ ); w = (P *)ALLOCA(n*sizeof(P)); for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ ) - w[i] = m->c; + w[i] = (P)m->c; gcdsf(CO,w,n,&g); if ( NUM(g) ) dp_monic_sf(p,rp); else { for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - NEXTMP(mr0,mr); divsp(CO,m->c,g,&mr->c); mr->dl = m->dl; + NEXTMP(mr0,mr); divsp(CO,(P)m->c,g,(P *)&mr->c); mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,p1); p1->sugar = p->sugar; dp_monic_sf(p1,rp); @@ -508,7 +577,7 @@ void dp_prim(DP p,DP *rp) for ( m = BDY(p), n = 0; m; m = NEXT(m), n++ ); if ( n == 1 ) { m = BDY(p); - NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0; + NEWMP(mr); mr->dl = m->dl; mr->c = (Obj)ONE; NEXT(mr) = 0; MKDP(p->nv,mr,*rp); (*rp)->sugar = p->sugar; return; } @@ -518,13 +587,13 @@ void dp_prim(DP p,DP *rp) if ( NUM(m->c) ) { c[i] = (Q)m->c; w[i] = (P)ONE; } else - ptozp(m->c,1,&c[i],&w[i]); + ptozp((P)m->c,1,&c[i],&w[i]); qltozl(c,n,&dvr); heu_nezgcdnpz(CO,w,n,&t); mulp(CO,t,(P)dvr,&g); if ( NUM(g) ) *rp = p; else { for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - NEXTMP(mr0,mr); divsp(CO,m->c,g,&mr->c); mr->dl = m->dl; + NEXTMP(mr0,mr); divsp(CO,(P)m->c,g,(P *)&mr->c); mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; add_denomlist(g); @@ -572,11 +641,11 @@ void dp_prim_mod(DP p,int mod,DP *rp) else if ( NoGCD ) *rp = p; else { - for ( m = BDY(p), g = m->c, m = NEXT(m); m; m = NEXT(m) ) { - gcdprsmp(CO,mod,g,m->c,&t); g = t; + for ( m = BDY(p), g = (P)m->c, m = NEXT(m); m; m = NEXT(m) ) { + gcdprsmp(CO,mod,g,(P)m->c,&t); g = t; } for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - NEXTMP(mr0,mr); divsmp(CO,mod,m->c,g,&mr->c); mr->dl = m->dl; + NEXTMP(mr0,mr); divsmp(CO,mod,(P)m->c,g,(P *)&mr->c); mr->dl = m->dl; } NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; } @@ -637,13 +706,13 @@ void dp_sp(DP p1,DP p2,DP *rp) } } - NEWMP(m); m->dl = d; m->c = (P)c2; NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)c2; NEXT(m) = 0; MKDP(n,m,s1); s1->sugar = d->td; muld(CO,s1,p1,&t); NEWDL(d,n); d->td = td - d2->td; for ( i = 0; i < n; i++ ) d->d[i] = w[i] - d2->d[i]; - NEWMP(m); m->dl = d; m->c = (P)c1; NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0; MKDP(n,m,s2); s2->sugar = d->td; muld(CO,s2,p2,&u); subd(CO,t,u,rp); @@ -651,17 +720,75 @@ void dp_sp(DP p1,DP p2,DP *rp) LIST hist; NODE node; - node = mknode(4,ONE,0,s1,ONE); + node = mknode(4,ONE,NULLP,s1,ONE); MKLIST(hist,node); MKNODE(TraceList,hist,0); - node = mknode(4,ONE,0,0,ONE); + node = mknode(4,ONE,NULLP,NULLP,ONE); chsgnd(s2,(DP *)&ARG2(node)); MKLIST(hist,node); MKNODE(node,hist,TraceList); TraceList = node; } } +void dpm_sp(DPM p1,DPM p2,DPM *rp) +{ + int i,n,td; + int *w; + DL d1,d2,d; + MP m; + DP s1,s2; + DPM t,u; + Q c,c1,c2; + N gn,tn; + + n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl; + if ( BDY(p1)->pos != BDY(p2)->pos ) { + *rp = 0; + return; + } + w = (int *)ALLOCA(n*sizeof(int)); + for ( i = 0, td = 0; i < n; i++ ) { + w[i] = MAX(d1->d[i],d2->d[i]); td += MUL_WEIGHT(w[i],i); + } + + NEWDL(d,n); d->td = td - d1->td; + for ( i = 0; i < n; i++ ) + d->d[i] = w[i] - d1->d[i]; + c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c; + if ( INT(c1) && INT(c2) ) { + gcdn(NM(c1),NM(c2),&gn); + if ( !UNIN(gn) ) { + divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c; + divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c; + } + } + + NEWMP(m); m->dl = d; m->c = (Obj)c2; NEXT(m) = 0; + MKDP(n,m,s1); s1->sugar = d->td; mulobjdpm(CO,(Obj)s1,p1,&t); + + NEWDL(d,n); d->td = td - d2->td; + for ( i = 0; i < n; i++ ) + d->d[i] = w[i] - d2->d[i]; + NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0; + MKDP(n,m,s2); s2->sugar = d->td; mulobjdpm(CO,(Obj)s2,p2,&u); + + subdpm(CO,t,u,rp); + if ( GenTrace ) { + LIST hist; + NODE node; + + node = mknode(4,ONE,NULLP,s1,ONE); + MKLIST(hist,node); + MKNODE(TraceList,hist,0); + + node = mknode(4,ONE,NULLP,NULLP,ONE); + chsgnd(s2,(DP *)&ARG2(node)); + MKLIST(hist,node); + MKNODE(node,hist,TraceList); TraceList = node; + } +} + void _dp_sp_dup(DP p1,DP p2,DP *rp) { int i,n,td; @@ -690,13 +817,13 @@ void _dp_sp_dup(DP p1,DP p2,DP *rp) } } - _NEWMP(m); m->dl = d; m->c = (P)c2; NEXT(m) = 0; + _NEWMP(m); m->dl = d; m->c = (Obj)c2; NEXT(m) = 0; _MKDP(n,m,s1); s1->sugar = d->td; _muld_dup(CO,s1,p1,&t); _free_dp(s1); _NEWDL(d,n); d->td = td - d2->td; for ( i = 0; i < n; i++ ) d->d[i] = w[i] - d2->d[i]; - _NEWMP(m); m->dl = d; chsgnp((P)c1,&m->c); NEXT(m) = 0; + _NEWMP(m); m->dl = d; chsgnp((P)c1,(P *)&m->c); NEXT(m) = 0; _MKDP(n,m,s2); s2->sugar = d->td; _muld_dup(CO,s2,p2,&u); _free_dp(s2); _addd_destructive(CO,t,u,rp); @@ -704,11 +831,11 @@ void _dp_sp_dup(DP p1,DP p2,DP *rp) LIST hist; NODE node; - node = mknode(4,ONE,0,s1,ONE); + node = mknode(4,ONE,NULLP,s1,ONE); MKLIST(hist,node); MKNODE(TraceList,hist,0); - node = mknode(4,ONE,0,0,ONE); + node = mknode(4,ONE,NULLP,NULLP,ONE); chsgnd(s2,(DP *)&ARG2(node)); MKLIST(hist,node); MKNODE(node,hist,TraceList); TraceList = node; @@ -731,12 +858,12 @@ void dp_sp_mod(DP p1,DP p2,int mod,DP *rp) NEWDL_NOINIT(d,n); d->td = td - d1->td; for ( i = 0; i < n; i++ ) d->d[i] = w[i] - d1->d[i]; - NEWMP(m); m->dl = d; m->c = (P)BDY(p2)->c; NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)BDY(p2)->c; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p1,s,&t); NEWDL_NOINIT(d,n); d->td = td - d2->td; for ( i = 0; i < n; i++ ) d->d[i] = w[i] - d2->d[i]; - NEWMP(m); m->dl = d; m->c = (P)BDY(p1)->c; NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)BDY(p1)->c; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,p2,s,&u); submd(CO,mod,t,u,rp); } @@ -762,7 +889,7 @@ void _dp_sp_mod_dup(DP p1,DP p2,int mod,DP *rp) _NEWDL(d,n); d->td = td - d2->td; for ( i = 0; i < n; i++ ) d->d[i] = w[i] - d2->d[i]; - _NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0; + _NEWMP(m); m->dl = d; m->c = (Obj)STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0; _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p2,&u); _free_dp(s); _addmd_destructive(mod,t,u,rp); } @@ -788,7 +915,7 @@ void _dp_sp_mod(DP p1,DP p2,int mod,DP *rp) NEWDL(d,n); d->td = td - d2->td; for ( i = 0; i < n; i++ ) d->d[i] = w[i] - d2->d[i]; - NEWMP(m); m->dl = d; m->c = STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)STOI(mod - ITOS(BDY(p1)->c)); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; mulmd_dup(mod,s,p2,&u); addmd_destructive(mod,t,u,rp); } @@ -833,13 +960,57 @@ void dp_red(DP p0,DP p1,DP p2,DP *head,DP *rest,P *dnp divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a; add_denomlist(g); } - NEWMP(m); m->dl = d; chsgnp((P)c1,&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; + NEWMP(m); m->dl = d; chsgnp((P)c1,(P *)&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; *multp = s; - muld(CO,s,p2,&t); muldc(CO,p1,(P)c2,&s); addd(CO,s,t,&r); - muldc(CO,p0,(P)c2,&h); + muld(CO,s,p2,&t); muldc(CO,p1,(Obj)c2,&s); addd(CO,s,t,&r); + muldc(CO,p0,(Obj)c2,&h); *head = h; *rest = r; *dnp = (P)c2; } +void dpm_red(DPM p0,DPM p1,DPM p2,DPM *head,DPM *rest,P *dnp,DP *multp) +{ + int i,n,pos; + DL d1,d2,d; + MP m; + DP s; + DPM t,r,h,u,w; + Q c,c1,c2; + N gn,tn; + P g,a; + P p[2]; + + n = p1->nv; d1 = BDY(p1)->dl; d2 = BDY(p2)->dl; pos = BDY(p1)->pos; + if ( pos != BDY(p2)->pos ) + error("dpm_red : cannot happen"); + NEWDL(d,n); d->td = d1->td - d2->td; + for ( i = 0; i < n; i++ ) + d->d[i] = d1->d[i]-d2->d[i]; + c1 = (Q)BDY(p1)->c; c2 = (Q)BDY(p2)->c; + if ( dp_fcoeffs == N_GFS ) { + p[0] = (P)c1; p[1] = (P)c2; + gcdsf(CO,p,2,&g); + divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a; + } else if ( dp_fcoeffs ) { + /* do nothing */ + } else if ( INT(c1) && INT(c2) ) { + gcdn(NM(c1),NM(c2),&gn); + if ( !UNIN(gn) ) { + divsn(NM(c1),gn,&tn); NTOQ(tn,SGN(c1),c); c1 = c; + divsn(NM(c2),gn,&tn); NTOQ(tn,SGN(c2),c); c2 = c; + } + } else { + ezgcdpz(CO,(P)c1,(P)c2,&g); + divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a; + add_denomlist(g); + } + NEWMP(m); m->dl = d; chsgnp((P)c1,(P *)&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; + *multp = s; + mulobjdpm(CO,(Obj)s,p2,&u); mulobjdpm(CO,(Obj)c2,p1,&w); adddpm(CO,u,w,&r); + mulobjdpm(CO,(Obj)c2,p0,&h); + *head = h; *rest = r; *dnp = (P)c2; +} + + /* * m-reduction by a marked poly * do content reduction over Z or Q(x,...) @@ -880,14 +1051,14 @@ void dp_red_marked(DP p0,DP p1,DP p2,DP hp2,DP *head,D ezgcdpz(CO,(P)c1,(P)c2,&g); divsp(CO,(P)c1,g,&a); c1 = (Q)a; divsp(CO,(P)c2,g,&a); c2 = (Q)a; } - NEWMP(m); m->dl = d; m->c = (P)c1; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; + NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; *multp = s; - muld(CO,s,p2,&t); muldc(CO,p1,(P)c2,&s); subd(CO,s,t,&r); - muldc(CO,p0,(P)c2,&h); + muld(CO,s,p2,&t); muldc(CO,p1,(Obj)c2,&s); subd(CO,s,t,&r); + muldc(CO,p0,(Obj)c2,&h); *head = h; *rest = r; *dnp = (P)c2; } -void dp_red_marked_mod(DP p0,DP p1,DP p2,DP hp2,int mod,DP *head,DP *rest,P *dnp) +void dp_red_marked_mod(DP p0,DP p1,DP p2,DP hp2,int mod,DP *head,DP *rest,P *dnp,DP *multp) { int i,n; DL d1,d2,d; @@ -905,12 +1076,14 @@ void dp_red_marked_mod(DP p0,DP p1,DP p2,DP hp2,int mo if ( NUM(c2) ) { divsmp(CO,mod,c1,c2,&u); c1 = u; c2 = (P)ONEM; } - NEWMP(m); m->dl = d; chsgnmp(mod,(P)c1,&m->c); NEXT(m) = 0; - MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,s,p2,&t); + NEWMP(m); m->dl = d; m->c = (Obj)c1; NEXT(m) = 0; + MKDP(n,m,s); s->sugar = d->td; + *multp = s; + mulmd(CO,mod,s,p2,&t); if ( NUM(c2) ) { - addmd(CO,mod,p1,t,&r); h = p0; + submd(CO,mod,p1,t,&r); h = p0; } else { - mulmdc(CO,mod,p1,c2,&s); addmd(CO,mod,s,t,&r); mulmdc(CO,mod,p0,c2,&h); + mulmdc(CO,mod,p1,c2,&s); submd(CO,mod,s,t,&r); mulmdc(CO,mod,p0,c2,&h); } *head = h; *rest = r; *dnp = c2; } @@ -934,12 +1107,37 @@ void dp_red_f(DP p1,DP p2,DP *rest) NEWMP(m); m->dl = d; divr(CO,(Obj)BDY(p1)->c,(Obj)BDY(p2)->c,&a); chsgnr(a,&b); - C(m) = (P)b; + C(m) = (Obj)b; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; muld(CO,s,p2,&t); addd(CO,p1,t,rest); } +void dpm_red_f(DPM p1,DPM p2,DPM *rest) +{ + int i,n; + DL d1,d2,d; + MP m; + DPM t; + DP s; + Obj a,b; + + n = p1->nv; + d1 = BDY(p1)->dl; d2 = BDY(p2)->dl; + + NEWDL(d,n); d->td = d1->td - d2->td; + for ( i = 0; i < n; i++ ) + d->d[i] = d1->d[i]-d2->d[i]; + + NEWMP(m); m->dl = d; + arf_div(CO,(Obj)BDY(p1)->c,(Obj)BDY(p2)->c,&a); arf_chsgn(a,&b); + C(m) = b; + NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; + + mulobjdpm(CO,(Obj)s,p2,&t); adddpm(CO,p1,t,rest); +} + + void dp_red_mod(DP p0,DP p1,DP p2,int mod,DP *head,DP *rest,P *dnp) { int i,n; @@ -958,7 +1156,7 @@ void dp_red_mod(DP p0,DP p1,DP p2,int mod,DP *head,DP if ( NUM(c2) ) { divsmp(CO,mod,c1,c2,&u); c1 = u; c2 = (P)ONEM; } - NEWMP(m); m->dl = d; chsgnmp(mod,(P)c1,&m->c); NEXT(m) = 0; + NEWMP(m); m->dl = d; chsgnmp(mod,(P)c1,(P *)&m->c); NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; mulmd(CO,mod,s,p2,&t); if ( NUM(c2) ) { addmd(CO,mod,p1,t,&r); h = p0; @@ -986,7 +1184,7 @@ void _dp_red_mod_destructive(DP p1,DP p2,int mod,DP *r c = invm(ITOS(BDY(p2)->c),mod); c2 = ITOS(BDY(p1)->c); DMAR(c,c2,0,mod,c1); - _NEWMP(m); m->dl = d; m->c = STOI(mod-c1); NEXT(m) = 0; + _NEWMP(m); m->dl = d; m->c = (Obj)STOI(mod-c1); NEXT(m) = 0; #if 0 _MKDP(n,m,s); s->sugar = d->td; _mulmd_dup(mod,s,p2,&t); _free_dp(s); @@ -1097,7 +1295,7 @@ void dp_removecont2(DP p1,DP p2,DP *r1p,DP *r2p,Q *con i = 0; if ( p1 ) { for ( m0 = 0, t = BDY(p1); i < n1; i++, t = NEXT(t) ) { - NEXTMP(m0,m); m->c = (P)w[i]; m->dl = t->dl; + NEXTMP(m0,m); m->c = (Obj)w[i]; m->dl = t->dl; } NEXT(m) = 0; MKDP(p1->nv,m0,*r1p); (*r1p)->sugar = p1->sugar; @@ -1105,7 +1303,7 @@ void dp_removecont2(DP p1,DP p2,DP *r1p,DP *r2p,Q *con *r1p = 0; if ( p2 ) { for ( m0 = 0, t = BDY(p2); i < n; i++, t = NEXT(t) ) { - NEXTMP(m0,m); m->c = (P)w[i]; m->dl = t->dl; + NEXTMP(m0,m); m->c = (Obj)w[i]; m->dl = t->dl; } NEXT(m) = 0; MKDP(p2->nv,m0,*r2p); (*r2p)->sugar = p2->sugar; @@ -1182,7 +1380,7 @@ last: void dp_true_nf_marked_mod(NODE b,DP g,DP *ps,DP *hps,int mod,DP *rp,P *dnp) { - DP hp,u,p,d,s,t; + DP hp,u,p,d,s,t,dmy; NODE l; MP m,mr; int i,n; @@ -1203,7 +1401,7 @@ void dp_true_nf_marked_mod(NODE b,DP g,DP *ps,DP *hps, for ( u = 0, i = 0; i < n; i++ ) { if ( dp_redble(g,hp = hps[wb[i]]) ) { p = ps[wb[i]]; - dp_red_marked_mod(d,g,p,hp,mod,&t,&u,&tdn); + dp_red_marked_mod(d,g,p,hp,mod,&t,&u,&tdn,&dmy); psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar; sugar = MAX(sugar,psugar); if ( !u ) { @@ -1247,7 +1445,7 @@ DP *dp_true_nf_and_quotient_marked (NODE b,DP g,DP *ps dn = (P)ONE; if ( !g ) { - *rp = 0; *dnp = dn; return; + *rp = 0; *dnp = dn; return 0; } for ( n = 0, l = b; l; l = NEXT(l), n++ ); wb = (int *)ALLOCA(n*sizeof(int)); @@ -1264,7 +1462,7 @@ DP *dp_true_nf_and_quotient_marked (NODE b,DP g,DP *ps psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar; sugar = MAX(sugar,psugar); for ( j = 0; j < n; j++ ) { - muldc(CO,q[j],(P)tdn,&dmy); q[j] = dmy; + muldc(CO,q[j],(Obj)tdn,&dmy); q[j] = dmy; } addd(CO,q[wb[i]],mult,&dmy); q[wb[i]] = dmy; mulp(CO,dn,tdn,&tdn1); dn = tdn1; @@ -1288,6 +1486,61 @@ last: return q; } +DP *dp_true_nf_and_quotient_marked_mod(NODE b,DP g,DP *ps,DP *hps,int mod,DP *rp,P *dnp) +{ + DP u,p,d,s,t,dmy,hp,mult; + DP *q; + NODE l; + MP m,mr; + int i,n,j; + int *wb; + int sugar,psugar; + P dn,tdn,tdn1; + + for ( n = 0, l = b; l; l = NEXT(l), n++ ); + q = (DP *)MALLOC(n*sizeof(DP)); + for ( i = 0; i < n; i++ ) q[i] = 0; + dn = (P)ONEM; + if ( !g ) { + *rp = 0; *dnp = dn; return 0; + } + wb = (int *)ALLOCA(n*sizeof(int)); + for ( i = 0, l = b; i < n; l = NEXT(l), i++ ) + wb[i] = QTOS((Q)BDY(l)); + sugar = g->sugar; + for ( d = 0; g; ) { + for ( u = 0, i = 0; i < n; i++ ) { + if ( dp_redble(g,hp = hps[wb[i]]) ) { + p = ps[wb[i]]; + dp_red_marked_mod(d,g,p,hp,mod,&t,&u,&tdn,&mult); + psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar; + sugar = MAX(sugar,psugar); + for ( j = 0; j < n; j++ ) { + mulmdc(CO,mod,q[j],(P)tdn,&dmy); q[j] = dmy; + } + addmd(CO,mod,q[wb[i]],mult,&dmy); q[wb[i]] = dmy; + mulmp(CO,mod,dn,tdn,&tdn1); dn = tdn1; + d = t; + if ( !u ) goto last; + break; + } + } + if ( u ) + g = u; + else { + m = BDY(g); NEWMP(mr); mr->dl = m->dl; mr->c = m->c; + NEXT(mr) = 0; MKDP(g->nv,mr,t); t->sugar = mr->dl->td; + addmd(CO,mod,d,t,&s); d = s; + dp_rest(g,&t); g = t; + } + } +last: + if ( d ) + d->sugar = sugar; + *rp = d; *dnp = dn; + return q; +} + /* nf computation over Z */ void dp_nf_z(NODE b,DP g,DP *ps,int full,int multiple,DP *rp) @@ -1359,6 +1612,75 @@ void dp_nf_z(NODE b,DP g,DP *ps,int full,int multiple, *rp = d; } +void dpm_nf_z(NODE b,DPM g,DPM *ps,int full,int multiple,DPM *rp) +{ + DPM u,p,d,s,t; + DP dmy1; + P dmy; + NODE l; + DMM m,mr; + int i,n; + int *wb; + int hmag; + int sugar,psugar; + + if ( !g ) { + *rp = 0; return; + } + for ( n = 0, l = b; l; l = NEXT(l), n++ ); + wb = (int *)ALLOCA(n*sizeof(int)); + for ( i = 0, l = b; i < n; l = NEXT(l), i++ ) + wb[i] = QTOS((Q)BDY(l)); + + hmag = multiple*HMAG(g); + sugar = g->sugar; + + for ( d = 0; g; ) { + for ( u = 0, i = 0; i < n; i++ ) { + if ( dpm_redble(g,p = ps[wb[i]]) ) { + dpm_red(d,g,p,&t,&u,&dmy,&dmy1); + psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar; + sugar = MAX(sugar,psugar); + if ( !u ) { + if ( d ) + d->sugar = sugar; + *rp = d; return; + } + d = t; + break; + } + } + if ( u ) { + g = u; + if ( d ) { + if ( multiple && HMAG(d) > hmag ) { + dpm_ptozp2(d,g,&t,&u); d = t; g = u; + hmag = multiple*HMAG(d); + } + } else { + if ( multiple && HMAG(g) > hmag ) { + dpm_ptozp(g,&t); g = t; + hmag = multiple*HMAG(g); + } + } + } + else if ( !full ) { + if ( g ) { + MKDPM(g->nv,BDY(g),t); t->sugar = sugar; g = t; + } + *rp = g; return; + } else { + m = BDY(g); NEWDMM(mr); mr->dl = m->dl; mr->c = m->c; mr->pos = m->pos; + NEXT(mr) = 0; MKDPM(g->nv,mr,t); t->sugar = mr->dl->td; + adddpm(CO,d,t,&s); d = s; + dpm_rest(g,&t); g = t; + } + } + if ( d ) + d->sugar = sugar; + *rp = d; +} + /* nf computation over a field */ void dp_nf_f(NODE b,DP g,DP *ps,int full,DP *rp) @@ -1412,6 +1734,57 @@ void dp_nf_f(NODE b,DP g,DP *ps,int full,DP *rp) *rp = d; } +void dpm_nf_f(NODE b,DPM g,DPM *ps,int full,DPM *rp) +{ + DPM u,p,d,s,t; + NODE l; + DMM m,mr; + int i,n; + int *wb; + int sugar,psugar; + + if ( !g ) { + *rp = 0; return; + } + for ( n = 0, l = b; l; l = NEXT(l), n++ ); + wb = (int *)ALLOCA(n*sizeof(int)); + for ( i = 0, l = b; i < n; l = NEXT(l), i++ ) + wb[i] = QTOS((Q)BDY(l)); + + sugar = g->sugar; + for ( d = 0; g; ) { + for ( u = 0, i = 0; i < n; i++ ) { + if ( dpm_redble(g,p = ps[wb[i]]) ) { + dpm_red_f(g,p,&u); + psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar; + sugar = MAX(sugar,psugar); + if ( !u ) { + if ( d ) + d->sugar = sugar; + *rp = d; return; + } + break; + } + } + if ( u ) + g = u; + else if ( !full ) { + if ( g ) { + MKDPM(g->nv,BDY(g),t); t->sugar = sugar; g = t; + } + *rp = g; return; + } else { + m = BDY(g); NEWDMM(mr); mr->dl = m->dl; mr->c = m->c; mr->pos = m->pos; + NEXT(mr) = 0; MKDPM(g->nv,mr,t); t->sugar = mr->dl->td; + adddpm(CO,d,t,&s); d = s; + dpm_rest(g,&t); g = t; + } + } + if ( d ) + d->sugar = sugar; + *rp = d; +} + /* nf computation over GF(mod) (only for internal use) */ void dp_nf_mod(NODE b,DP g,DP *ps,int mod,int full,DP *rp) @@ -1591,8 +1964,8 @@ void dp_lnf_f(DP p1,DP p2,NODE g,DP *r1p,DP *r2p) b2 = (DP)BDY(NEXT(b)); divr(CO,(Obj)ONE,(Obj)BDY(b1)->c,&c1); mulr(CO,c1,(Obj)BDY(r1)->c,&c2); chsgnr(c2,&c); - muldc(CO,b1,(P)c,&t); addd(CO,r1,t,&s); r1 = s; - muldc(CO,b2,(P)c,&t); addd(CO,r2,t,&s); r2 = s; + muldc(CO,b1,(Obj)c,&t); addd(CO,r1,t,&s); r1 = s; + muldc(CO,b2,(Obj)c,&t); addd(CO,r2,t,&s); r2 = s; } } *r1p = r1; *r2p = r2; @@ -1643,7 +2016,7 @@ void dp_nf_tab_mod(DP p,LIST *tab,int mod,DP *rp) h = m->dl; while ( !dl_equal(n,h,BDY((DP)BDY(BDY(tab[i])))->dl ) ) i++; - mulmdc(CO,mod,(DP)BDY(NEXT(BDY(tab[i]))),m->c,&t); + mulmdc(CO,mod,(DP)BDY(NEXT(BDY(tab[i]))),(P)m->c,&t); addmd(CO,mod,s,t,&u); s = u; } *rp = s; @@ -1678,10 +2051,11 @@ void dp_nf_tab_f(DP p,LIST *tab,DP *rp) int create_order_spec(VL vl,Obj obj,struct order_spec **specp) { - int i,j,n,s,row,col,ret; + int i,j,n,s,row,col,ret,wlen; struct order_spec *spec; struct order_pair *l; - NODE node,t,tn; + Obj wp,wm; + NODE node,t,tn,wpair; MAT m; VECT v; pointer **b,*bv; @@ -1700,18 +2074,71 @@ int create_order_spec(VL vl,Obj obj,struct order_spec spec->ord.simple = QTOS((Q)obj); return 1; } else if ( OID(obj) == O_LIST ) { + /* module order; obj = [0|1,w,ord] or [0|1,ord] */ node = BDY((LIST)obj); - for ( n = 0, t = node; t; t = NEXT(t), n++ ); - l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair)); - for ( i = 0, t = node, s = 0; i < n; t = NEXT(t), i++ ) { - tn = BDY((LIST)BDY(t)); l[i].order = QTOS((Q)BDY(tn)); - tn = NEXT(tn); l[i].length = QTOS((Q)BDY(tn)); - s += l[i].length; - } - spec->id = 1; spec->obj = obj; - spec->ord.block.order_pair = l; - spec->ord.block.length = n; spec->nv = s; - return 1; + if ( !BDY(node) || NUM(BDY(node)) ) { + switch ( length(node) ) { + case 2: + create_order_spec(0,(Obj)BDY(NEXT(node)),&spec); + spec->id += 256; spec->obj = obj; + spec->top_weight = 0; + spec->module_rank = 0; + spec->module_top_weight = 0; + spec->ispot = (BDY(node)!=0); + if ( spec->ispot ) { + n = QTOS((Q)BDY(node)); + if ( n < 0 ) + spec->pot_nelim = -n; + else + spec->pot_nelim = 0; + } + break; + + case 3: + create_order_spec(0,(Obj)BDY(NEXT(NEXT(node))),&spec); + spec->id += 256; spec->obj = obj; + spec->ispot = (BDY(node)!=0); + node = NEXT(node); + if ( !BDY(node) || OID(BDY(node)) != O_LIST ) + error("create_order_spec : [weight_for_poly,weight_for_modlue] must be specified as a module topweight"); + wpair = BDY((LIST)BDY(node)); + if ( length(wpair) != 2 ) + error("create_order_spec : [weight_for_poly,weight_for_modlue] must be specified as a module topweight"); + + wp = BDY(wpair); + wm = BDY(NEXT(wpair)); + if ( !wp || OID(wp) != O_LIST || !wm || OID(wm) != O_LIST ) + error("create_order_spec : [weight_for_poly,weight_for_modlue] must be specified as a module topweight"); + spec->nv = length(BDY((LIST)wp)); + spec->top_weight = (int *)MALLOC_ATOMIC(spec->nv*sizeof(int)); + for ( i = 0, t = BDY((LIST)wp); i < spec->nv; t = NEXT(t), i++ ) + spec->top_weight[i] = QTOS((Q)BDY(t)); + + spec->module_rank = length(BDY((LIST)wm)); + spec->module_top_weight = (int *)MALLOC_ATOMIC(spec->module_rank*sizeof(int)); + for ( i = 0, t = BDY((LIST)wm); i < spec->module_rank; t = NEXT(t), i++ ) + spec->module_top_weight[i] = QTOS((Q)BDY(t)); + break; + default: + error("create_order_spec : invalid arguments for module order"); + } + + *specp = spec; + return 1; + } else { + /* block order in polynomial ring */ + for ( n = 0, t = node; t; t = NEXT(t), n++ ); + l = (struct order_pair *)MALLOC_ATOMIC(n*sizeof(struct order_pair)); + for ( i = 0, t = node, s = 0; i < n; t = NEXT(t), i++ ) { + tn = BDY((LIST)BDY(t)); l[i].order = QTOS((Q)BDY(tn)); + tn = NEXT(tn); l[i].length = QTOS((Q)BDY(tn)); + s += l[i].length; + } + spec->id = 1; spec->obj = obj; + spec->ord.block.order_pair = l; + spec->ord.block.length = n; spec->nv = s; + return 1; + } } else if ( OID(obj) == O_MAT ) { m = (MAT)obj; row = m->row; col = m->col; b = BDY(m); w = almat(row,col); @@ -1722,14 +2149,6 @@ int create_order_spec(VL vl,Obj obj,struct order_spec spec->nv = col; spec->ord.matrix.row = row; spec->ord.matrix.matrix = w; return 1; - } else if ( OID(obj) == O_VECT ) { - v = (VECT)obj; bv = BDY(v); - if ( v->len < 2 ) error("create_order_spec : invalid argument"); - create_order_spec(0,(Obj)bv[1],&spec); - spec->id += 256; spec->obj = obj; - spec->istop = bv[0]==0; - *specp = spec; - return 1; } else return 0; } @@ -2124,13 +2543,13 @@ void dp_mod(DP p,int mod,NODE subst,DP *rp) *rp = 0; else { for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - for ( tn = subst, s = m->c; tn; tn = NEXT(tn) ) { + for ( tn = subst, s = (P)m->c; tn; tn = NEXT(tn) ) { v = VR((P)BDY(tn)); tn = NEXT(tn); substp(CO,s,v,(P)BDY(tn),&s1); s = s1; } ptomp(mod,s,&t); if ( t ) { - NEXTMP(mr0,mr); mr->c = t; mr->dl = m->dl; + NEXTMP(mr0,mr); mr->c = (Obj)t; mr->dl = m->dl; } } if ( mr0 ) { @@ -2148,7 +2567,7 @@ void dp_rat(DP p,DP *rp) *rp = 0; else { for ( mr0 = 0, m = BDY(p); m; m = NEXT(m) ) { - NEXTMP(mr0,mr); mptop(m->c,&mr->c); mr->dl = m->dl; + NEXTMP(mr0,mr); mptop((P)m->c,(P *)&mr->c); mr->dl = m->dl; } if ( mr0 ) { NEXT(mr) = 0; MKDP(p->nv,mr0,*rp); (*rp)->sugar = p->sugar; @@ -2193,17 +2612,18 @@ void homogenize_order(struct order_spec *old,int n,str error("homogenize_order : invalid input"); } break; - case 1: + case 1: case 257: length = old->ord.block.length; l = (struct order_pair *) MALLOC_ATOMIC((length+1)*sizeof(struct order_pair)); bcopy((char *)old->ord.block.order_pair,(char *)l,length*sizeof(struct order_pair)); l[length].order = 2; l[length].length = 1; - new->id = 1; new->nv = n+1; + new->id = old->id; new->nv = n+1; new->ord.block.order_pair = l; new->ord.block.length = length+1; + new->ispot = old->ispot; break; - case 2: + case 2: case 258: nv = old->nv; row = old->ord.matrix.row; oldm = old->ord.matrix.matrix; newm = almat(row+1,nv+1); for ( i = 0; i <= nv; i++ ) @@ -2213,10 +2633,11 @@ void homogenize_order(struct order_spec *old,int n,str newm[i+1][j] = oldm[i][j]; newm[i+1][j] = 0; } - new->id = 2; new->nv = nv+1; + new->id = old->id; new->nv = nv+1; new->ord.matrix.row = row+1; new->ord.matrix.matrix = newm; + new->ispot = old->ispot; break; - case 3: + case 3: case 259: onv = old->nv; nnv = onv+1; olen = old->ord.composite.length; @@ -2250,12 +2671,33 @@ void homogenize_order(struct order_spec *old,int n,str (struct sparse_weight *)MALLOC(sizeof(struct sparse_weight)); nwb[i].body.sparse_weight[0].pos = onv; nwb[i].body.sparse_weight[0].value = 1; - new->id = 3; + new->id = old->id; new->nv = nnv; new->ord.composite.length = nlen; new->ord.composite.w_or_b = nwb; + new->ispot = old->ispot; print_composite_order_spec(new); break; + case 256: /* simple module order */ + switch ( old->ord.simple ) { + case 0: + new->id = 256; new->ord.simple = 0; break; + case 1: + l = (struct order_pair *) + MALLOC_ATOMIC(2*sizeof(struct order_pair)); + l[0].length = n; l[0].order = old->ord.simple; + l[1].length = 1; l[1].order = 2; + new->id = 257; + new->ord.block.order_pair = l; + new->ord.block.length = 2; new->nv = n+1; + break; + case 2: + new->id = 256; new->ord.simple = 1; break; + default: + error("homogenize_order : invalid input"); + } + new->ispot = old->ispot; + break; default: error("homogenize_order : invalid input"); } @@ -2276,9 +2718,12 @@ void qltozl(Q *w,int n,Q *dvr) v.id = O_VECT; v.len = n; v.body = (pointer *)w; igcdv(&v,dvr); return; } - c = w[0]; nm = NM(c); dn = INT(c) ? ONEN : DN(c); - for ( i = 1; i < n; i++ ) { - c = w[i]; l1 = INT(c) ? ONEN : DN(c); + for ( i = 0; !w[i]; i++ ); + c = w[i]; nm = NM(c); dn = INT(c) ? ONEN : DN(c); + for ( i++; i < n; i++ ) { + c = w[i]; + if ( !c ) continue; + l1 = INT(c) ? ONEN : DN(c); gcdn(nm,NM(c),&g); nm = g; gcdn(dn,l1,&l2); muln(dn,l1,&l3); divsn(l3,l2,&dn); } @@ -2321,6 +2766,24 @@ int dp_redble(DP p1,DP p2) } } +int dpm_redble(DPM p1,DPM p2) +{ + int i,n; + DL d1,d2; + + if ( BDY(p1)->pos != BDY(p2)->pos ) return 0; + d1 = BDY(p1)->dl; d2 = BDY(p2)->dl; + if ( d1->td < d2->td ) + return 0; + else { + for ( i = 0, n = p1->nv; i < n; i++ ) + if ( d1->d[i] < d2->d[i] ) + return 0; + return 1; + } +} + + void dp_subd(DP p1,DP p2,DP *rp) { int i,n; @@ -2332,7 +2795,7 @@ void dp_subd(DP p1,DP p2,DP *rp) NEWDL(d,n); d->td = d1->td - d2->td; for ( i = 0; i < n; i++ ) d->d[i] = d1->d[i]-d2->d[i]; - NEWMP(m); m->dl = d; m->c = (P)ONE; NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)ONE; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; *rp = s; } @@ -2342,7 +2805,7 @@ void dltod(DL d,int n,DP *rp) MP m; DP s; - NEWMP(m); m->dl = d; m->c = (P)ONE; NEXT(m) = 0; + NEWMP(m); m->dl = d; m->c = (Obj)ONE; NEXT(m) = 0; MKDP(n,m,s); s->sugar = d->td; *rp = s; } @@ -2368,11 +2831,38 @@ void dp_ht(DP p,DP *rp) *rp = 0; else { m = BDY(p); - NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0; + NEWMP(mr); mr->dl = m->dl; mr->c = (Obj)ONE; NEXT(mr) = 0; MKDP(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */ } } +void dpm_hm(DPM p,DPM *rp) +{ + DMM m,mr; + + if ( !p ) + *rp = 0; + else { + m = BDY(p); + NEWDMM(mr); mr->dl = m->dl; mr->c = m->c; mr->pos = m->pos; NEXT(mr) = 0; + MKDPM(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */ + } +} + +void dpm_ht(DPM p,DPM *rp) +{ + DMM m,mr; + + if ( !p ) + *rp = 0; + else { + m = BDY(p); + NEWDMM(mr); mr->dl = m->dl; mr->pos = m->pos; mr->c = (Obj)ONE; NEXT(mr) = 0; + MKDPM(p->nv,mr,*rp); (*rp)->sugar = mr->dl->td; /* XXX */ + } +} + + void dp_rest(DP p,DP *rp) { MP m; @@ -2387,6 +2877,20 @@ void dp_rest(DP p,DP *rp) } } +void dpm_rest(DPM p,DPM *rp) +{ + DMM m; + + m = BDY(p); + if ( !NEXT(m) ) + *rp = 0; + else { + MKDPM(p->nv,NEXT(m),*rp); + if ( *rp ) + (*rp)->sugar = p->sugar; + } +} + DL lcm_of_DL(int nv,DL dl1,DL dl2,DL dl) { register int i, *d1, *d2, *d, td; @@ -2678,7 +3182,7 @@ int dpv_ht(DPV p,DP *h) return -1; } else { m = BDY(e[maxi]); - NEWMP(mr); mr->dl = m->dl; mr->c = (P)ONE; NEXT(mr) = 0; + NEWMP(mr); mr->dl = m->dl; mr->c = (Obj)ONE; NEXT(mr) = 0; MKDP(e[maxi]->nv,mr,*h); (*h)->sugar = mr->dl->td; /* XXX */ return maxi; } @@ -2951,6 +3455,90 @@ NODE compute_essential_df(DP *g,DP *gh,int ng) } MKNODE(r1,0,ri); MKLIST(l,r1); BDY(rt) = (pointer)l; + } + return r; +} + +int comp_bits_divisible(int *a,int *b,int n) +{ + int bpi,i,wi,bi; + + bpi = (sizeof(int)/sizeof(char))*8; + for ( i = 0; i < n; i++ ) { + wi = i/bpi; bi = i%bpi; + if ( !(a[wi]&(1< bb ) return 1; + else if ( ba < bb ) return -1; + } + return 0; +} + +NODE mono_raddec(NODE ideal) +{ + DP p; + int nv,w,i,bpi,di,c,len; + int *d,*s,*u,*new; + NODE t,t1,v,r,rem,prev; + + if( !ideal ) return 0; + p = (DP)BDY(ideal); + nv = NV(p); + bpi = (sizeof(int)/sizeof(char))*8; + w = (nv+(bpi-1))/bpi; + d = p->body->dl->d; + if ( !NEXT(ideal) ) { + for ( t = 0, i = nv-1; i >= 0; i-- ) { + if ( d[i] ) { + s = (int *)CALLOC(w,sizeof(int)); + s[i/bpi] |= 1<<(i%bpi); + MKNODE(t1,s,t); + t = t1; + } + } + return t; + } + rem = mono_raddec(NEXT(ideal)); + r = 0; + len = w*sizeof(int); + u = (int *)CALLOC(w,sizeof(int)); + for ( i = nv-1; i >= 0; i-- ) { + if ( d[i] ) { + for ( t = rem; t; t = NEXT(t) ) { + bcopy((char *)BDY(t),(char *)u,len); + u[i/bpi] |= 1<<(i%bpi); + for ( v = r; v; v = NEXT(v) ) { + if ( comp_bits_divisible(u,(int *)BDY(v),nv) ) break; + } + if ( v ) continue; + for ( v = r, prev = 0; v; v = NEXT(v) ) { + if ( comp_bits_divisible((int *)BDY(v),u,nv) ) { + if ( prev ) NEXT(prev) = NEXT(v); + else r = NEXT(r); + } else prev =v; + } + for ( v = r, prev = 0; v; prev = v, v = NEXT(v) ) { + if ( comp_bits_lex(u,(int *)BDY(v),nv) < 0 ) break; + } + new = (int *)CALLOC(w,sizeof(int)); + bcopy((char *)u,(char *)new,len); + MKNODE(t1,new,v); + if ( prev ) NEXT(prev) = t1; + else r = t1; + } + } } return r; }