=================================================================== RCS file: /home/cvs/OpenXM/src/kan96xx/Kan/ecart.c,v retrieving revision 1.11 retrieving revision 1.16 diff -u -p -r1.11 -r1.16 --- OpenXM/src/kan96xx/Kan/ecart.c 2003/08/21 02:30:23 1.11 +++ OpenXM/src/kan96xx/Kan/ecart.c 2003/08/27 03:11:12 1.16 @@ -1,4 +1,4 @@ -/* $OpenXM: OpenXM/src/kan96xx/Kan/ecart.c,v 1.10 2003/08/20 05:18:35 takayama Exp $ */ +/* $OpenXM: OpenXM/src/kan96xx/Kan/ecart.c,v 1.15 2003/08/26 05:52:43 takayama Exp $ */ #include #include "datatype.h" #include "extern2.h" @@ -36,6 +36,8 @@ static POLY reduction_ecart1(POLY r,struct gradedPolyS static POLY reduction_ecart1_mod(POLY r,struct gradedPolySet *gset); static POLY ecartCheckSyz0(POLY cf,POLY r_0,POLY syz, struct gradedPolySet *gg,POLY r); +static void ecartCheckSyz0_printinfo(POLY cf,POLY r_0,POLY syz, + struct gradedPolySet *gg,POLY r); extern int DebugReductionRed; extern int TraceLift; @@ -67,8 +69,20 @@ static POLY ecartDivideSv(POLY r,int *d) { } r = r->next; } + if (k > 0) { *d = k; + return ppMult(cxx(1,0,-k,f->m->ringp),f); + }else{ + return f; + } + + /* Do not do the below. It caused a bug. cf. misc-2003/07/ecart/b4.sm1 test2. + Note. 2003.8.26 + */ + /* + if (k > 0) { + *d = k; f = r; while (r != POLYNULL) { r->m->e[0].x -= k; @@ -76,6 +90,7 @@ static POLY ecartDivideSv(POLY r,int *d) { } } return f; + */ } static int ecartGetEll(POLY f,POLY g) { @@ -240,7 +255,36 @@ static POLY ecartCheckSyz0(POLY cf,POLY r_0,POLY syz, return f; } +static void ecartCheckSyz0_printinfo(POLY cf,POLY r_0,POLY syz, + struct gradedPolySet *gg,POLY r) +{ + POLY f; + int grd,i; + POLY q; + struct coeff *c; + fprintf(stderr,"cf=%s\n",POLYToString(cf,'*',1)); + fprintf(stderr,"r_0=%s\n",POLYToString(r_0,'*',1)); + fprintf(stderr,"syz=%s\n",POLYToString(syz,'*',1)); + fprintf(stderr,"r=%s\n",POLYToString(r,'*',1)); + f = ppMult(cf,r_0); + while (syz != POLYNULL) { + grd = syz->m->e[0].x; + i = syz->m->e[0].D; + c = syz->coeffp; + if (c->tag == POLY_COEFF) { + q = c->val.f; + }else{ + q = POLYNULL; + } + fprintf(stderr,"[grd,idx]=[%d,%d], %s\n",grd,i, + POLYToString(((gg->polys[grd])->g)[i],'*',1)); + /* f = ppAdd(f,ppMult(q,((gg->polys[grd])->g)[i])); */ + syz = syz->next; + } + /* f = ppSub(f,r); */ +} + POLY reduction_ecart(r,gset,needSyz,syzp) POLY r; struct gradedPolySet *gset; @@ -248,13 +292,17 @@ POLY reduction_ecart(r,gset,needSyz,syzp) struct syz0 *syzp; /* set */ { POLY rn; + if (TraceLift && needSyz) { + warningGradedSet("TraceLift cannot be used to get syzygy. TraceLift is turned off.\n"); + TraceLift = 0; + } if (TraceLift) { if (EcartAutomaticHomogenization) { if (TraceLift_ringmod == NULL) { warningPoly("reduction_ecart: TraceLift_ringmod is not set.\n"); return reduction_ecart1(r,gset,needSyz,syzp); } - rn = reduction_ecart1_mod(r,gset); /* BUG: syzygy is not obtained. */ + rn = reduction_ecart1_mod(r,gset); if (rn == POLYNULL) return rn; else return reduction_ecart1(r,gset,needSyz,syzp); }else{ @@ -293,6 +341,8 @@ static POLY reduction_ecart0(r,gset,needSyz,syzp) POLY cf_o; POLY syz_o; POLY r_0; + int se; + struct coeff *cont; extern struct ring *CurrentRingp; struct ring *rp; @@ -312,6 +362,8 @@ static POLY reduction_ecart0(r,gset,needSyz,syzp) } } + if (DoCancel && (r != POLYNULL)) shouldReduceContent(r,1); + if (DebugReductionEcart&1) printf("--------------------------------------\n"); do { if (DebugReductionRed) printf("r=%s\n",POLYToString(r,'*',1)); @@ -333,8 +385,18 @@ static POLY reduction_ecart0(r,gset,needSyz,syzp) if (DebugReductionEcart & 4) printf("#"); pp = (gg->pa)[ells.ggi]; } - if (ell > 0) r = mpMult(cxx(1,0,ell,rp),r); /* r = s^ell r */ + if (ell > 0) r = ppMult(cxx(1,0,ell,rp),r); /* r = s^ell r */ r = (*reduction1)(r,pp,needSyz,&cc,&cg); + + if (DoCancel && (r != POLYNULL)) { + if (shouldReduceContent(r,0)) { + r = reduceContentOfPoly(r,&cont); + shouldReduceContent(r,1); + if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("CONT=%s ",coeffToString(cont)); + } + } + + if (needSyz) { if (ells.first) { if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp)); @@ -356,7 +418,8 @@ static POLY reduction_ecart0(r,gset,needSyz,syzp) tp = ecartCheckSyz0(cf,r_0,syz,gset,r); if (tp != POLYNULL) { fprintf(stderr,"reduction_ecart0(): sygyzy is broken. Return the Current values.\n"); - fprintf(stderr,"%s\n",POLYToString(tp,'*',1)); + fprintf(stderr,"tp=%s\n",POLYToString(tp,'*',1)); + ecartCheckSyz0_printinfo(cf,r_0,syz,gset,r); syzp->cf = cf; syzp->syz = syz; return r; @@ -364,16 +427,32 @@ static POLY reduction_ecart0(r,gset,needSyz,syzp) } } if (r ISZERO) goto ss; - /*r = ecartDivideSv(r,&se); r = r/s^? Don't do this. */ + + /* r = r/s^? Don't do this?? */ + r = ecartDivideSv(r,&se); + if (needSyz && (se > 0)) { + POLY tt; + tt = cxx(1,0,-se,rp); + cf = ppMult(tt,cf); + syz = cpMult(toSyzCoeff(tt),syz); + } + } }while (ell >= 0); ss: ; if (needSyz) { - syzp->cf = cf; /* cf is in the CurrentRingp */ - syzp->syz = syz; /* syz is in the SyzRingp */ + syzp->cf = cf; /* cf is in the ring of r */ + syzp->syz = syz; /* syz is in the syzRing of r */ } + if (DoCancel && (r != POLYNULL)) { + if (r->m->ringp->p == 0) { + r = reduceContentOfPoly(r,&cont); + if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("cont=%s ",coeffToString(cont)); + } + } + return(r); } @@ -400,6 +479,9 @@ static POLY reduction_ecart1(r,gset,needSyz,syzp) struct ecartPolyArray *gg; POLY pp; int ell; + POLY cf_o; + POLY syz_o; + POLY r_0; int se; struct coeff *cont; @@ -407,6 +489,7 @@ static POLY reduction_ecart1(r,gset,needSyz,syzp) struct ring *rp; extern struct ring *SmallRingp; + r_0 = r; gg = NULL; if (needSyz) { if (r ISZERO) { rp = CurrentRingp; } else { rp = r->m->ringp; } @@ -434,7 +517,11 @@ static POLY reduction_ecart1(r,gset,needSyz,syzp) ell = ells.ell; if (ell > 0) { if (DebugReductionEcart & 2) printf("%"); - gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL); + if (needSyz) { + gg = ecartPutPolyInG(r,gg,cf,syz); + }else{ + gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL); + } } if (ell >= 0) { if (ells.first) { @@ -443,10 +530,10 @@ static POLY reduction_ecart1(r,gset,needSyz,syzp) if (DebugReductionEcart & 4) {printf("+"); fflush(NULL);} pp = (gg->pa)[ells.ggi]; } - if (ell > 0) r = mpMult(cxx(1,0,ell,rp),r); /* r = s^ell r */ + if (ell > 0) r = ppMult(cxx(1,0,ell,rp),r); /* r = s^ell r */ r = (*reduction1)(r,pp,needSyz,&cc,&cg); - if (DoCancel && (r != POLYNULL)) { /* BUG: syzygy should be corrected. */ + if (DoCancel && (r != POLYNULL)) { if (shouldReduceContent(r,0)) { r = reduceContentOfPoly(r,&cont); shouldReduceContent(r,1); @@ -456,37 +543,96 @@ static POLY reduction_ecart1(r,gset,needSyz,syzp) if (needSyz) { if (ells.first) { + if (ell > 0) cc = ppMult(cc,cxx(1,0,ell,rp)); cf = ppMult(cc,cf); syz = cpMult(toSyzCoeff(cc),syz); - syz = ppAddv(syz,toSyzPoly(cg,ells.grade,ells.gseti)); + syz = ppAdd(syz,toSyzPoly(cg,ells.grade,ells.gseti)); }else{ - fprintf(stderr,"It has not yet implemented.\n"); - exit(10); - /* BUG: not yet */ + if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp)); + cf_o = (gg->cf)[ells.ggi]; + syz_o = (gg->syz)[ells.ggi]; + cf = ppMult(cc,cf); + cf = ppAdd(cf,ppMult(cg,cf_o)); + syz = cpMult(toSyzCoeff(cc),syz); + syz = ppAdd(syz,cpMult(toSyzCoeff(cg),syz_o)); + /* Note. 2003.07.19 */ } } + + if (DebugReductionRed) { + POLY tp; + tp = ecartCheckSyz0(cf,r_0,syz,gset,r); + tp = goDeHomogenizeS(tp); + if (tp != POLYNULL) { + fprintf(stderr,"Error: reduction_ecart1(): sygyzy is broken. Return the Current values.\n"); + fprintf(stderr,"%s\n",POLYToString(tp,'*',1)); + syzp->cf = cf; + syzp->syz = syz; + return r; + } + } + if (r ISZERO) goto ss1; r = ecartDivideSv(r,&se); /* r = r/s^? */ + + if (needSyz && (se > 0)) { /* It may not necessary because of dehomo*/ + POLY tt; + /*printf("!1/H!"); fflush(NULL);*/ /* misc-2003/ecart/t1.sm1 foo4 */ + tt = cxx(1,0,-se,rp); + cf = ppMult(tt,cf); + syz = cpMult(toSyzCoeff(tt),syz); + } + + /* For debug */ + if (DebugReductionRed && needSyz) { + POLY tp; + tp = ecartCheckSyz0(cf,r_0,syz,gset,r); + tp = goDeHomogenizeS(tp); + if (tp != POLYNULL) { + fprintf(stderr,"Error: reduction_ecart1() after divide: sygyzy is broken. Return the Current values.\n"); + fprintf(stderr,"%s\n",POLYToString(tp,'*',1)); + syzp->cf = cf; + syzp->syz = syz; + return r; + } + } + } }while (ell >= 0); ss1: ; if (needSyz) { + /* dehomogenize the syzygy. BUG, this may be inefficient. */ + cf = goDeHomogenizeS(cf); + syz = goDeHomogenizeS(syz); + /*printf("cf=%s\n",POLYToString(cf,'*',1)); + printf("syz=%s\n",POLYToString(syz,'*',1));*/ syzp->cf = cf; /* cf is in the CurrentRingp */ syzp->syz = syz; /* syz is in the SyzRingp */ - /* BUG: dehomogenize the syzygy */ - fprintf(stderr,"It has not yet implemented.\n"); - exit(10); } r = goDeHomogenizeS(r); - if (DoCancel && (r != POLYNULL)) { /* BUG: syzygy should be corrected. */ + if (DoCancel && (r != POLYNULL)) { if (r->m->ringp->p == 0) { - r = reduceContentOfPoly(r,&cont); - if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("cont=%s ",coeffToString(cont)); + r = reduceContentOfPoly(r,&cont); + if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("cont=%s ",coeffToString(cont)); } } + /* For debug */ + if (DebugReductionRed && needSyz) { + POLY tp; + tp = ecartCheckSyz0(cf,r_0,syz,gset,r); + tp = goDeHomogenizeS(tp); + if (tp != POLYNULL) { + fprintf(stderr,"Error: reduction_ecart1() last step: sygyzy is broken. Return the Current values.\n"); + fprintf(stderr,"%s\n",POLYToString(tp,'*',1)); + syzp->cf = cf; + syzp->syz = syz; + return r; + } + } + return(r); } @@ -621,7 +767,7 @@ static POLY reduction_ecart1_mod(r,gset) if (DebugReductionEcart & 4) {printf("+"); fflush(NULL);} pp = (gg->pa)[ells.ggi]; } - if (ell > 0) r = mpMult(cxx(1,0,ell,rp),r); /* r = s^ell r */ + if (ell > 0) r = ppMult(cxx(1,0,ell,rp),r); /* r = s^ell r */ r = (*reduction1)(r,pp,0,&cc,&cg); if (r ISZERO) goto ss1; r = ecartDivideSv(r,&se); /* r = r/s^? */