===================================================================
RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/gr.c,v
retrieving revision 1.14
retrieving revision 1.21
diff -u -p -r1.14 -r1.21
--- OpenXM_contrib2/asir2000/builtin/gr.c	2000/12/08 02:39:05	1.14
+++ OpenXM_contrib2/asir2000/builtin/gr.c	2001/09/04 08:48:19	1.21
@@ -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.13 2000/12/05 08:29:43 noro Exp $ 
+ * $OpenXM: OpenXM_contrib2/asir2000/builtin/gr.c,v 1.20 2001/07/23 05:03:22 noro Exp $ 
 */
 #include "ca.h"
 #include "parse.h"
@@ -105,7 +105,8 @@ DP_pairs criterion_M(DP_pairs);
 DP_pairs criterion_B(DP_pairs,int);
 DP_pairs newpairs(NODE,int);
 DP_pairs updpairs(DP_pairs,NODE,int);
-void _dp_nf(NODE,DP,DP *,int,int,DP *);
+void _dp_nf(NODE,DP,DP *,int,DP *);
+void _dp_nf_z(NODE,DP,DP *,int,int,DP *);
 NODE gb_mod(NODE,int);
 NODE gbd(NODE,int,NODE,NODE);
 NODE gb(NODE,int,NODE);
@@ -134,8 +135,8 @@ void pltovl(LIST,VL *);
 void printdl(DL);
 int DPPlength(DP_pairs);
 void dp_gr_mod_main(LIST,LIST,Num,int,struct order_spec *,LIST *);
-void dp_gr_main(LIST,LIST,Num,int,struct order_spec *,LIST *);
-void dp_f4_main(LIST,LIST,struct order_spec *,LIST *);
+void dp_gr_main(LIST,LIST,Num,int,int,struct order_spec *,LIST *);
+void dp_f4_main(LIST,LIST,int,struct order_spec *,LIST *);
 void dp_f4_mod_main(LIST,LIST,int,struct order_spec *,LIST *);
 double get_rtime();
 void _dpmod_to_vect(DP,DL *,int *);
@@ -150,17 +151,18 @@ static P *psc;
 
 static int *pss;
 static int psn,pslen;
-static int NVars,CNVars,PCoeffs;
+static int NVars,CNVars;
 static VL VC;
 
+int PCoeffs;
 int DP_Print = 0;
 int DP_Multiple = 0;
+int DP_NFStat = 0;
 LIST Dist = 0;
 int NoGCD = 0;
 int GenTrace = 0;
 int OXCheck = -1;
 
-static DP_NFStat = 0;
 static int NoSugar = 0;
 static int NoCriB = 0;
 static int NoGC = 0;
@@ -255,10 +257,10 @@ NODE f;
 	printf("\n");
 }
 
-void dp_gr_main(f,v,homo,modular,ord,rp)
+void dp_gr_main(f,v,homo,modular,field,ord,rp)
 LIST f,v;
 Num homo;
-int modular;
+int modular,field;
 struct order_spec *ord;
 LIST *rp;
 {
@@ -267,7 +269,7 @@ LIST *rp;
 	VL fv,vv,vc;
 	NODE fd,fd0,fi,fi0,r,r0,t,subst,x,s,xx;
 
-	mindex = 0; nochk = 0; dp_fcoeffs = 0;
+	mindex = 0; nochk = 0; dp_fcoeffs = field;
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
 	NVars = length((NODE)vv); PCoeffs = vc ? 1 : 0; VC = vc;
 	CNVars = homo ? NVars+1 : NVars;
@@ -294,7 +296,7 @@ LIST *rp;
 		modular = -modular; nochk = 1;
 	}
 	if ( modular )
-		m = modular > 1 ? modular : lprime[mindex];
+		m = modular > 1 ? modular : get_lprime(mindex);
 	else
 		m = 0;
 	makesubst(vc,&subst);
@@ -324,7 +326,7 @@ LIST *rp;
 			if ( modular > 1 ) {
 				*rp = 0; return;
 			} else
