=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/nd.c,v retrieving revision 1.28 retrieving revision 1.31 diff -u -p -r1.28 -r1.31 --- OpenXM_contrib2/asir2000/engine/nd.c 2003/08/11 06:58:01 1.28 +++ OpenXM_contrib2/asir2000/engine/nd.c 2003/08/13 03:13:22 1.31 @@ -1,4 +1,4 @@ -/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.27 2003/08/10 01:31:24 noro Exp $ */ +/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.30 2003/08/12 09:07:19 noro Exp $ */ #include "ca.h" #include "inline.h" @@ -23,14 +23,15 @@ typedef struct oPGeoBucket { typedef struct oND { struct oNM *body; int nv; + int len; int sugar; } *ND; typedef struct oNDV { struct oNMV *body; int nv; - int sugar; int len; + int sugar; } *NDV; typedef struct oNM { @@ -109,7 +110,7 @@ extern int Top,Reverse; ((r)=(RHist)MALLOC(sizeof(struct oRHist)+(nd_wpd-1)*sizeof(unsigned int))) #define NEWND_pairs(m) if(!_ndp_free_list)_NDP_alloc(); (m)=_ndp_free_list; _ndp_free_list = NEXT(_ndp_free_list) #define NEWNM(m) if(!_nm_free_list)_NM_alloc(); (m)=_nm_free_list; _nm_free_list = NEXT(_nm_free_list) -#define MKND(n,m,d) if(!_nd_free_list)_ND_alloc(); (d)=_nd_free_list; _nd_free_list = (ND)BDY(_nd_free_list); (d)->nv=(n); BDY(d)=(m) +#define MKND(n,m,len,d) if(!_nd_free_list)_ND_alloc(); (d)=_nd_free_list; _nd_free_list = (ND)BDY(_nd_free_list); NV(d)=(n); LEN(d)=(len); BDY(d)=(m) #define NEXTRHist(r,c) \ if(!(r)){NEWRHist(r);(c)=(r);}else{NEWRHist(NEXT(c));(c)=NEXT(c);} @@ -134,7 +135,7 @@ void nd_mul_c_q(ND p,Q mul); ND normalize_pbucket(int mod,PGeoBucket g); int head_pbucket(int mod,PGeoBucket g); int head_pbucket_q(PGeoBucket g); -void add_pbucket(int mod,PGeoBucket g,ND d,int l); +void add_pbucket(int mod,PGeoBucket g,ND d); void free_pbucket(PGeoBucket b); void mulq_pbucket(PGeoBucket g,Q c); PGeoBucket create_pbucket(); @@ -175,6 +176,9 @@ int nd_sp(int mod,ND_pairs p,ND *nf); int nd_find_reducer(ND g); int nd_find_reducer_direct(ND g,NDV *ps,int len); int nd_nf(int mod,ND g,int full,ND *nf); +int nd_nf_pbucket(int mod,ND g,int full,ND *nf); +int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full,ND *rp); +int nd_nf_direct_pbucket(int mod,ND g,NDV *ps,int len,int full,ND *rp); ND nd_reduce(ND p1,ND p2); ND nd_reduce_special(ND p1,ND p2); NODE nd_reduceall(int m,NODE f); @@ -257,7 +261,7 @@ void _NDP_alloc() } } -INLINE nd_length(ND p) +INLINE int nd_length(ND p) { NM m; int i; @@ -340,7 +344,7 @@ int ndl_dehomogenize(unsigned int *d) { unsigned int mask; unsigned int h; - int i; + int i,bits; if ( is_rlex ) { if ( nd_bpe == 32 ) { @@ -350,7 +354,8 @@ int ndl_dehomogenize(unsigned int *d) d[i-1] = 0; return h; } else { - mask = (1<<(nd_epw*nd_bpe))-1; + bits = nd_epw*nd_bpe; + mask = bits==32?0xffffffff:((1<<(nd_epw*nd_bpe))-1); h = (d[0]>>((nd_epw-1)*nd_bpe))&nd_mask0; for ( i = 0; i < nd_wpd; i++ ) d[i] = ((d[i]<index; } } - if ( Reverse ) for ( i = nd_psn-1; i >= 0; i-- ) { r = nd_psh[i]; @@ -781,8 +785,9 @@ INLINE int nd_find_reducer(ND g) INLINE int nd_find_reducer_direct(ND g,NDV *ps,int len) { - int i; NDV r; + RHist s; + int d,k,i; if ( Reverse ) for ( i = len-1; i >= 0; i-- ) { @@ -806,7 +811,7 @@ INLINE int nd_find_reducer_direct(ND g,NDV *ps,int len ND nd_add(int mod,ND p1,ND p2) { int n,c; - int t; + int t,can; ND r; NM m1,m2,mr0,mr,s; @@ -817,6 +822,7 @@ ND nd_add(int mod,ND p1,ND p2) else if ( !mod ) return nd_add_q(p1,p2); else { + can = 0; for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; ) { if ( TD(m1) > TD(m2) ) c = 1; @@ -831,8 +837,10 @@ ND nd_add(int mod,ND p1,ND p2) t += mod; s = m1; m1 = NEXT(m1); if ( t ) { + can++; NEXTNM2(mr0,mr,s); CM(mr) = (t); } else { + can += 2; FREENM(s); } s = m2; m2 = NEXT(m2); FREENM(s); @@ -860,6 +868,7 @@ ND nd_add(int mod,ND p1,ND p2) NEXT(mr) = 0; BDY(p1) = mr0; SG(p1) = MAX(SG(p1),SG(p2)); + LEN(p1) = LEN(p1)+LEN(p2)-can; FREEND(p2); return p1; } @@ -867,7 +876,7 @@ ND nd_add(int mod,ND p1,ND p2) ND nd_add_q(ND p1,ND p2) { - int n,c; + int n,c,can; ND r; NM m1,m2,mr0,mr,s; Q t; @@ -877,6 +886,7 @@ ND nd_add_q(ND p1,ND p2) else if ( !p2 ) return p1; else { + can = 0; for ( n = NV(p1), m1 = BDY(p1), m2 = BDY(p2), mr0 = 0; m1 && m2; ) { if ( TD(m1) > TD(m2) ) c = 1; @@ -889,8 +899,10 @@ ND nd_add_q(ND p1,ND p2) addq(CQ(m1),CQ(m2),&t); s = m1; m1 = NEXT(m1); if ( t ) { + can++; NEXTNM2(mr0,mr,s); CQ(mr) = (t); } else { + can += 2; FREENM(s); } s = m2; m2 = NEXT(m2); FREENM(s); @@ -918,12 +930,12 @@ ND nd_add_q(ND p1,ND p2) NEXT(mr) = 0; BDY(p1) = mr0; SG(p1) = MAX(SG(p1),SG(p2)); + LEN(p1) = LEN(p1)+LEN(p2)-can; FREEND(p2); return p1; } } -#if !USE_GEOBUCKET /* ret=1 : success, ret=0 : overflow */ int nd_nf(int mod,ND g,int full,ND *rp) { @@ -931,7 +943,7 @@ int nd_nf(int mod,ND g,int full,ND *rp) NM m,mrd,tail; NM mul; int n,sugar,psugar,sugar0,stat,index; - int c,c1,c2; + int c,c1,c2,dummy; RHist h; NDV p,red; Q cg,cred,gcd; @@ -986,14 +998,16 @@ afo: m = BDY(g); if ( NEXT(m) ) { BDY(g) = NEXT(m); NEXT(m) = 0; + LEN(g)--; } else { FREEND(g); g = 0; } if ( d ) { NEXT(tail)=m; tail=m; + LEN(d)++; } else { - MKND(n,m,d); + MKND(n,m,1,d); tail = BDY(d); } } @@ -1004,89 +1018,13 @@ afo: return 1; } -int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full,ND *rp) +int nd_nf_pbucket(int mod,ND g,int full,ND *rp) { - ND d; - NM m,mrd,tail; - NM mul; - int n,sugar,psugar,sugar0,stat,index; - int c,c1,c2; - RHist h; - NDV p,red; - Q cg,cred,gcd; - double hmag; - - if ( !g ) { - *rp = 0; - return 1; - } -#if 0 - if ( !mod ) - hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; -#else - /* XXX */ - hmag = 0; -#endif - - sugar0 = sugar = SG(g); - n = NV(g); - mul = (NM)ALLOCA(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int)); - for ( d = 0; g; ) { - index = nd_find_reducer_direct(g,ps,len); - if ( index >= 0 ) { - p = ps[index]; - ndl_sub(HDL(g),HDL(p),DL(mul)); - TD(mul) = HTD(g)-HTD(p); - if ( ndl_check_bound2_direct(HDL(p),DL(mul)) ) { - nd_free(g); nd_free(d); - return 0; - } - if ( mod ) { - c1 = invm(HCM(p),mod); c2 = mod-HCM(g); - DMAR(c1,c2,0,mod,c); CM(mul) = c; - } else { - igcd_cofactor(HCQ(g),HCQ(p),&gcd,&cg,&cred); - chsgnq(cg,&CQ(mul)); - nd_mul_c_q(d,cred); nd_mul_c_q(g,cred); - } - g = nd_add(mod,g,ndv_mul_nm(mod,p,mul)); - sugar = MAX(sugar,SG(p)+TD(mul)); - if ( !mod && hmag && g && ((double)(p_mag((P)HCQ(g))) > hmag) ) { - nd_removecont2(d,g); - hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; - } - } else if ( !full ) { - *rp = g; - return 1; - } else { - m = BDY(g); - if ( NEXT(m) ) { - BDY(g) = NEXT(m); NEXT(m) = 0; - } else { - FREEND(g); g = 0; - } - if ( d ) { - NEXT(tail)=m; - tail=m; - } else { - MKND(n,m,d); - tail = BDY(d); - } - } - } - if ( d ) - SG(d) = sugar; - *rp = d; - return 1; -} -#else -int nd_nf(int mod,ND g,int full,ND *rp) -{ int hindex,index; NDV p; ND u,d,red; NODE l; - NM mul,m,mrd; + NM mul,m,mrd,tail; int sugar,psugar,n,h_reducible; PGeoBucket bucket; int c,c1,c2; @@ -1103,7 +1041,7 @@ int nd_nf(int mod,ND g,int full,ND *rp) if ( !mod ) hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; bucket = create_pbucket(); - add_pbucket(mod,bucket,g,nd_length(g)); + add_pbucket(mod,bucket,g); d = 0; mul = (NM)ALLOCA(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int)); while ( 1 ) { @@ -1142,7 +1080,7 @@ int nd_nf(int mod,ND g,int full,ND *rp) red = ndv_mul_nm(mod,p,mul); bucket->body[hindex] = nd_remove_head(g); red = nd_remove_head(red); - add_pbucket(mod,bucket,red,LEN(p)-1); + add_pbucket(mod,bucket,red); sugar = MAX(sugar,SG(p)+TD(mul)); if ( !mod && hmag && (gmag > hmag) ) { g = normalize_pbucket(mod,bucket); @@ -1154,7 +1092,7 @@ int nd_nf(int mod,ND g,int full,ND *rp) } nd_removecont2(d,g); hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; - add_pbucket(mod,bucket,g,nd_length(g)-1); + add_pbucket(mod,bucket,g); } } else if ( !full ) { g = normalize_pbucket(mod,bucket); @@ -1166,16 +1104,19 @@ int nd_nf(int mod,ND g,int full,ND *rp) m = BDY(g); if ( NEXT(m) ) { BDY(g) = NEXT(m); NEXT(m) = 0; + LEN(g)--; } else { FREEND(g); g = 0; } bucket->body[hindex] = g; NEXT(m) = 0; if ( d ) { - for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) ); - NEXT(mrd) = m; + NEXT(tail)=m; + tail=m; + LEN(d)++; } else { - MKND(n,m,d); + MKND(n,m,1,d); + tail = BDY(d); } } } @@ -1183,11 +1124,89 @@ int nd_nf(int mod,ND g,int full,ND *rp) int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full,ND *rp) { + ND d; + NM m,mrd,tail; + NM mul; + int n,sugar,psugar,sugar0,stat,index; + int c,c1,c2; + RHist h; + NDV p,red; + Q cg,cred,gcd; + double hmag; + + if ( !g ) { + *rp = 0; + return 1; + } +#if 0 + if ( !mod ) + hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; +#else + /* XXX */ + hmag = 0; +#endif + + sugar0 = sugar = SG(g); + n = NV(g); + mul = (NM)ALLOCA(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int)); + for ( d = 0; g; ) { + index = nd_find_reducer_direct(g,ps,len); + if ( index >= 0 ) { + p = ps[index]; + ndl_sub(HDL(g),HDL(p),DL(mul)); + TD(mul) = HTD(g)-HTD(p); + if ( ndl_check_bound2_direct(HDL(p),DL(mul)) ) { + nd_free(g); nd_free(d); + return 0; + } + if ( mod ) { + c1 = invm(HCM(p),mod); c2 = mod-HCM(g); + DMAR(c1,c2,0,mod,c); CM(mul) = c; + } else { + igcd_cofactor(HCQ(g),HCQ(p),&gcd,&cg,&cred); + chsgnq(cg,&CQ(mul)); + nd_mul_c_q(d,cred); nd_mul_c_q(g,cred); + } + g = nd_add(mod,g,ndv_mul_nm(mod,p,mul)); + sugar = MAX(sugar,SG(p)+TD(mul)); + if ( !mod && hmag && g && ((double)(p_mag((P)HCQ(g))) > hmag) ) { + nd_removecont2(d,g); + hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; + } + } else if ( !full ) { + *rp = g; + return 1; + } else { + m = BDY(g); + if ( NEXT(m) ) { + BDY(g) = NEXT(m); NEXT(m) = 0; + LEN(g)--; + } else { + FREEND(g); g = 0; + } + if ( d ) { + NEXT(tail)=m; + tail=m; + LEN(d)++; + } else { + MKND(n,m,1,d); + tail = BDY(d); + } + } + } + if ( d ) + SG(d) = sugar; + *rp = d; + return 1; +} + +int nd_nf_direct_pbucket(int mod,ND g,NDV *ps,int len,int full,ND *rp) +{ int hindex,index; NDV p; ND u,d,red; NODE l; - NM mul,m,mrd; + NM mul,m,mrd,tail; int sugar,psugar,n,h_reducible; PGeoBucket bucket; int c,c1,c2; @@ -1209,7 +1228,7 @@ int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full hmag = 0; #endif bucket = create_pbucket(); - add_pbucket(mod,bucket,g,nd_length(g)); + add_pbucket(mod,bucket,g); d = 0; mul = (NM)ALLOCA(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int)); while ( 1 ) { @@ -1246,7 +1265,7 @@ int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full red = ndv_mul_nm(mod,p,mul); bucket->body[hindex] = nd_remove_head(g); red = nd_remove_head(red); - add_pbucket(mod,bucket,red,LEN(p)-1); + add_pbucket(mod,bucket,red); sugar = MAX(sugar,SG(p)+TD(mul)); if ( !mod && hmag && (gmag > hmag) ) { g = normalize_pbucket(mod,bucket); @@ -1258,7 +1277,7 @@ int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full } nd_removecont2(d,g); hmag = ((double)p_mag((P)HCQ(g)))*nd_scale; - add_pbucket(mod,bucket,g,nd_length(g)-1); + add_pbucket(mod,bucket,g); } } else if ( !full ) { g = normalize_pbucket(mod,bucket); @@ -1270,21 +1289,23 @@ int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full m = BDY(g); if ( NEXT(m) ) { BDY(g) = NEXT(m); NEXT(m) = 0; + LEN(g)--; } else { FREEND(g); g = 0; } bucket->body[hindex] = g; NEXT(m) = 0; if ( d ) { - for ( mrd = BDY(d); NEXT(mrd); mrd = NEXT(mrd) ); - NEXT(mrd) = m; + NEXT(tail)=m; + tail=m; + LEN(d)++; } else { - MKND(n,m,d); + MKND(n,m,1,d); + tail = BDY(d); } } } } -#endif /* input : list of DP, cand : list of DP */ @@ -1295,6 +1316,7 @@ int nd_check_candidate(NODE input,NODE cand) NODE t; nd_setup(0,cand); + n = length(cand); /* membercheck : list is a subset of Id(cand) ? */ for ( t = input; t; t = NEXT(t) ) { @@ -1320,8 +1342,10 @@ ND nd_remove_head(ND p) if ( !NEXT(m) ) { FREEND(p); p = 0; - } else + } else { BDY(p) = NEXT(m); + LEN(p)--; + } FREENM(m); return p; } @@ -1346,14 +1370,17 @@ void free_pbucket(PGeoBucket b) { GC_free(b); } -void add_pbucket(int mod,PGeoBucket g,ND d,int l) +void add_pbucket(int mod,PGeoBucket g,ND d) { - int k,m; + int l,i,k,m; - for ( k = 0, m = 1; l > m; k++, m <<= 2 ); - /* 4^(k-1) < l <= 4^k */ + if ( !d ) + return; + l = LEN(d); + for ( k = 0, m = 1; l > m; k++, m <<= 1 ); + /* 2^(k-1) < l <= 2^k (=m) */ d = nd_add(mod,g->body[k],d); - for ( ; d && nd_length(d) > 1<<(2*k); k++ ) { + for ( ; d && LEN(d) > m; k++, m <<= 1 ) { g->body[k] = 0; d = nd_add(mod,g->body[k+1],d); } @@ -1519,7 +1546,7 @@ again: d = nd_reconstruct(m,0,d); goto again; } - stat = nd_nf(m,h,!Top,&nf); + stat = m?nd_nf_pbucket(m,h,!Top,&nf):nd_nf(m,h,!Top,&nf); if ( !stat ) { NEXT(l) = d; d = l; d = nd_reconstruct(m,0,d); @@ -1572,7 +1599,7 @@ again: d = nd_reconstruct(m,1,d); goto again; } - stat = nd_nf(m,h,!Top,&nf); + stat = nd_nf_pbucket(m,h,!Top,&nf); if ( !stat ) { NEXT(l) = d; d = l; d = nd_reconstruct(m,1,d); @@ -1914,13 +1941,6 @@ ND_pairs nd_minp( ND_pairs d, ND_pairs *prest ) c = ndl_compare(p->lcm,lcm); if ( c < 0 ) goto find; -#if 0 - else if ( c == 0 ) { - tlen = nd_psl[p->i1]+nd_psl[p->i2]; - if ( tlen < len ) - goto find; - } -#endif } } continue; @@ -2274,13 +2294,13 @@ ND dptond(int mod,DP p) ND d; NM m0,m; MP t; - int n; + int n,l; if ( !p ) return 0; n = NV(p); m0 = 0; - for ( t = BDY(p); t; t = NEXT(t) ) { + for ( t = BDY(p), l = 0; t; t = NEXT(t), l++ ) { NEXTNM(m0,m); if ( mod ) CM(m) = ITOS(C(t)); @@ -2290,8 +2310,7 @@ ND dptond(int mod,DP p) dltondl(n,DL(t),DL(m)); } NEXT(m) = 0; - MKND(n,m0,d); - NV(d) = n; + MKND(n,m0,l,d); SG(d) = SG(p); return d; } @@ -2746,6 +2765,17 @@ void nd_reconstruct_direct(int mod,NDV *ps,int len) if ( mod != 0 ) for ( i = len-1; i >= 0; i-- ) ndv_realloc(ps[i],obpe,oadv); + for ( i = 0; i < REDTAB_LEN; i++ ) { + for ( mr0 = 0, r = nd_red[i]; r; r = NEXT(r) ) { + NEXTRHist(mr0,mr); + mr->index = r->index; + SG(mr) = SG(r); + TD(mr) = TD(r); + ndl_dup(obpe,DL(r),DL(mr)); + } + if ( mr0 ) NEXT(mr) = 0; + nd_red[i] = mr0; + } prev_nm_free_list = 0; prev_ndp_free_list = 0; GC_gcollect(); @@ -2809,7 +2839,7 @@ ND nd_copy(ND p) ndl_copy(DL(m),DL(mr)); } NEXT(mr) = 0; - MKND(NV(p),mr0,r); + MKND(NV(p),mr0,LEN(p),r); SG(r) = SG(p); return r; } @@ -2916,7 +2946,7 @@ ND ndv_mul_nm(int mod,NDV p,NM m0) } } NEXT(mr) = 0; - MKND(NV(p),mr0,r); + MKND(NV(p),mr0,len,r); SG(r) = SG(p) + td; return r; } @@ -2961,7 +2991,7 @@ NDV ndtondv(int mod,ND p) if ( !p ) return 0; - len = nd_length(p); + len = LEN(p); if ( mod ) m0 = m = (NMV)MALLOC_ATOMIC(len*nmv_adv); else @@ -2994,7 +3024,7 @@ ND ndvtond(int mod,NDV p) CQ(m) = CQ(t); } NEXT(m) = 0; - MKND(NV(p),m0,d); + MKND(NV(p),m0,len,d); SG(d) = SG(p); return d; }