=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/nd.c,v retrieving revision 1.157 retrieving revision 1.163 diff -u -p -r1.157 -r1.163 --- OpenXM_contrib2/asir2000/engine/nd.c 2009/01/04 05:44:51 1.157 +++ OpenXM_contrib2/asir2000/engine/nd.c 2009/01/05 06:29:46 1.163 @@ -1,4 +1,4 @@ -/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.156 2008/05/23 01:24:21 noro Exp $ */ +/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.162 2009/01/05 02:47:39 noro Exp $ */ #include "nd.h" @@ -46,7 +46,7 @@ static int nd_worb_len; static int nd_found,nd_create,nd_notfirst; static int nmv_adv; static int nd_demand; -static int nd_module,nd_istop,nd_mpos; +static int nd_module,nd_ispot,nd_mpos; NumberField get_numberfield(); UINT *nd_det_compute_bound(NDV **dm,int n,int j); @@ -57,7 +57,8 @@ NDV plain_vect_to_ndv_q(Q *mat,int col,UINT *s0vect); LIST ndvtopl(int mod,VL vl,VL dvl,NDV p,int rank); NDV pltondv(VL vl,VL dvl,LIST p); void pltozpl(LIST l,Q *cont,LIST *pp); -void ndl_lcm_nocheck(UINT *d1,unsigned *d2,UINT *d); +void ndl_max(UINT *d1,unsigned *d2,UINT *d); +pointer GC_malloc_atomic_ignore_off_page(int); extern int Denominator,DP_Multiple; @@ -229,9 +230,9 @@ void ndl_homogenize(UINT *d,UINT *r,int obpe,EPOS oepo } w = TD(d); PUT_EXP(r,nd_nvar-1,weight-w); + if ( nd_module ) MPOS(r) = MPOS(d); TD(r) = weight; if ( nd_blockmask ) ndl_weight_mask(r); - if ( nd_module ) MPOS(r) = MPOS(d); } void ndl_dehomogenize(UINT *d) @@ -363,94 +364,16 @@ void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) d[i] = u; } #endif + if ( nd_module ) MPOS(d) = MPOS(d1); TD(d) = ndl_weight(d); if ( nd_blockmask ) ndl_weight_mask(d); - if ( nd_module ) MPOS(d) = MPOS(d1); } -void ndl_lcm_nocheck(UINT *d1,unsigned *d2,UINT *d) +void ndl_max(UINT *d1,unsigned *d2,UINT *d) { UINT t1,t2,u,u1,u2; int i,j,l; -#if USE_UNROLL - switch ( nd_bpe ) { - case 3: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - t1 = (u1&0x38000000); t2 = (u2&0x38000000); u = t1>t2?t1:t2; - t1 = (u1& 0x7000000); t2 = (u2& 0x7000000); u |= t1>t2?t1:t2; - t1 = (u1& 0xe00000); t2 = (u2& 0xe00000); u |= t1>t2?t1:t2; - t1 = (u1& 0x1c0000); t2 = (u2& 0x1c0000); u |= t1>t2?t1:t2; - t1 = (u1& 0x38000); t2 = (u2& 0x38000); u |= t1>t2?t1:t2; - t1 = (u1& 0x7000); t2 = (u2& 0x7000); u |= t1>t2?t1:t2; - t1 = (u1& 0xe00); t2 = (u2& 0xe00); u |= t1>t2?t1:t2; - t1 = (u1& 0x1c0); t2 = (u2& 0x1c0); u |= t1>t2?t1:t2; - t1 = (u1& 0x38); t2 = (u2& 0x38); u |= t1>t2?t1:t2; - t1 = (u1& 0x7); t2 = (u2& 0x7); u |= t1>t2?t1:t2; - d[i] = u; - } - break; - case 4: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - t1 = (u1&0xf0000000); t2 = (u2&0xf0000000); u = t1>t2?t1:t2; - t1 = (u1& 0xf000000); t2 = (u2& 0xf000000); u |= t1>t2?t1:t2; - t1 = (u1& 0xf00000); t2 = (u2& 0xf00000); u |= t1>t2?t1:t2; - t1 = (u1& 0xf0000); t2 = (u2& 0xf0000); u |= t1>t2?t1:t2; - t1 = (u1& 0xf000); t2 = (u2& 0xf000); u |= t1>t2?t1:t2; - t1 = (u1& 0xf00); t2 = (u2& 0xf00); u |= t1>t2?t1:t2; - t1 = (u1& 0xf0); t2 = (u2& 0xf0); u |= t1>t2?t1:t2; - t1 = (u1& 0xf); t2 = (u2& 0xf); u |= t1>t2?t1:t2; - d[i] = u; - } - break; - case 6: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - t1 = (u1&0x3f000000); t2 = (u2&0x3f000000); u = t1>t2?t1:t2; - t1 = (u1& 0xfc0000); t2 = (u2& 0xfc0000); u |= t1>t2?t1:t2; - t1 = (u1& 0x3f000); t2 = (u2& 0x3f000); u |= t1>t2?t1:t2; - t1 = (u1& 0xfc0); t2 = (u2& 0xfc0); u |= t1>t2?t1:t2; - t1 = (u1& 0x3f); t2 = (u2& 0x3f); u |= t1>t2?t1:t2; - d[i] = u; - } - break; - case 8: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - t1 = (u1&0xff000000); t2 = (u2&0xff000000); u = t1>t2?t1:t2; - t1 = (u1& 0xff0000); t2 = (u2& 0xff0000); u |= t1>t2?t1:t2; - t1 = (u1& 0xff00); t2 = (u2& 0xff00); u |= t1>t2?t1:t2; - t1 = (u1& 0xff); t2 = (u2& 0xff); u |= t1>t2?t1:t2; - d[i] = u; - } - break; - case 16: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - t1 = (u1&0xffff0000); t2 = (u2&0xffff0000); u = t1>t2?t1:t2; - t1 = (u1& 0xffff); t2 = (u2& 0xffff); u |= t1>t2?t1:t2; - d[i] = u; - } - break; - case 32: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - d[i] = u1>u2?u1:u2; - } - break; - default: - for ( i = nd_exporigin; i < nd_wpd; i++ ) { - u1 = d1[i]; u2 = d2[i]; - for ( j = 0, u = 0; j < nd_epw; j++ ) { - t1 = (u1&nd_mask[j]); t2 = (u2&nd_mask[j]); u |= t1>t2?t1:t2; - } - d[i] = u; - } - break; - } -#else for ( i = nd_exporigin; i < nd_wpd; i++ ) { u1 = d1[i]; u2 = d2[i]; for ( j = 0, u = 0; j < nd_epw; j++ ) { @@ -458,10 +381,6 @@ void ndl_lcm_nocheck(UINT *d1,unsigned *d2,UINT *d) } d[i] = u; } -#endif - TD(d) = ndl_weight(d); - if ( nd_blockmask ) ndl_weight_mask(d); - if ( nd_module ) MPOS(d) = MPOS(d1); } int ndl_weight(UINT *d) @@ -480,6 +399,8 @@ int ndl_weight(UINT *d) for ( j = 0; j < nd_epw; j++, u>>=nd_bpe ) t += (u&nd_mask0); } + if ( nd_module && current_module_weight_vector && MPOS(d) ) + t += current_module_weight_vector[MPOS(d)]; return t; } @@ -649,7 +570,7 @@ int ndl_module_grlex_compare(UINT *d1,UINT *d2) { int i; - if ( nd_istop ) { + if ( nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -658,7 +579,7 @@ int ndl_module_grlex_compare(UINT *d1,UINT *d2) for ( i = nd_nvar-1; i >= 0; i-- ) if ( GET_EXP(d1,i) < GET_EXP(d2,i) ) return 1; else if ( GET_EXP(d1,i) > GET_EXP(d2,i) ) return -1; - if ( !nd_istop ) { + if ( !nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -669,7 +590,7 @@ int ndl_module_glex_compare(UINT *d1,UINT *d2) { int i; - if ( nd_istop ) { + if ( nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -678,7 +599,7 @@ int ndl_module_glex_compare(UINT *d1,UINT *d2) for ( i = 0; i < nd_nvar; i++ ) if ( GET_EXP(d1,i) > GET_EXP(d2,i) ) return 1; else if ( GET_EXP(d1,i) < GET_EXP(d2,i) ) return -1; - if ( !nd_istop ) { + if ( !nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -689,14 +610,14 @@ int ndl_module_lex_compare(UINT *d1,UINT *d2) { int i; - if ( nd_istop ) { + if ( nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } for ( i = 0; i < nd_nvar; i++ ) if ( GET_EXP(d1,i) > GET_EXP(d2,i) ) return 1; else if ( GET_EXP(d1,i) < GET_EXP(d2,i) ) return -1; - if ( !nd_istop ) { + if ( !nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -707,12 +628,12 @@ int ndl_module_block_compare(UINT *d1,UINT *d2) { int i,c; - if ( nd_istop ) { + if ( nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } if ( c = ndl_block_compare(d1,d2) ) return c; - if ( !nd_istop ) { + if ( !nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -723,12 +644,12 @@ int ndl_module_matrix_compare(UINT *d1,UINT *d2) { int i,c; - if ( nd_istop ) { + if ( nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } if ( c = ndl_matrix_compare(d1,d2) ) return c; - if ( !nd_istop ) { + if ( !nd_ispot ) { if ( MPOS(d1) < MPOS(d2) ) return 1; else if ( MPOS(d1) > MPOS(d2) ) return -1; } @@ -739,12 +660,12 @@ int ndl_module_composite_compare(UINT *d1,UINT *d2) { int i,c; - if ( nd_istop ) { + if ( nd_ispot ) { if ( MPOS(d1) > MPOS(d2) ) return 1; else if ( MPOS(d1) < MPOS(d2) ) return -1; } if ( c = ndl_composite_compare(d1,d2) ) return c; - if ( !nd_istop ) { + if ( !nd_ispot ) { if ( MPOS(d1) > MPOS(d2) ) return 1; else if ( MPOS(d1) < MPOS(d2) ) return -1; } @@ -806,6 +727,10 @@ INLINE void ndl_add(UINT *d1,UINT *d2,UINT *d) { int i; + if ( nd_module ) { + if ( MPOS(d1) && MPOS(d2) && (MPOS(d1) != MPOS(d2)) ) + error("ndl_add : invalid operation"); + } #if 1 switch ( nd_wpd ) { case 2: @@ -824,11 +749,6 @@ INLINE void ndl_add(UINT *d1,UINT *d2,UINT *d) #else for ( i = 0; i < nd_wpd; i++ ) d[i] = d1[i]+d2[i]; #endif - if ( nd_module ) { - if ( MPOS(d1) && MPOS(d2) && (MPOS(d1) != MPOS(d2)) ) - error("ndl_add : invalid operation"); - MPOS(d) = MPOS(d1); - } } /* d1 += d2 */ @@ -836,6 +756,10 @@ INLINE void ndl_addto(UINT *d1,UINT *d2) { int i; + if ( nd_module ) { + if ( MPOS(d1) && MPOS(d2) && (MPOS(d1) != MPOS(d2)) ) + error("ndl_addto : invalid operation"); + } #if 1 switch ( nd_wpd ) { case 2: @@ -854,11 +778,6 @@ INLINE void ndl_addto(UINT *d1,UINT *d2) #else for ( i = 0; i < nd_wpd; i++ ) d1[i] += d2[i]; #endif - if ( nd_module ) { - if ( MPOS(d1) && MPOS(d2) && (MPOS(d1) != MPOS(d2)) ) - error("ndl_addto : invalid operation"); - MPOS(d1) = MPOS(d2); - } } INLINE void ndl_sub(UINT *d1,UINT *d2,UINT *d) @@ -1972,10 +1891,10 @@ again: } if ( nd_demand ) for ( t = g; t; t = NEXT(t) ) - BDY(t) = (pointer)ndv_load((int)BDY(t)); + BDY(t) = (pointer)ndv_load((long)BDY(t)); else for ( t = g; t; t = NEXT(t) ) - BDY(t) = (pointer)nd_ps[(int)BDY(t)]; + BDY(t) = (pointer)nd_ps[(long)BDY(t)]; if ( !checkonly && DP_Print ) { printf("nd_gb done.\n"); fflush(stdout); } return g; } @@ -2152,10 +2071,10 @@ again: } if ( nd_demand ) for ( t = g; t; t = NEXT(t) ) - BDY(t) = (pointer)ndv_load((int)BDY(t)); + BDY(t) = (pointer)ndv_load((long)BDY(t)); else for ( t = g; t; t = NEXT(t) ) - BDY(t) = (pointer)nd_ps_trace[(int)BDY(t)]; + BDY(t) = (pointer)nd_ps_trace[(long)BDY(t)]; if ( nd_nalg ) { print_eg("monic",&eg_monic); print_eg("invdalg",&eg_invdalg); @@ -2262,10 +2181,10 @@ ND_pairs nd_newpairs( NODE g, int t ) dl = DL(nd_psh[t]); ts = SG(nd_psh[t]) - TD(dl); for ( r0 = 0, h = g; h; h = NEXT(h) ) { - if ( nd_module && (MPOS(DL(nd_psh[(int)BDY(h)])) != MPOS(dl)) ) + if ( nd_module && (MPOS(DL(nd_psh[(long)BDY(h)])) != MPOS(dl)) ) continue; NEXTND_pairs(r0,r); - r->i1 = (int)BDY(h); + r->i1 = (long)BDY(h); r->i2 = t; ndl_lcm(DL(nd_psh[r->i1]),dl,r->lcm); s = SG(nd_psh[r->i1])-TD(DL(nd_psh[r->i1])); @@ -2290,16 +2209,24 @@ ND_pairs crit_B( ND_pairs d, int s ) lcm = (UINT *)ALLOCA(nd_wpd*sizeof(UINT)); while ( cur ) { tl = cur->lcm; - if ( ndl_reducible(tl,t) - && (ndl_lcm(DL(nd_psh[cur->i1]),t,lcm),!ndl_equal(lcm,tl)) - && (ndl_lcm(DL(nd_psh[cur->i2]),t,lcm),!ndl_equal(lcm,tl)) ) { - remove = cur; - if ( !prev ) { - head = cur = NEXT(cur); - } else { - cur = NEXT(prev) = NEXT(cur); - } - FREENDP(remove); + if ( ndl_reducible(tl,t) ) { + ndl_lcm(DL(nd_psh[cur->i1]),t,lcm); + if ( !ndl_equal(lcm,tl) ) { + ndl_lcm(DL(nd_psh[cur->i2]),t,lcm); + if (!ndl_equal(lcm,tl)) { + remove = cur; + if ( !prev ) { + head = cur = NEXT(cur); + } else { + cur = NEXT(prev) = NEXT(cur); + } + FREENDP(remove); + } else { + prev = cur; cur = NEXT(cur); + } + } else { + prev = cur; cur = NEXT(cur); + } } else { prev = cur; cur = NEXT(cur); } @@ -2423,7 +2350,7 @@ NODE update_base(NODE nd,int ndp) dl = DL(nd_psh[ndp]); for ( head = last = 0, p = nd; p; ) { - dln = DL(nd_psh[(int)BDY(p)]); + dln = DL(nd_psh[(long)BDY(p)]); if ( ndl_reducible( dln, dl ) ) { p = NEXT(p); if ( last ) NEXT(last) = p; @@ -2718,6 +2645,7 @@ void nd_gr(LIST f,LIST v,int m,int f4,struct order_spe NumberField nf; struct order_spec *ord1; + nd_module = 0; if ( !m && Demand ) nd_demand = 1; else nd_demand = 0; @@ -2871,7 +2799,7 @@ void nd_gr_postproc(LIST f,LIST v,int m,struct order_s } } else { for ( t = x; t; t = NEXT(t) ) - BDY(t) = (pointer)nd_ps[(int)BDY(t)]; + BDY(t) = (pointer)nd_ps[(long)BDY(t)]; } x = ndv_reducebase(x); x = ndv_reduceall(m,x); @@ -2895,15 +2823,16 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,int P p,zp; Q dmy; EPOS oepos; - int obpe,oadv,wmax,i,len,cbpe,ishomo,nalg; + int obpe,oadv,wmax,i,len,cbpe,ishomo,nalg,mrank,trank; Alg alpha,dp; P poly; - LIST f1,f2; + LIST f1,f2,zpl; Obj obj; NumberField nf; struct order_spec *ord1; struct oEGT eg_check,eg0,eg1; + nd_module = 0; if ( DP_Multiple ) nd_scale = ((double)DP_Multiple)/(double)(Denominator?Denominator:1); @@ -2942,18 +2871,34 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,int nocheck = 1; } m = trace > 1 ? trace : get_lprime(mindex); + nd_init_ord(ord); + mrank = 0; for ( t = BDY(f), max = 0; t; t = NEXT(t) ) for ( tv = vv; tv; tv = NEXT(tv) ) { - e = getdeg(tv->v,(P)BDY(t)); - max = MAX(e,max); + if ( nd_module ) { + s = BDY((LIST)BDY(t)); + trank = length(s); + mrank = MAX(mrank,trank); + for ( ; s; s = NEXT(s) ) { + e = getdeg(tv->v,(P)BDY(s)); + max = MAX(e,max); + } + } else { + e = getdeg(tv->v,(P)BDY(t)); + max = MAX(e,max); + } } - nd_init_ord(ord); nd_setup_parameters(nvar,max); obpe = nd_bpe; oadv = nmv_adv; oepos = nd_epos; ishomo = 1; for ( in0 = 0, fd0 = 0, t = BDY(f); t; t = NEXT(t) ) { - ptozp((P)BDY(t),1,&dmy,&zp); - c = (pointer)ptondv(CO,vv,zp); + if ( nd_module ) { + pltozpl((LIST)BDY(t),&dmy,&zpl); + c = (pointer)pltondv(CO,vv,zpl); + } else { + ptozp((P)BDY(t),1,&dmy,&zp); + c = (pointer)ptondv(CO,vv,zp); + } if ( ishomo ) ishomo = ishomo && ndv_ishomo(c); if ( c ) { @@ -3024,8 +2969,10 @@ void nd_gr_trace(LIST f,LIST v,int trace,int homo,int /* dp->p */ nd_bpe = cbpe; nd_setup_parameters(nd_nvar,0); - for ( r = cand; r; r = NEXT(r) ) - BDY(r) = (pointer)ndvtop(0,CO,vv,BDY(r)); + for ( r = cand; r; r = NEXT(r) ) { + if ( nd_module ) BDY(r) = ndvtopl(0,CO,vv,BDY(r),mrank); + else BDY(r) = (pointer)ndvtop(0,CO,vv,BDY(r)); + } if ( nalg ) cand = postprocess_algcoef(av,alist,cand); MKLIST(*rp,cand); @@ -3456,7 +3403,7 @@ UINT *ndv_compute_bound(NDV p) len = LEN(p); m = BDY(p); ndl_copy(DL(m),d1); NMV_ADV(m); for ( i = 1; i < len; i++, NMV_ADV(m) ) { - ndl_lcm_nocheck(DL(m),d1,d2); + ndl_max(DL(m),d1,d2); t = d1; d1 = d2; d2 = t; } l = nd_nvar+31; @@ -3836,6 +3783,7 @@ void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta } if ( nd_module ) { mpos = MPOS(d1); + TD(d1) = ndl_weight(d1); if ( MPOS(d0) ) error("weyl_mul_nm_nmv : invalid operation"); } tab[0] = m; @@ -3872,8 +3820,11 @@ void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta TD(d) = s; PUT_EXP(d,n-1,s-h); } else TD(d) = h; + if ( nd_module ) { + MPOS(d) = mpos; + TD(d) = ndl_weight(d); + } if ( nd_blockmask ) ndl_weight_mask(d); - if ( nd_module ) MPOS(d) = mpos; if ( mod ) c = ctab[j]; else q = ctab_q[j]; p = tab+curlen*j; @@ -4213,9 +4164,12 @@ NDV pltondv(VL vl,VL dvl,LIST p) r = 0; for ( i = 1, t = BDY(p); t; t = NEXT(t), i++ ) { ri = ptond(vl,dvl,(P)BDY(t)); - if ( ri ) - for ( m = BDY(ri); m; m = NEXT(m) ) - MPOS(DL(m)) = i; + if ( ri ) + for ( m = BDY(ri); m; m = NEXT(m) ) { + MPOS(DL(m)) = i; + TD(DL(m)) = ndl_weight(DL(m)); + if ( nd_blockmask ) ndl_weight_mask(DL(m)); + } r = nd_add(0,r,ri); } return ndtondv(0,r); @@ -4516,7 +4470,7 @@ void nd_init_ord(struct order_spec *ord) /* module order */ case 256: - nd_istop = ord->istop; + nd_ispot = ord->ispot; nd_dcomp = -1; nd_isrlex = 0; switch ( ord->ord.simple ) { @@ -5226,7 +5180,7 @@ NODE nd_f4(int m) g = update_base(g,nh); } } - for ( r = g; r; r = NEXT(r) ) BDY(r) = (pointer)nd_ps[(int)BDY(r)]; + for ( r = g; r; r = NEXT(r) ) BDY(r) = (pointer)nd_ps[(long)BDY(r)]; #if 0 fprintf(asir_out,"ndv_alloc=%d\n",ndv_alloc); #endif @@ -5326,7 +5280,7 @@ NODE nd_f4_trace(int m) g = update_base(g,nh); } } - for ( r = g; r; r = NEXT(r) ) BDY(r) = (pointer)nd_ps_trace[(int)BDY(r)]; + for ( r = g; r; r = NEXT(r) ) BDY(r) = (pointer)nd_ps_trace[(long)BDY(r)]; #if 0 fprintf(asir_out,"ndv_alloc=%d\n",ndv_alloc); #endif @@ -6479,9 +6433,9 @@ int nd_monic(int mod,ND *p) PUT_EXP(DL(mr),i+nd_ntrans,e); td += MUL_WEIGHT(e,i+nd_ntrans); } + if ( nd_module ) MPOS(DL(mr)) = MPOS(DL(m)); TD(DL(mr)) = td; if ( nd_blockmask) ndl_weight_mask(DL(mr)); - if ( nd_module ) MPOS(DL(mr)) = MPOS(DL(m)); } } NEXT(mr) = 0;