-				m = lprime[++mindex];
+				m = get_lprime(++mindex);
 		makesubst(vc,&subst);
 		psn = length(s);
 		for ( i = psn; i < pslen; i++ ) {
@@ -353,6 +355,7 @@ LIST *rp;
 	VL fv,vv,vc;
 	NODE fd,fd0,r,r0,t,x,s,xx;
 	DP a,b,c;
+extern struct oEGT eg_red_mod;
 
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
 	NVars = length((NODE)vv); PCoeffs = vc ? 1 : 0; VC = vc;
@@ -390,7 +393,9 @@ LIST *rp;
 	if ( homo ) {
 		initd(&ord1); CNVars = NVars+1;
 	}
+/* init_eg(&eg_red_mod); */
 	x = gb_mod(s,m);
+/* print_eg("Red_mod",&eg_red_mod); */
 	if ( homo ) {
 		reducebase_dehomo(x,&xx); x = xx;
 		initd(ord); CNVars = NVars;
@@ -409,8 +414,9 @@ LIST *rp;
 	MKLIST(*rp,r0);
 }
 
-void dp_f4_main(f,v,ord,rp)
+void dp_f4_main(f,v,field,ord,rp)
 LIST f,v;
+int field;
 struct order_spec *ord;
 LIST *rp;
 {
@@ -419,7 +425,7 @@ LIST *rp;
 	VL fv,vv,vc;
 	NODE fd,fd0,fi,fi0,r,r0,t,subst,x,s,xx;
 
-	dp_fcoeffs = 0;
+	dp_fcoeffs = field;
 	get_vars((Obj)f,&fv); pltovl(v,&vv); vlminus(fv,vv,&vc);
 	NVars = length((NODE)vv); PCoeffs = vc ? 1 : 0; VC = vc;
 	CNVars = NVars;
@@ -726,11 +732,13 @@ int m;
 			for ( j = 0; j < spcol; j++ )
 				if ( spmat[i][j] )
 					nonzero++;
-		if ( DP_Print )
+		if ( DP_Print && nsp )
 			fprintf(asir_out,"spmat : %d x %d (nonzero=%f%%)...",
 				nsp,spcol,((double)nonzero*100)/(nsp*spcol));
-		rank = generic_gauss_elim_mod(spmat,nsp,spcol,m,colstat);
-
+		if ( nsp )
+			rank = generic_gauss_elim_mod(spmat,nsp,spcol,m,colstat);
+		else
+			rank = 0;
 		get_eg(&tmp1); add_eg(&eg_elim2,&tmp0,&tmp1);
 		init_eg(&eg_split_elim2); add_eg(&eg_split_elim2,&tmp0,&tmp1);
 
@@ -895,7 +903,11 @@ int m;
 	int i;
 	NODE s,s0,f0;
 
+#if 1
 	f0 = f = NODE_sortb(f,1);
+#else
+	f0 = f;
+#endif
 	psn = length(f); pslen = 2*psn;
 	ps = (DP *)MALLOC(pslen*sizeof(DP));
 	psh = (DL *)MALLOC(pslen*sizeof(DL));
@@ -1060,7 +1072,7 @@ NODE *h;
 			MKLIST(hist,node);
 			MKNODE(TraceList,hist,0);
 		}
-		_dp_nf(top,ps[w[i]],ps,1,PtozpRA?DP_Multiple:0,&g);
+		_dp_nf(top,ps[w[i]],ps,1,&g);
 		prim_part(g,0,&g1);
 		get_eg(&tmp1); add_eg(&eg_ra,&tmp0,&tmp1);
 		if ( DP_Print || DP_PrintShort ) {
@@ -1405,7 +1417,10 @@ NODE subst;
 				new_sugar = h->sugar;
 			get_eg(&tnf0);
 			t_0 = get_rtime();
-			_dp_nf(gall,h,ps,!Top,DP_Multiple,&nf);
+			if ( PCoeffs || dp_fcoeffs )
+				_dp_nf(gall,h,ps,!Top,&nf);
+			else
+				_dp_nf_z(gall,h,ps,!Top,DP_Multiple,&nf);
 			if ( DP_Print )
 				fprintf(asir_out,"(%.3g)",get_rtime()-t_0);
 			get_eg(&tnf1); add_eg(&eg_nf,&tnf0,&tnf1);
@@ -1803,8 +1818,10 @@ NODE f;
 	while ( d ) {
 		l = d; d = NEXT(d);
 		get_eg(&tmp0);
-		dp_load(l->dp1,&dp1); dp_load(l->dp2,&dp2); dp_sp(dp1,dp2,&h);
-		_dp_nf(gall,h,ps,1,0,&nf);
+		dp_load(l->dp1,&dp1); dp_load(l->dp2,&dp2);
+		dp_sp(dp1,dp2,&h);
+/* fprintf(stderr,"{%d,%d}",l->dp1,l->dp2); */
+		_dp_nf(gall,h,ps,1,&nf);
 		get_eg(&tmp1); add_eg(&eg_gc,&tmp0,&tmp1);
 		if ( DP_Print || DP_PrintShort ) {
 			fprintf(asir_out,"."); fflush(asir_out);
@@ -1830,7 +1847,7 @@ NODE f,x;
 	}
 	for ( ; f; f = NEXT(f) ) {
 		get_eg(&tmp0);
-		_dp_nf(x,(DP)BDY(f),ps,1,0,&g);
+		_dp_nf(x,(DP)BDY(f),ps,1,&g);
 		get_eg(&tmp1); add_eg(&eg_mc,&tmp0,&tmp1);
 		if ( DP_Print ) {
 			print_split_eg(&tmp0,&tmp1); fflush(asir_out);
@@ -1914,6 +1931,7 @@ LIST *list;
 	STOQ(Reverse,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"Reverse"); MKNODE(n1,name,n); n = n1;
 	STOQ(Stat,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"Stat"); MKNODE(n1,name,n); n = n1;
 	STOQ(DP_Print,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"Print"); MKNODE(n1,name,n); n = n1;
+	STOQ(DP_PrintShort,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"PrintShort"); MKNODE(n1,name,n); n = n1;
 	STOQ(DP_NFStat,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"NFStat"); MKNODE(n1,name,n); n = n1;
 	STOQ(OXCheck,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"OXCheck"); MKNODE(n1,name,n); n = n1;
 	STOQ(GenTrace,v); MKNODE(n1,v,n); n = n1; MKSTR(name,"GenTrace"); MKNODE(n1,name,n); n = n1;
@@ -2006,17 +2024,95 @@ void print_stat() {
  *
  */
 
-#define SAFENM(q) ((q)?NM(q):0)
-
 double pz_t_e, pz_t_d, pz_t_d1, pz_t_c, im_t_s, im_t_r;
 
 extern int GenTrace;
 extern NODE TraceList;
+extern int mpi_mag;
 
-void _dp_nf(b,g,ps,full,multiple,r)
+void dp_mulc_d(p,c,r)
+DP p;
+P c;
+DP *r;
+{
+	if ( Dist && BDY(Dist)
+		&& HMAG(p) > mpi_mag
+		&& p_mag((P)c) > mpi_mag ) {
+		if ( DP_NFStat ) fprintf(asir_out,"~");
+		dp_imul_d(p,(Q)c,r);
+	} else {
+		if ( DP_NFStat ) fprintf(asir_out,"_");
+		muldc(CO,p,c,r);
+	}
+}
+
+void _dp_nf(b,g,ps,full,rp)
 NODE b;
 DP g;
 DP *ps;
+int full;
+DP *rp;
+{
+	DP u,p,d,s,t,mult;
+	P coef;
+	NODE l;
+	MP m,mr;
+	int sugar,psugar;
+
+	if ( !g ) {
+		*rp = 0; return;
+	}
+	sugar = g->sugar;
+	for ( d = 0; g; ) {
+		for ( u = 0, l = b; l; l = NEXT(l) ) {
+			if ( dl_redble(BDY(g)->dl,psh[(int)BDY(l)]) ) {
+				dp_load((int)BDY(l),&p);
+				/* t+u = coef*(d+g) - mult*p (t = coef*d) */
+				dp_red(d,g,p,&t,&u,&coef,&mult);
+				psugar = (BDY(g)->dl->td - BDY(p)->dl->td) + p->sugar;
+				sugar = MAX(sugar,psugar);
+				if ( GenTrace ) {
+					LIST hist;
+					Q cq;
+					NODE node,node0;
+
+					STOQ((int)BDY(l),cq);
+					node0 = mknode(4,coef,cq,mult,ONE);
+					MKLIST(hist,node0);
+					MKNODE(node,hist,TraceList); TraceList = node;
+				}
+				if ( !u ) {
+					if ( d )
+						d->sugar = sugar;
+					*rp = d; return;
+				}
+				d = t;
+				break;
+			}
+		}
+		if ( u )
+			g = u;
+		else if ( !full ) {
+			if ( g ) {
+				MKDP(g->nv,BDY(g),t); t->sugar = sugar; g = t;
+			}
+			*rp = g; return;
+		} 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;
+			addd(CO,d,t,&s); d = s;
+			dp_rest(g,&t); g = t;
+		}
+	}
+	if ( d )
+		d->sugar = sugar;
+	*rp = d;
+}
+
+void _dp_nf_z(b,g,ps,full,multiple,r)
+NODE b;
+DP g;
+DP *ps;
 int full,multiple;
 DP *r;
 {
@@ -2029,11 +2125,9 @@ DP *r;
 	int sugar,psugar;
 	NODE dist;
 	STRING imul;
-	int ndist;
 	int kara_bit;
-	extern int mpi_mag;
 	double get_rtime();	
-	double t_0,t_00,tt,ttt,t_p,t_m,t_m1,t_m2,t_s,t_g,t_a;
+	double t_0,t_00,tt,ttt,t_p,t_m,t_g,t_a;
 	LIST hist;
 	NODE node;
 	Q rcred,mrcred;
@@ -2042,14 +2136,10 @@ DP *r;
 		*r = 0; return;
 	}
 	pz_t_e = pz_t_d = pz_t_d1 = pz_t_c = 0;
-	t_p = t_m = t_m1 = t_m2 = t_s = t_g = t_a = 0;
+	t_p = t_m = t_g = t_a = 0;
 
 	denom = Denominator?Denominator:1;
 	hmag = multiple*HMAG(g)/denom;
-	if ( Dist ) {
-		dist = BDY(Dist);
-		ndist = length(dist);
-	}
 	sugar = g->sugar;
 
 	dc = 0; dp = 0; rc = ONE; rp = g;
@@ -2062,49 +2152,20 @@ DP *r;
 				t_0 = get_rtime();
 				dp_load((int)BDY(l),&red);
 				hr = (Q)BDY(rp)->c; hred = (Q)BDY(red)->c;
-				/*
-				 * hr = HC(rp), hred = HC(red)
-				 * cred = hr/GCD(hr,hred), cr = hred/GCD(hr,hred)
-				 */
 				igcd_cofactor((Q)BDY(rp)->c,(Q)BDY(red)->c,&gcd,&cred,&cr);
-				tt = get_rtime()-t_0; t_p += tt;
+				tt = get_rtime(); t_p += tt-t_0;
 
-				t_0 = get_rtime();
 				dp_subd(rp,red,&shift);
-
-				t_00 = get_rtime();
-				if ( Dist && ndist
-					&& HMAG(red) > mpi_mag
-					&& p_mag((P)cr) > mpi_mag ) {
-					if ( DP_NFStat ) fprintf(asir_out,"~");
-					dp_imul_d(rp,cr,&t);
-				} else {
-					if ( DP_NFStat ) fprintf(asir_out,"_");
-					muldc(CO,rp,(P)cr,&t);
-				}
-				ttt = get_rtime()-t_00; t_m1 += ttt/dp_nt(rp);
-
-				t_00 = get_rtime();
+				dp_mulc_d(rp,cr,&t);
 				chsgnp((P)cred,(P *)&mcred);
-				if ( Dist && ndist
-					&& HMAG(red) > mpi_mag
-					&& p_mag((P)mcred) > mpi_mag ) {
-					if ( DP_NFStat ) fprintf(asir_out,"=");
-					dp_imul_d(red,mcred,&t1);
-				} else {
-					if ( DP_NFStat ) fprintf(asir_out,"_");
-					muldc(CO,red,(P)mcred,&t1);
-				}
-				ttt = get_rtime()-t_00; t_m2 += ttt/dp_nt(red);
+				dp_mulc_d(red,mcred,&t1);
+				muld(CO,shift,t1,&t1);
+				addd(CO,t,t1,&u);
+				t_m += get_rtime()-tt;
 
-				t_00 = get_rtime();
-				muld(CO,shift,t1,&t2);
-				addd(CO,t,t2,&u);
-				tt = get_rtime(); t_m += tt-t_0;
-				ttt = get_rtime(); t_s += ttt-t_00;
-
 				psugar = (BDY(rp)->dl->td - BDY(red)->dl->td) + red->sugar;
 				sugar = MAX(sugar,psugar);
+
 				if ( GenTrace ) {
 					/* u = cr*rp + (-cred)*shift*red */ 
 					STOQ((int)BDY(l),cq);
@@ -2114,6 +2175,7 @@ DP *r;
 					muldc(CO,shift,(P)mrcred,(DP *)&ARG2(node));
 					MKLIST(hist,node);
 				}
+
 				if ( !u ) {
 					if ( dp )
 						dp->sugar = sugar;
@@ -2130,37 +2192,29 @@ DP *r;
 		if ( u ) {
 			if ( multiple && HMAG(u) > hmag ) {
 				t_0 = get_rtime();
-				if ( Dist && HMAG(u) > mpi_mag ) {
-					if ( DP_NFStat )
-						fprintf(asir_out,"D");
-					dp_ptozp_d(dist,ndist,u,&t);
-				} else {
-					if ( DP_NFStat )
-						fprintf(asir_out,"L");
-					dp_ptozp_d(0,0,u,&t);
-				}
-				tt = get_rtime()-t_0; t_g += tt;
-				t_0 = get_rtime();
-				divsn(NM((Q)BDY(u)->c),NM((Q)BDY(t)->c),&cn); NTOQ(cn,1,cont);
+				dp_ptozp_d(u,&rp);
+				tt = get_rtime(); t_g += tt-t_0;
+
+				divsq((Q)BDY(u)->c,(Q)BDY(rp)->c,&cont);
 				if ( !dp_fcoeffs && DP_NFStat ) {
-					fprintf(asir_out,"(%d)",p_mag((P)cont)*100/p_mag((P)BDY(u)->c));
+					fprintf(asir_out,
+						"(%d)",p_mag((P)cont)*100/p_mag((P)BDY(u)->c));
 					fflush(asir_out);
 				}
-				mulq(cr,dc,&dcq);
-				mulq(cont,rc,&rcq);
+				mulq(cr,dc,&dcq); mulq(cont,rc,&rcq);
 				igcd_cofactor(dcq,rcq,&gcd,&dc,&rc);
-				tt = get_rtime()-t_0; t_a += tt;
-				rp = t;
+				t_a = get_rtime()-tt;
+
 				hmag = multiple*HMAG(rp)/denom;
 				if ( GenTrace ) {
 					ARG3(BDY(hist)) = (pointer)gcd;
 					MKNODE(node,hist,TraceList); TraceList = node;
 				}
 			} else {
-				t_0 = get_rtime();
-				mulq(cr,dc,&dcq); dc = dcq;
-				tt = get_rtime()-t_0; t_a += tt;
 				rp = u;
+				t_0 = get_rtime();
+				mulq(cr,dc,&dc);
+				t_a += get_rtime()-t_0;
 				if ( GenTrace ) {
 					ARG3(BDY(hist)) = (pointer)ONE;
 					MKNODE(node,hist,TraceList); TraceList = node;
@@ -2175,15 +2229,11 @@ DP *r;
 		} else {
 			t_0 = get_rtime();
 			mulq((Q)BDY(rp)->c,rc,&c);
-			igcd_cofactor(dc,c,&gcd,&dcq,&cq);
-			muldc(CO,dp,(P)dcq,&t1);
-			m = BDY(rp); NEWMP(mr); mr->dl = m->dl; mr->c = (P)cq;
-			NEXT(mr) = 0; MKDP(rp->nv,mr,t); t->sugar = mr->dl->td;
-			addd(CO,t,t1,&dp);
-			dc = gcd;
-
-			dp_rest(rp,&t); rp = t;
-			tt = get_rtime()-t_0; t_a += tt;
+			igcd_cofactor(dc,c,&dc,&dcq,&cq);
+			muldc(CO,dp,(P)dcq,&t);
+			dp_hm(rp,&t1); BDY(t1)->c = (P)cq;  addd(CO,t,t1,&dp);
+			dp_rest(rp,&rp);
+			t_a += get_rtime()-t_0;
 		}
 	}
 	if ( GenTrace ) {
@@ -2195,9 +2245,9 @@ DP *r;
 	*r = dp;
 final:
 	if ( DP_NFStat )
-		fprintf(asir_out,"(%.3g %.3g %.3g %.3g %.3g %.3g %.3g %.3g %.3g %.3g %.3g)",
-			t_p,t_m,t_m1,t_m2,t_s,
-			t_g,t_a,
+		fprintf(asir_out,
+			"(%.3g %.3g %.3g %.3g %.3g %.3g %.3g %.3g)",
+			t_p,t_m,t_g,t_a,
 			pz_t_e, pz_t_d, pz_t_d1, pz_t_c);
 }
 
@@ -2218,7 +2268,6 @@ DP *rp;
 	NODE tn,dist,n0,n1,n2;
 	Obj dmy;
 	STRING imul;
-
 	extern LIST Dist;
 
 	if ( !p || !q ) {