=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/dp-supp.c,v retrieving revision 1.51 retrieving revision 1.60 diff -u -p -r1.51 -r1.60 --- OpenXM_contrib2/asir2000/builtin/dp-supp.c 2009/01/05 01:47:30 1.51 +++ OpenXM_contrib2/asir2000/builtin/dp-supp.c 2015/08/08 14:19:41 1.60 @@ -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.50 2009/01/04 10:02:00 noro Exp $ + * $OpenXM: OpenXM_contrib2/asir2000/builtin/dp-supp.c,v 1.59 2013/11/05 02:55:03 noro Exp $ */ #include "ca.h" #include "base.h" @@ -399,6 +399,9 @@ void dp_ptozp_d(DP p,DP *rp) if ( 0 ) fprintf(stderr,"(%d,%d)",p_mag((P)d0)-blen,blen); } +#if defined(__MINGW32__) || defined(__MINGW64__) + fflush(stderr); +#endif } void dp_ptozp2_d(DP p0,DP p1,DP *hp,DP *rp) @@ -651,11 +654,11 @@ 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; @@ -704,11 +707,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; @@ -887,7 +890,7 @@ void dp_red_marked(DP p0,DP p1,DP p2,DP hp2,DP *head,D *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 +908,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 = (P)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; } @@ -1182,7 +1187,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 +1208,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 +1252,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)); @@ -1288,6 +1293,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) @@ -1701,6 +1761,23 @@ int create_order_spec(VL vl,Obj obj,struct order_spec return 1; } else if ( OID(obj) == O_LIST ) { node = BDY((LIST)obj); + if ( !BDY(node) || NUM(BDY(node)) ) { + if ( length(node) < 2 ) + error("create_order_spec : invalid argument"); + create_order_spec(0,(Obj)BDY(NEXT(node)),&spec); + spec->id += 256; spec->obj = obj; + 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; + } + *specp = spec; + return 1; + } + 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++ ) { @@ -1722,14 +1799,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->ispot = (bv[0]==1); - *specp = spec; - return 1; } else return 0; } @@ -2299,9 +2368,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); } @@ -2482,6 +2554,9 @@ void _print_mp(int nv,MP m) fprintf(stderr,">",C(m)); } fprintf(stderr,"\n"); +#if defined(__MINGW32__) || defined(__MINGW64__) + fflush(stderr); +#endif } static int cmp_mp_nvar; @@ -2974,6 +3049,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; }