===================================================================
RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/gr.c,v
retrieving revision 1.50
retrieving revision 1.56
diff -u -p -r1.50 -r1.56
--- OpenXM_contrib2/asir2000/builtin/gr.c	2003/07/01 08:12:38	1.50
+++ OpenXM_contrib2/asir2000/builtin/gr.c	2003/12/26 02:38:10	1.56
@@ -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/gr.c,v 1.49 2003/06/21 02:09:15 noro Exp $
+ * $OpenXM: OpenXM_contrib2/asir2000/builtin/gr.c,v 1.55 2003/12/25 08:46:19 noro Exp $
 */
 #include "ca.h"
 #include "parse.h"
@@ -97,7 +97,7 @@ int NoGCD = 0;
 int GenTrace = 0;
 int OXCheck = -1;
 
-static int NoSugar = 0;
+int NoSugar = 0;
 static int NoCriB = 0;
 static int NoGC = 0;
 static int NoMC = 0;
@@ -105,11 +105,11 @@ static int NoRA = 0;
 static int ShowMag = 0;
 static int Stat = 0;
 static int Denominator = 1;
-static int Top = 0;
-static int Reverse = 0;
+int Top = 0;
+int Reverse = 0;
 static int Max_mag = 0;
 static int Max_coef = 0;
-static char *Demand = 0;
+char *Demand = 0;
 static int PtozpRA = 0;
 
 int doing_f4;
@@ -123,17 +123,70 @@ void Pox_pop_local(NODE,Obj *);
 INLINE int eqdl(int nv,DL dl1,DL dl2)
 {
 	int i;
-	int *b1,*b2;
+	int *p1,*p2;
 
 	if ( dl1->td != dl2->td )
 		return 0;
-	for ( i = 0, b1 = dl1->d, b2 = dl2->d; i < nv; i++ )
-		if ( b1[i] != b2[i] )
-			break;
-	if ( i == nv )
-		return 1;
-	else
-		return 0;
+	i = nv-1;
+	p1 = dl1->d;
+	p2 = dl2->d;
+	while ( i >= 7 ) {
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		if ( *p1++ != *p2++ ) return 0;
+		i -= 8;
+	}
+	switch ( i ) {
+		case 6:
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		case 5:
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		case 4:
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		case 3:
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		case 2:
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		case 1:
+			if ( *p1++ != *p2++ ) return 0;
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		case 0:
+			if ( *p1++ != *p2++ ) return 0;
+			return 1;
+		default:
+			return 1;
+	}
 }
 
 /* b[] should be cleared */
