===================================================================
RCS file: /home/cvs/OpenXM_contrib2/asir2000/engine/nd.c,v
retrieving revision 1.26
retrieving revision 1.32
diff -u -p -r1.26 -r1.32
--- OpenXM_contrib2/asir2000/engine/nd.c 2003/08/07 09:47:08 1.26
+++ OpenXM_contrib2/asir2000/engine/nd.c 2003/08/19 04:38:09 1.32
@@ -1,4 +1,4 @@
-/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.25 2003/08/07 08:46:50 noro Exp $ */
+/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.31 2003/08/13 03:13:22 noro Exp $ */
#include "ca.h"
#include "inline.h"
@@ -11,6 +11,8 @@
#define INLINE
#endif
+#define USE_GEOBUCKET 1
+
#define REDTAB_LEN 32003
typedef struct oPGeoBucket {
@@ -21,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 {
@@ -64,10 +67,13 @@ typedef struct oND_pairs {
unsigned int lcm[1];
} *ND_pairs;
+int (*nm_compare_generic)(NM a,NM b);
+int (*nmv_compare_generic)(NMV a,NMV b);
+
double nd_scale=2;
static unsigned int **nd_bound;
int nd_nvar;
-int is_rlex;
+int nd_isrlex;
int nd_epw,nd_bpe,nd_wpd;
unsigned int nd_mask[32];
unsigned int nd_mask0,nd_mask1;
@@ -88,6 +94,7 @@ int nd_red_len;
int nd_found,nd_create,nd_notfirst;
int nm_adv;
int nmv_adv;
+int nd_degcompat;
extern int Top,Reverse;
@@ -102,12 +109,16 @@ extern int Top,Reverse;
#define SG(a) ((a)->sugar)
#define LEN(a) ((a)->len)
+#define NM_COMPARE(m1,m2)\
+nd_degcompat?(TD(m1)>TD(m2)?1\
+ :(TD(m1)
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);}
@@ -122,6 +133,8 @@ if(!(r)){(c)=(r)=(s);}else{NEXT(c)=(s);(c)=(s);}
#define NEXTND_pairs(r,c) \
if(!(r)){NEWND_pairs(r);(c)=(r);}else{NEWND_pairs(NEXT(c));(c)=NEXT(c);}
+void nd_init_ord(struct order_spec *spec);
+int nd_check_candidate(NODE input,NODE cand);
void nd_removecont(int mod,ND p);
void nd_removecont2(ND p1,ND p2);
void ndv_removecont(int mod,NDV p);
@@ -131,7 +144,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();
@@ -139,6 +152,7 @@ ND nd_remove_head(ND p);
void GC_gcollect();
NODE append_one(NODE,int);
+NODE nd_reducebase(NODE x);
void removecont_array(Q *c,int n);
ND_pairs crit_B( ND_pairs d, int s );
@@ -156,7 +170,7 @@ ND_pairs crit_F( ND_pairs d1 );
ND_pairs crit_M( ND_pairs d1 );
ND_pairs nd_newpairs( NODE g, int t );
ND_pairs update_pairs( ND_pairs d, NODE /* of index */ g, int t);
-NODE nd_gb(int m);
+NODE nd_gb(int m,int checkonly);
NODE nd_gb_trace(int m);
void nd_free_private_storage();
void _NM_alloc();
@@ -171,9 +185,14 @@ 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);
+int nd_gbcheck(int m,NODE f);
+int nd_membercheck(int m,NODE f);
void nd_free(ND p);
void ndv_free(NDV p);
void ndl_print(unsigned int *dl);
@@ -251,7 +270,7 @@ void _NDP_alloc()
}
}
-INLINE nd_length(ND p)
+INLINE int nd_length(ND p)
{
NM m;
int i;
@@ -334,9 +353,9 @@ int ndl_dehomogenize(unsigned int *d)
{
unsigned int mask;
unsigned int h;
- int i;
+ int i,bits;
- if ( is_rlex ) {
+ if ( nd_isrlex ) {
if ( nd_bpe == 32 ) {
h = d[0];
for ( i = 1; i < nd_wpd; i++ )
@@ -344,7 +363,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]< *d2 )
- return is_rlex ? -1 : 1;
+ return nd_isrlex ? -1 : 1;
else if ( *d1 < *d2 )
- return is_rlex ? 1 : -1;
+ return nd_isrlex ? 1 : -1;
return 0;
}
@@ -751,7 +771,6 @@ INLINE int nd_find_reducer(ND g)
return r->index;
}
}
-
if ( Reverse )
for ( i = nd_psn-1; i >= 0; i-- ) {
r = nd_psh[i];
@@ -775,8 +794,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-- ) {
@@ -800,7 +820,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;
@@ -811,13 +831,9 @@ 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;
- else if ( TD(m1) < TD(m2) )
- c = -1;
- else
- c = ndl_compare(DL(m1),DL(m2));
+ c = NM_COMPARE(m1,m2);
switch ( c ) {
case 0:
t = ((CM(m1))+(CM(m2))) - mod;
@@ -825,8 +841,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);
@@ -854,6 +872,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;
}
@@ -861,7 +880,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;
@@ -871,20 +890,18 @@ 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;
- else if ( TD(m1) < TD(m2) )
- c = -1;
- else
- c = ndl_compare(DL(m1),DL(m2));
+ c = NM_COMPARE(m1,m2);
switch ( c ) {
case 0:
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);
@@ -912,12 +929,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 1
/* ret=1 : success, ret=0 : overflow */
int nd_nf(int mod,ND g,int full,ND *rp)
{
@@ -925,7 +942,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;
@@ -980,14 +997,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);
}
}
@@ -997,19 +1016,20 @@ afo:
*rp = d;
return 1;
}
-#else
-int nd_nf(int mod,ND g,int full,ND *rp)
+
+int nd_nf_pbucket(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;
Q cg,cred,gcd,zzz;
RHist h;
+ double hmag,gmag;
if ( !g ) {
*rp = 0;
@@ -1017,8 +1037,10 @@ int nd_nf(int mod,ND g,int full,ND *rp)
}
sugar = SG(g);
n = NV(g);
+ 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 ) {
@@ -1052,12 +1074,25 @@ int nd_nf(int mod,ND g,int full,ND *rp)
nd_mul_c_q(d,cred);
mulq_pbucket(bucket,cred);
g = bucket->body[hindex];
+ gmag = (double)p_mag((P)HCQ(g));
}
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);
+ if ( !g ) {
+ if ( d )
+ SG(d) = sugar;
+ *rp = d;
+ return 1;
+ }
+ nd_removecont2(d,g);
+ hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;
+ add_pbucket(mod,bucket,g);
+ }
} else if ( !full ) {
g = normalize_pbucket(mod,bucket);
if ( g )
@@ -1068,21 +1103,23 @@ 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);
}
}
}
}
-#endif
int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full,ND *rp)
{
@@ -1142,14 +1179,16 @@ 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;
}
if ( d ) {
NEXT(tail)=m;
tail=m;
+ LEN(d)++;
} else {
- MKND(n,m,d);
+ MKND(n,m,1,d);
tail = BDY(d);
}
}
@@ -1160,6 +1199,140 @@ int nd_nf_direct(int mod,ND g,NDV *ps,int len,int full
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,tail;
+ int sugar,psugar,n,h_reducible;
+ PGeoBucket bucket;
+ int c,c1,c2;
+ Q cg,cred,gcd,zzz;
+ RHist h;
+ double hmag,gmag;
+
+ if ( !g ) {
+ *rp = 0;
+ return 1;
+ }
+ sugar = SG(g);
+ n = NV(g);
+#if 0
+ if ( !mod )
+ hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;
+#else
+ /* XXX */
+ hmag = 0;
+#endif
+ bucket = create_pbucket();
+ add_pbucket(mod,bucket,g);
+ d = 0;
+ mul = (NM)ALLOCA(sizeof(struct oNM)+(nd_wpd-1)*sizeof(unsigned int));
+ while ( 1 ) {
+ hindex = mod?head_pbucket(mod,bucket):head_pbucket_q(bucket);
+ if ( hindex < 0 ) {
+ if ( d )
+ SG(d) = sugar;
+ *rp = d;
+ return 1;
+ }
+ g = bucket->body[hindex];
+ 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(d);
+ free_pbucket(bucket);
+ *rp = 0;
+ 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);
+ mulq_pbucket(bucket,cred);
+ g = bucket->body[hindex];
+ gmag = (double)p_mag((P)HCQ(g));
+ }
+ red = ndv_mul_nm(mod,p,mul);
+ bucket->body[hindex] = nd_remove_head(g);
+ red = nd_remove_head(red);
+ add_pbucket(mod,bucket,red);
+ sugar = MAX(sugar,SG(p)+TD(mul));
+ if ( !mod && hmag && (gmag > hmag) ) {
+ g = normalize_pbucket(mod,bucket);
+ if ( !g ) {
+ if ( d )
+ SG(d) = sugar;
+ *rp = d;
+ return 1;
+ }
+ nd_removecont2(d,g);
+ hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;
+ add_pbucket(mod,bucket,g);
+ }
+ } else if ( !full ) {
+ g = normalize_pbucket(mod,bucket);
+ if ( g )
+ SG(g) = sugar;
+ *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;
+ }
+ bucket->body[hindex] = g;
+ NEXT(m) = 0;
+ if ( d ) {
+ NEXT(tail)=m;
+ tail=m;
+ LEN(d)++;
+ } else {
+ MKND(n,m,1,d);
+ tail = BDY(d);
+ }
+ }
+ }
+}
+
+/* input : list of DP, cand : list of DP */
+
+int nd_check_candidate(NODE input,NODE cand)
+{
+ int n,i,stat;
+ ND nf,d;
+ NODE t;
+
+ nd_setup(0,cand);
+ n = length(cand);
+
+ /* membercheck : list is a subset of Id(cand) ? */
+ for ( t = input; t; t = NEXT(t) ) {
+ d = dptond(0,(DP)BDY(t));
+ stat = nd_nf_direct(0,d,nd_psq,n,0,&nf);
+ if ( !stat )
+ nd_reconstruct_direct(0,nd_psq,n);
+ else if ( nf )
+ return 0;
+ }
+ /* gbcheck : cand is a GB of Id(cand) ? */
+ if ( !nd_gb(0,1) )
+ return 0;
+ /* XXX */
+ return 1;
+}
+
ND nd_remove_head(ND p)
{
NM m;
@@ -1168,8 +1341,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;
}
@@ -1194,14 +1369,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);
}
@@ -1238,12 +1416,7 @@ int head_pbucket(int mod,PGeoBucket g)
} else {
di = HDL(gi);
nv = NV(gi);
- if ( HTD(gi) > HTD(gj) )
- c = 1;
- else if ( HTD(gi) < HTD(gj) )
- c = -1;
- else
- c = ndl_compare(di,dj);
+ c = NM_COMPARE(BDY(gi),BDY(gj));
if ( c > 0 ) {
if ( sum )
HCM(gj) = sum;
@@ -1292,12 +1465,7 @@ int head_pbucket_q(PGeoBucket g)
} else {
di = HDL(gi);
nv = NV(gi);
- if ( HTD(gi) > HTD(gj) )
- c = 1;
- else if ( HTD(gi) < HTD(gj) )
- c = -1;
- else
- c = ndl_compare(di,dj);
+ c = NM_COMPARE(BDY(gi),BDY(gj));
if ( c > 0 ) {
if ( sum )
HCQ(gj) = sum;
@@ -1330,12 +1498,17 @@ ND normalize_pbucket(int mod,PGeoBucket g)
ND r,t;
r = 0;
- for ( i = 0; i <= g->m; i++ )
+ for ( i = 0; i <= g->m; i++ ) {
r = nd_add(mod,r,g->body[i]);
+ g->body[i] = 0;
+ }
+ g->m = -1;
return r;
}
-NODE nd_gb(int m)
+/* return value = 0 => input is not a GB */
+
+NODE nd_gb(int m,int checkonly)
{
int i,nh,sugar,stat;
NODE r,g,t;
@@ -1362,12 +1535,13 @@ 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);
goto again;
} else if ( nf ) {
+ if ( checkonly ) return 0;
printf("+"); fflush(stdout);
nh = nd_newps(m,nf);
d = update_pairs(d,g,nh);
@@ -1414,7 +1588,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);
@@ -1426,6 +1600,8 @@ again:
if ( nfq ) {
printf("+"); fflush(stdout);
nh = nd_newps_trace(m,nf,nfq);
+ /* failure; m|HC(nfq) */
+ if ( nf < 0 ) return 0;
d = update_pairs(d,g,nh);
g = update_base(g,nh);
} else {
@@ -1445,10 +1621,13 @@ int ndv_compare(NDV *p1,NDV *p2)
{
int td1,td2;
- td1 = HTD(*p1); td2 = HTD(*p2);
- if ( td1 > td2 ) return 1;
- else if ( td1 < td2 ) return -1;
- else return ndl_compare(HDL(*p1),HDL(*p2));
+ if ( nd_degcompat ) {
+ td1 = HTD(*p1); td2 = HTD(*p2);
+ if ( td1 > td2 ) return 1;
+ else if ( td1 < td2 ) return -1;
+ else return ndl_compare(HDL(*p1),HDL(*p2));
+ } else
+ return (*nmv_compare_generic)(BDY(*p1),BDY(*p2));
}
int ndv_compare_rev(NDV *p1,NDV *p2)
@@ -1754,13 +1933,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;
@@ -1826,6 +1998,7 @@ int nd_newps_trace(int mod,ND nf,ND nfq)
nd_bound = (unsigned int **)
REALLOC((char *)nd_bound,nd_pslen*sizeof(unsigned int *));
}
+ if ( !rem(NM(HCQ(nfq)),mod) ) return -1;
nd_removecont(mod,nf);
nd_ps[nd_psn] = ndtondv(mod,nf);
@@ -1958,18 +2131,8 @@ void nd_gr(LIST f,LIST v,int m,struct order_spec *ord,
get_vars((Obj)f,&fv); pltovl(v,&vv);
nd_nvar = length(vv);
- if ( ord->id )
- error("nd_gr : unsupported order");
- switch ( ord->ord.simple ) {
- case 0:
- is_rlex = 1;
- break;
- case 1:
- is_rlex = 0;
- break;
- default:
- error("nd_gr : unsupported order");
- }
+ nd_init_ord(ord);
+ /* XXX for DP */
initd(ord);
for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
ptod(CO,vv,(P)BDY(t),&b);
@@ -1977,7 +2140,7 @@ void nd_gr(LIST f,LIST v,int m,struct order_spec *ord,
}
if ( fd0 ) NEXT(fd) = 0;
nd_setup(m,fd0);
- x = nd_gb(m);
+ x = nd_gb(m,0);
fprintf(asir_out,"found=%d,notfirst=%d,create=%d\n",
nd_found,nd_notfirst,nd_create);
x = nd_reduceall(m,x);
@@ -1999,59 +2162,31 @@ void nd_gr_trace(LIST f,LIST v,int m,int homo,struct o
{
struct order_spec ord1;
VL fv,vv,vc;
- NODE fd,fd0,r,r0,t,x,s,xx,t0;
+ NODE fd,fd0,in0,in,r,r0,t,s,cand;
DP a,b,c,h;
- NDV *w;
- int len,i,j;
+ P p;
get_vars((Obj)f,&fv); pltovl(v,&vv);
nd_nvar = length(vv);
- if ( ord->id )
- error("nd_gr : unsupported order");
initd(ord);
if ( homo ) {
homogenize_order(ord,nd_nvar,&ord1);
- switch ( ord1.ord.simple ) {
- case 0: is_rlex = 1; break;
- case 1: is_rlex = 0; break;
- default: error("nd_gr : unsupported order");
- }
- for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
+ for ( fd0 = 0, in0 = 0, t = BDY(f); t; t = NEXT(t) ) {
ptod(CO,vv,(P)BDY(t),&c);
if ( c ) {
dp_homo(c,&h); NEXTNODE(fd0,fd); BDY(fd) = (pointer)h;
+ NEXTNODE(in0,in); BDY(in) = (pointer)c;
}
}
if ( fd0 ) NEXT(fd) = 0;
+ if ( in0 ) NEXT(in) = 0;
+ nd_init_ord(&ord1);
initd(&ord1);
nd_nvar++;
- nd_setup_trace(m,fd0);
- x = nd_gb_trace(m);
-
- /* dehomogenization */
- for ( t = x; t; t = NEXT(t) )
- ndv_dehomogenize((NDV)BDY(t));
- nd_nvar--;
- nd_setup_parameters();
- initd(ord);
- len = length(x);
- w = (NDV *)ALLOCA(len*sizeof(NDV));
- for ( i = 0, t = x; i < len; i++, t = NEXT(t) ) w[i] = BDY(t);
- for ( i = 0; i < len; i++ ) {
- for ( j = 0; j < i; j++ ) {
- if ( w[i] && w[j] )
- if ( ndl_reducible(HDL(w[i]),HDL(w[j])) ) w[i] = 0;
- else if ( ndl_reducible(HDL(w[j]),HDL(w[i])) ) w[j] = 0;
- }
- }
- for ( i = len-1, t0 = 0; i >= 0; i-- ) {
- if ( w[i] ) { NEXTNODE(t0,t); BDY(t) = (pointer)w[i]; }
- }
- NEXT(t) = 0; x = t0;
} else {
switch ( ord->ord.simple ) {
- case 0: is_rlex = 1; break;
- case 1: is_rlex = 0; break;
+ case 0: nd_isrlex = 1; nd_degcompat = 1; break;
+ case 1: nd_isrlex = 0; nd_degcompat = 1; break;
default: error("nd_gr : unsupported order");
}
for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
@@ -2061,20 +2196,38 @@ void nd_gr_trace(LIST f,LIST v,int m,int homo,struct o
}
}
if ( fd0 ) NEXT(fd) = 0;
- /* setup over GF(m) */
+ in0 = fd0;
+ }
+ do {
nd_setup_trace(m,fd0);
- x = nd_gb_trace(m);
+ cand = nd_gb_trace(m);
+ if ( !cand ) continue;
+ if ( homo ) {
+ /* dehomogenization */
+ for ( t = cand; t; t = NEXT(t) )
+ ndv_dehomogenize((NDV)BDY(t));
+ nd_nvar--;
+ nd_setup_parameters();
+ initd(ord);
+ cand = nd_reducebase(cand);
+ }
+ fprintf(asir_out,"found=%d,notfirst=%d,create=%d\n",
+ nd_found,nd_notfirst,nd_create);
+ cand = nd_reduceall(0,cand);
+ initd(ord);
+ for ( r0 = 0; cand; cand = NEXT(cand) ) {
+ NEXTNODE(r0,r);
+ BDY(r) = (pointer)ndvtodp(0,(NDV)BDY(cand));
+ }
+ if ( r0 ) NEXT(r) = 0;
+ cand = r0;
+ } while ( !nd_check_candidate(in0,cand) );
+ /* dp->p */
+ for ( r = cand; r; r = NEXT(r) ) {
+ dtop(CO,vv,BDY(r),&p);
+ BDY(r) = (pointer)p;
}
- fprintf(asir_out,"found=%d,notfirst=%d,create=%d\n",
- nd_found,nd_notfirst,nd_create);
- x = nd_reduceall(0,x);
- for ( r0 = 0; x; x = NEXT(x) ) {
- NEXTNODE(r0,r);
- a = ndvtodp(0,(NDV)BDY(x));
- dtop(CO,vv,a,(P *)&BDY(r));
- }
- if ( r0 ) NEXT(r) = 0;
- MKLIST(*rp,r0);
+ MKLIST(*rp,cand);
}
void dltondl(int n,DL dl,unsigned int *r)
@@ -2084,7 +2237,7 @@ void dltondl(int n,DL dl,unsigned int *r)
d = dl->d;
bzero(r,nd_wpd*sizeof(unsigned int));
- if ( is_rlex )
+ if ( nd_isrlex )
for ( i = 0; i < n; i++ )
r[(n-1-i)/nd_epw] |= (d[i]<<((nd_epw-((n-1-i)%nd_epw)-1)*nd_bpe));
else
@@ -2101,7 +2254,7 @@ DL ndltodl(int n,int td,unsigned int *ndl)
NEWDL(dl,n);
TD(dl) = td;
d = dl->d;
- if ( is_rlex )
+ if ( nd_isrlex )
for ( i = 0; i < n; i++ )
d[i] = (ndl[(n-1-i)/nd_epw]>>((nd_epw-((n-1-i)%nd_epw)-1)*nd_bpe))
&((1<>((nd_epw-((n-1-i)%nd_epw)-1)*nd_bpe))
@@ -2589,6 +2741,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();
@@ -2604,7 +2767,7 @@ void ndl_dup(int obpe,unsigned int *d,unsigned int *r)
cbpe = nd_bpe;
for ( i = 0; i < nd_wpd; i++ )
r[i] = 0;
- if ( is_rlex )
+ if ( nd_isrlex )
for ( i = 0; i < n; i++ ) {
ei = (d[(n-1-i)/oepw]>>((oepw-((n-1-i)%oepw)-1)*obpe))
&((1<= 0; i-- ) {
+ if ( w[i] ) { NEXTNODE(t0,t); BDY(t) = (pointer)w[i]; }
+ }
+ NEXT(t) = 0; x = t0;
+ return x;
+}
+
+void nd_init_ord(struct order_spec *ord)
+{
+ if ( ord->id )
+ error("nd_gr : unsupported order");
+ nd_degcompat = 0;
+ switch ( ord->ord.simple ) {
+ case 0:
+ nd_degcompat = 1;
+ nd_isrlex = 1;
+ break;
+ case 1:
+ nd_degcompat = 1;
+ nd_isrlex = 0;
+ break;
+ default:
+ error("nd_gr : unsupported order");
+ }
+}
+
|