@@ -251,6 +304,7 @@ void dp_gr_main(LIST f,LIST v,Num homo,int modular,int
 	NODE fd,fd0,fi,fi0,r,r0,t,subst,x,s,xx;
 	NODE ind,ind0;
 	LIST trace,gbindex;
+	int input_is_dp = 0;
 
 	mindex = 0; nochk = 0; dp_fcoeffs = field;
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
@@ -263,14 +317,22 @@ void dp_gr_main(LIST f,LIST v,Num homo,int modular,int
 		homogenize_order(ord,NVars,&ord1);
 		for ( fd0 = fi0 = 0, t = BDY(f); t; t = NEXT(t) ) {
 			NEXTNODE(fd0,fd); NEXTNODE(fi0,fi);
-			ptod(CO,vv,(P)BDY(t),(DP *)&BDY(fi)); dp_homo((DP)BDY(fi),(DP *)&BDY(fd));
+			if ( BDY(t) && OID(BDY(t)) == O_DP ) {
+				dp_sort((DP)BDY(t),(DP *)&BDY(fi)); input_is_dp = 1;
+			} else 
+				ptod(CO,vv,(P)BDY(t),(DP *)&BDY(fi));
+			dp_homo((DP)BDY(fi),(DP *)&BDY(fd));
 		}
 		if ( fd0 ) NEXT(fd) = 0;
 		if ( fi0 ) NEXT(fi) = 0;
 		initd(&ord1);
 	} else {
 		for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
-			NEXTNODE(fd0,fd); ptod(CO,vv,(P)BDY(t),(DP *)&BDY(fd));
+			NEXTNODE(fd0,fd);
+			if ( BDY(t) && OID(BDY(t)) == O_DP ) {
+				dp_sort((DP)BDY(t),(DP *)&BDY(fd)); input_is_dp = 1;
+			} else
+				ptod(CO,vv,(P)BDY(t),(DP *)&BDY(fd));
 		}
 		if ( fd0 ) NEXT(fd) = 0;
 		fi0 = fd0;
@@ -318,7 +380,10 @@ void dp_gr_main(LIST f,LIST v,Num homo,int modular,int
 	}
 	for ( r0 = 0, ind0 = 0; x; x = NEXT(x) ) {
 		NEXTNODE(r0,r); dp_load((int)BDY(x),&ps[(int)BDY(x)]);
-		dtop(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
+		if ( input_is_dp )
+			BDY(r) = (pointer)ps[(int)BDY(x)];
+		else
+			dtop(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
 		NEXTNODE(ind0,ind);
 		STOQ((int)BDY(x),q); BDY(ind) = q;
 	}
@@ -348,7 +413,8 @@ void dp_gr_mod_main(LIST f,LIST v,Num homo,int m,struc
 	VL fv,vv,vc;
 	NODE fd,fd0,r,r0,t,x,s,xx;
 	DP a,b,c;
-extern struct oEGT eg_red_mod;
+	extern struct oEGT eg_red_mod;
+	int input_is_dp = 0;
 
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
 	NVars = length((NODE)vv); PCoeffs = vc ? 1 : 0; VC = vc;
@@ -358,7 +424,11 @@ extern struct oEGT eg_red_mod;
 	initd(ord);
 	if ( homo ) {
 		for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
-			ptod(CO,vv,(P)BDY(t),&a); dp_homo(a,&b);
+			if ( BDY(t) && OID(BDY(t)) == O_DP ) {
+				dp_sort((DP)BDY(t),&a); input_is_dp = 1;
+			} else 
+				ptod(CO,vv,(P)BDY(t),&a);
+			dp_homo(a,&b);
 			if ( PCoeffs )
 				dp_mod(b,m,0,&c);
 			else
@@ -370,7 +440,10 @@ extern struct oEGT eg_red_mod;
 		homogenize_order(ord,NVars,&ord1); initd(&ord1);
 	} else {
 		for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {
-			ptod(CO,vv,(P)BDY(t),&b);
+			if ( BDY(t) && OID(BDY(t)) == O_DP ) {
+				dp_sort((DP)BDY(t),&b); input_is_dp = 1;
+			} else 
+				ptod(CO,vv,(P)BDY(t),&b);
 			if ( PCoeffs )
 				dp_mod(b,m,0,&c);
 			else
@@ -396,11 +469,19 @@ extern struct oEGT eg_red_mod;
 	reduceall_mod(x,m,&xx); x = xx;
 	if ( PCoeffs )
 		for ( r0 = 0; x; x = NEXT(x) ) {
-			NEXTNODE(r0,r); mdtop(CO,m,vv,ps[(int)BDY(x)],(P *)&BDY(r));
+			NEXTNODE(r0,r); 
+			if ( input_is_dp )
+				mdtodp(ps[(int)BDY(x)],(DP *)&BDY(r));
+			else
+				mdtop(CO,m,vv,ps[(int)BDY(x)],(P *)&BDY(r));
 		}
 	else
 		for ( r0 = 0; x; x = NEXT(x) ) {
-			NEXTNODE(r0,r); _dtop_mod(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
+			NEXTNODE(r0,r); 
+			if ( input_is_dp )
+				_mdtodp(ps[(int)BDY(x)],(DP *)&BDY(r));
+			else
+				_dtop_mod(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
 		}
 	print_stat();
 	if ( r0 ) NEXT(r) = 0;
@@ -412,6 +493,7 @@ void dp_f4_main(LIST f,LIST v,struct order_spec *ord,L
 	int homogen;
 	VL fv,vv,vc;
 	NODE fd,fd0,r,r0,t,x,s,xx;
+	int input_is_dp = 0;
 
 	dp_fcoeffs = 0;
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
@@ -421,7 +503,11 @@ void dp_f4_main(LIST f,LIST v,struct order_spec *ord,L
 		error("dp_f4_main : invalid order specification");
 	initd(ord);
 	for ( fd0 = 0, t = BDY(f), homogen = 1; t; t = NEXT(t) ) {
-		NEXTNODE(fd0,fd); ptod(CO,vv,(P)BDY(t),(DP *)&BDY(fd));
+		NEXTNODE(fd0,fd); 
+		if ( BDY(t) && OID(BDY(t)) == O_DP ) {
+			dp_sort((DP)BDY(t),(DP *)&BDY(fd)); input_is_dp = 1;
+		} else 
+			ptod(CO,vv,(P)BDY(t),(DP *)&BDY(fd));
 		if ( homogen )
 			homogen = dp_homogeneous(BDY(fd));
 	}
@@ -433,7 +519,10 @@ void dp_f4_main(LIST f,LIST v,struct order_spec *ord,L
 	}
 	for ( r0 = 0; x; x = NEXT(x) ) {
 		NEXTNODE(r0,r); dp_load((int)BDY(x),&ps[(int)BDY(x)]);
-		dtop(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
+		if ( input_is_dp )
+			BDY(r) = (pointer)ps[(int)BDY(x)];
+		else
+			dtop(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
 	}
 	if ( r0 ) NEXT(r) = 0;
 	MKLIST(*rp,r0);
@@ -445,6 +534,7 @@ void dp_f4_mod_main(LIST f,LIST v,int m,struct order_s
 	VL fv,vv,vc;
 	DP b,c,c1;
 	NODE fd,fd0,r,r0,t,x,s,xx;
+	int input_is_dp = 0;
 
 	dp_fcoeffs = 0;
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
@@ -454,7 +544,10 @@ void dp_f4_mod_main(LIST f,LIST v,int m,struct order_s
 		error("dp_f4_mod_main : invalid order specification");
 	initd(ord);
 	for ( fd0 = 0, t = BDY(f), homogen = 1; t; t = NEXT(t) ) {
-		ptod(CO,vv,(P)BDY(t),&b);
+		if ( BDY(t) && OID(BDY(t)) == O_DP ) {
+			dp_sort((DP)BDY(t),&b); input_is_dp = 1;
+		} else 
+			ptod(CO,vv,(P)BDY(t),&b);
 		if ( homogen )
 			homogen = dp_homogeneous(b);
 		_dp_mod(b,m,0,&c);
@@ -474,7 +567,11 @@ void dp_f4_mod_main(LIST f,LIST v,int m,struct order_s
 		reduceall_mod(x,m,&xx); x = xx;
 	}
 	for ( r0 = 0; x; x = NEXT(x) ) {
-		NEXTNODE(r0,r); _dtop_mod(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
+		NEXTNODE(r0,r);
+		if ( input_is_dp )
+			_mdtodp(ps[(int)BDY(x)],(DP *)&BDY(r));
+		else
+			_dtop_mod(CO,vv,ps[(int)BDY(x)],(P *)&BDY(r));
 	}
 	if ( r0 ) NEXT(r) = 0;
 	MKLIST(*rp,r0);
@@ -601,6 +698,8 @@ NODE gb_f4(NODE f)
 /* initial bases are monic */
 
 unsigned int **psca;
+GeoBucket create_bucket();
+DL remove_head_bucket(GeoBucket,int);
 
 NODE gb_f4_mod(NODE f,int m)
 {
@@ -619,6 +718,8 @@ NODE gb_f4_mod(NODE f,int m)
 	int *indred,*isred;
 	CDP ri;
 	int pscalen;
+	GeoBucket bucket;
+	DL head;
 	struct oEGT tmp0,tmp1,eg_split_symb,eg_split_conv,eg_split_elim1,eg_split_elim2;
 	extern struct oEGT eg_symb,eg_conv,eg_elim1,eg_elim2;
 
@@ -641,37 +742,49 @@ NODE gb_f4_mod(NODE f,int m)
 		minsugar(d,&dm,&dr); d = dr;
 		if ( DP_Print )
 			fprintf(asir_out,"sugar=%d\n",dm->sugar);
-		blist = 0; s0 = 0;
+		blist = 0;
+		bucket = create_bucket();
 		/* asph : sum of all head terms of spoly */
 		for ( t = dm; t; t = NEXT(t) ) {
 			_dp_sp_mod(ps[t->dp1],ps[t->dp2],m,&sp);
 /*			fprintf(stderr,"splen=%d-",dp_nt(sp)); */
 			if ( sp ) {
 				MKNODE(bt,sp,blist); blist = bt;
-				s0 = symb_merge(s0,dp_dllist(sp),nv);
+				add_bucket(bucket,dp_dllist(sp),nv);
 /*				fprintf(stderr,"%d-",length(s0)); */
 			}
 		}
+#if 0
 		if ( DP_Print )
 			fprintf(asir_out,"initial spmat : %d x %d ",length(blist),length(s0));
+#endif
 		/* s0 : all the terms appeared in symbolic reduction */
-		for ( s = s0, nred = 0; s; s = NEXT(s) ) {
+		nred = 0;
+		s0 = 0;
+		while ( 1 ) {
+			head = remove_head_bucket(bucket,nv);
+			if ( !head ) break;
+			else {
+				NEXTNODE(s0,s);
+				BDY(s) = (pointer)head;
+			}
 			for ( r = gall;	r; r = NEXT(r) )
-				if ( _dl_redble(BDY(ps[(int)BDY(r)])->dl,BDY(s),nv) )
+				if ( _dl_redble(BDY(ps[(int)BDY(r)])->dl,head,nv) )
 					break;
 			if ( r ) {
-				dltod(BDY(s),nv,&tdp);
+				dltod(head,nv,&tdp);
 				dp_subd(tdp,ps[(int)BDY(r)],&sd);
 				dt = mul_dllist(BDY(sd)->dl,ps[(int)BDY(r)]);
+				add_bucket(bucket,NEXT(dt),nv);
 /*				fprintf(stderr,"[%d]",length(dt)); */
 				/* list of [t,f] */
 				bt1 = mknode(2,BDY(sd)->dl,BDY(r));
 				MKNODE(bt,bt1,blist); blist = bt;
-				symb_merge(s,dt,nv);
 /*				fprintf(stderr,"%d-",length(s0));  */
 				nred++;
 			}
 		}
+		if ( s0 ) NEXT(s) = 0;
 /*		fprintf(stderr,"\n"); */
 		get_eg(&tmp1); add_eg(&eg_symb,&tmp0,&tmp1);
 		init_eg(&eg_split_symb); add_eg(&eg_split_symb,&tmp0,&tmp1);
@@ -1905,7 +2018,7 @@ DP_pairs criterion_F( DP_pairs d1 )
 	for ( head = last = 0, p = d1; NEXT(p); ) {
 		s = (r = w = collect_pairs_of_hdlcm( p, &rest ))->sugar;
 		while ( w = NEXT(w) )
-			if ( criterion_2( w->dp1, w->dp2 ) ) {
+			if ( !do_weyl && criterion_2( w->dp1, w->dp2 ) ) {
 				r = w;
 				break;
 			} else if ( w->sugar < s ) s = (r = w)->sugar;