Annotation of OpenXM/src/kan96xx/Kan/ecart.c, Revision 1.21
1.21 ! takayama 1: /* $OpenXM: OpenXM/src/kan96xx/Kan/ecart.c,v 1.20 2003/09/21 02:19:43 takayama Exp $ */
1.1 takayama 2: #include <stdio.h>
3: #include "datatype.h"
4: #include "extern2.h"
5: #include "gradedset.h"
6:
7: #define outofmemory(a) if (a == NULL) {fprintf(stderr,"ecart.c: out of memory.\n"); exit(10);}
8: /* Local structures */
9: struct ecartReducer {
10: int ell; /* s^\ell */
11: int first; /* first =1 ==> gset else gg */
12: int grade; int gseti; /* if first==1, gradedPolySet */
13: int ggi; /* if first==0 ecartPoly [] */
14: };
15: struct ecartPolyArray {
16: int size;
17: int limit;
18: POLY *pa;
1.5 takayama 19: POLY *cf;
20: POLY *syz;
1.1 takayama 21: };
22:
23: static struct ecartReducer ecartFindReducer(POLY r,struct gradedPolySet *gset,struct ecartPolyArray *epa);
1.8 takayama 24: static struct ecartReducer ecartFindReducer_mod(POLY r,struct gradedPolySet *gset,struct ecartPolyArray *epa);
1.1 takayama 25: static int ecartCheckPoly(POLY f); /* check if it does not contain s-variables */
26: static int ecartCheckEnv(); /* check if the environment is OK for ecart div*/
1.5 takayama 27: static struct ecartPolyArray *ecartPutPolyInG(POLY g,struct ecartPolyArray *eparray,POLY cf, POLY syz);
1.1 takayama 28: static int ecartGetEll(POLY r,POLY g);
1.21 ! takayama 29: static int ecartGetEllPartial(POLY r,POLY g);
1.5 takayama 30: static POLY ecartDivideSv(POLY r,int *d);
31: /* No automatic homogenization and s is used as a standart var. */
32: static POLY reduction_ecart0(POLY r,struct gradedPolySet *gset,
33: int needSyz, struct syz0 *syzp);
34: /* Automatic homogenization and s->1 */
35: static POLY reduction_ecart1(POLY r,struct gradedPolySet *gset,
36: int needSyz, struct syz0 *syzp);
1.8 takayama 37: static POLY reduction_ecart1_mod(POLY r,struct gradedPolySet *gset);
1.5 takayama 38: static POLY ecartCheckSyz0(POLY cf,POLY r_0,POLY syz,
39: struct gradedPolySet *gg,POLY r);
1.16 takayama 40: static void ecartCheckSyz0_printinfo(POLY cf,POLY r_0,POLY syz,
41: struct gradedPolySet *gg,POLY r);
1.1 takayama 42:
1.8 takayama 43: extern int DebugReductionRed;
44: extern int TraceLift;
45: struct ring *TraceLift_ringmod;
1.9 takayama 46: extern DoCancel;
1.6 takayama 47: int DebugReductionEcart = 0;
1.11 takayama 48: extern DebugContentReduction;
1.19 takayama 49: extern int Sugar;
1.1 takayama 50:
51: /* This is used for goHomogenization */
52: extern int DegreeShifto_size;
53: extern int *DegreeShifto_vec;
1.19 takayama 54: int Ecart_sugarGrade;
1.1 takayama 55:
1.3 takayama 56: /* It is used reduction_ecart() and ecartFindReducer()
57: to determine if we homogenize in this function */
58: extern int EcartAutomaticHomogenization;
59:
1.1 takayama 60: #define LARGE 0x7fffffff
61:
62:
1.5 takayama 63: static POLY ecartDivideSv(POLY r,int *d) {
1.1 takayama 64: POLY f;
65: int k;
1.5 takayama 66: *d = 0;
1.1 takayama 67: if (r == POLYNULL) return r;
68: f = r; k = LARGE;
69: while (r != POLYNULL) {
70: if (r->m->e[0].x < k) {
71: k = r->m->e[0].x;
72: }
73: r = r->next;
74: }
1.16 takayama 75:
76: if (k > 0) {
77: *d = k;
78: return ppMult(cxx(1,0,-k,f->m->ringp),f);
79: }else{
80: return f;
81: }
82:
83: /* Do not do the below. It caused a bug. cf. misc-2003/07/ecart/b4.sm1 test2.
84: Note. 2003.8.26
85: */
86: /*
1.1 takayama 87: if (k > 0) {
1.5 takayama 88: *d = k;
1.1 takayama 89: f = r;
90: while (r != POLYNULL) {
91: r->m->e[0].x -= k;
92: r = r->next;
93: }
94: }
95: return f;
1.16 takayama 96: */
1.1 takayama 97: }
98:
99: static int ecartGetEll(POLY f,POLY g) {
100: int n,i,p;
101: MONOMIAL tf;
102: MONOMIAL tg;
103:
104: if (f ISZERO) return(-1);
105: if (g ISZERO) return(-1);
106:
107: checkRingIsR(f,g);
108:
109: if (!isSameComponent_x(f,g)) return(-1);
110: tf = f->m; tg = g->m; n = tf->ringp->n;
111: for (i=1; i<n; i++) {
112: if (tf->e[i].x < tg->e[i].x) return(-1);
113: if (tf->e[i].D < tg->e[i].D) return(-1);
114: }
115: if (tf->e[0].D < tg->e[0].D) return(-1); /* h */
1.2 takayama 116: p = tf->e[0].x - tg->e[0].x; /* H, s */
1.1 takayama 117: if (p >=0 ) return 0;
118: else return(-p);
119: }
120:
1.21 ! takayama 121: static int ecartGetEllPartial(POLY f,POLY g) {
! 122: int n,i,p;
! 123: MONOMIAL tf;
! 124: MONOMIAL tg;
! 125: int nglob; int *glob; int k;
! 126:
! 127: if (f ISZERO) return(-1);
! 128: if (g ISZERO) return(-1);
! 129:
! 130: checkRingIsR(f,g);
! 131:
! 132: if (!isSameComponent_x(f,g)) return(-1);
! 133: tf = f->m; tg = g->m; n = tf->ringp->n;
! 134: for (i=1; i<n; i++) {
! 135: if (tf->e[i].x < tg->e[i].x) return(-1);
! 136: if (tf->e[i].D < tg->e[i].D) return(-1);
! 137: }
! 138: if (tf->e[0].D < tg->e[0].D) return(-1); /* h */
! 139:
! 140: /* Find u and ell s.t. s^\ell f = u g.
! 141: When u contains variables in glob (index x var), then returns -1.
! 142: x' = glob, x'' = complement of glob. Division in Q(x'')[x'].
! 143: */
! 144: nglob = tf->ringp->partialEcart;
! 145: glob = tf->ringp->partialEcartGlobalVarX;
! 146: for (i=0; i<nglob; i++) {
! 147: k = glob[i];
! 148: if (tf->e[k].x > tg->e[k].x) return(-1);
! 149: }
! 150:
! 151: p = tf->e[0].x - tg->e[0].x; /* H, s */
! 152: if (p >=0 ) return 0;
! 153: else return(-p);
! 154: }
! 155:
1.1 takayama 156:
157: #define EP_SIZE 10
1.5 takayama 158: static struct ecartPolyArray *ecartPutPolyInG(POLY g,struct ecartPolyArray *eparray,POLY cf,POLY syz)
1.1 takayama 159: {
160: struct ecartPolyArray *a;
161: POLY *ep;
162: int i;
163: if (eparray == (struct ecartPolyArray *)NULL) {
164: a = (struct ecartPolyArray *) sGC_malloc(sizeof(struct ecartPolyArray));
165: outofmemory(a);
166: ep = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
167: outofmemory(ep);
1.5 takayama 168: a->cf = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
169: outofmemory(a->cf);
170: a->syz = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
171: outofmemory(a->syz);
1.1 takayama 172: a->limit = EP_SIZE;
173: a->size = 0;
174: a->pa = ep;
175: eparray = a;
176: }
177: if (eparray->size >= eparray->limit) {
178: a = (struct ecartPolyArray *) sGC_malloc(sizeof(struct ecartPolyArray));
179: outofmemory(a);
180: ep = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
181: outofmemory(ep);
1.5 takayama 182: a->cf = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
183: outofmemory(a->cf);
184: a->syz = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
185: outofmemory(a->syz);
1.1 takayama 186: a->limit = (eparray->limit)*2;
187: a->size = eparray->size;
188: a->pa = ep;
189: for (i=0; i<eparray->size; i++) {
190: (a->pa)[i] = (eparray->pa)[i];
1.5 takayama 191: (a->cf)[i] = (eparray->cf)[i];
192: (a->syz)[i] = (eparray->syz)[i];
1.1 takayama 193: }
194: eparray = a;
195: }
196: i = eparray->size;
197: (eparray->pa)[i] = g;
1.5 takayama 198: (eparray->cf)[i] = cf;
199: (eparray->syz)[i] = syz;
1.1 takayama 200: (eparray->size)++;
201: return eparray;
202: }
203:
1.20 takayama 204: #ifndef EXPERIMENT
1.1 takayama 205: static struct ecartReducer ecartFindReducer(POLY r,struct gradedPolySet *gset,
206: struct ecartPolyArray *epa)
207: {
208: int grd;
209: struct polySet *set;
210: int minGrade = 0;
211: int minGseti = 0;
212: int minGgi = 0;
213: int ell1 = LARGE;
214: int ell2 = LARGE;
215: int ell;
216: int i;
217: struct ecartReducer er;
218: /* Try to find a reducer in gset; */
219: grd = 0;
220: while (grd < gset->maxGrade) {
221: set = gset->polys[grd];
222: for (i=0; i<set->size; i++) {
1.2 takayama 223: if (set->gh[i] == POLYNULL) {
224: /* goHomogenize set->gh[i] */
1.4 takayama 225: if (EcartAutomaticHomogenization) {
226: set->gh[i] = goHomogenize11(set->g[i],DegreeShifto_vec,DegreeShifto_size,-1,1);
227: }else{
228: set->gh[i] = set->g[i];
229: }
1.2 takayama 230: }
231: ell = ecartGetEll(r,set->gh[i]);
1.1 takayama 232: if ((ell>=0) && (ell < ell1)) {
233: ell1 = ell;
234: minGrade = grd; minGseti=i;
235: }
236: }
237: grd++;
238: }
239: if (epa != NULL) {
240: /* Try to find in the second group. */
241: for (i=0; i< epa->size; i++) {
1.21 ! takayama 242: ell = ecartGetEllPartial(r,(epa->pa)[i]);
1.1 takayama 243: if ((ell>=0) && (ell < ell2)) {
244: ell2 = ell;
245: minGgi = i;
246: }
247: }
248: }
249:
1.18 takayama 250: if ((DebugReductionRed&1) || (DebugReductionEcart&1)) {
1.1 takayama 251: printf("ecartFindReducer(): ell1=%d, ell2=%d, minGrade=%d, minGseti=%d, minGgi=%d\n",ell1,ell2,minGrade,minGseti,minGgi);
252: }
1.19 takayama 253:
1.1 takayama 254: if (ell1 <= ell2) {
255: if (ell1 == LARGE) {
256: er.ell = -1;
257: return er;
258: }else{
259: er.ell = ell1;
260: er.first = 1;
261: er.grade = minGrade;
262: er.gseti = minGseti;
263: return er;
264: }
265: }else{
266: er.ell = ell2;
267: er.first = 0;
268: er.ggi = minGgi;
269: return er;
270: }
271: }
1.20 takayama 272: #endif
1.1 takayama 273:
1.5 takayama 274: static POLY ecartCheckSyz0(POLY cf,POLY r_0,POLY syz,
275: struct gradedPolySet *gg,POLY r)
276: {
277: POLY f;
278: int grd,i;
279: POLY q;
280: struct coeff *c;
281: f = ppMult(cf,r_0);
282: while (syz != POLYNULL) {
283: grd = syz->m->e[0].x;
284: i = syz->m->e[0].D;
285: c = syz->coeffp;
286: if (c->tag == POLY_COEFF) {
287: q = c->val.f;
288: }else{
289: q = POLYNULL;
290: }
291: f = ppAdd(f,ppMult(q,((gg->polys[grd])->g)[i]));
292: /* not gh, works only for _ecart0 */
293: syz = syz->next;
294: }
295: f = ppSub(f,r);
296: return f;
297: }
298:
1.16 takayama 299: static void ecartCheckSyz0_printinfo(POLY cf,POLY r_0,POLY syz,
300: struct gradedPolySet *gg,POLY r)
301: {
302: POLY f;
303: int grd,i;
304: POLY q;
305: struct coeff *c;
306: fprintf(stderr,"cf=%s\n",POLYToString(cf,'*',1));
307: fprintf(stderr,"r_0=%s\n",POLYToString(r_0,'*',1));
308: fprintf(stderr,"syz=%s\n",POLYToString(syz,'*',1));
309: fprintf(stderr,"r=%s\n",POLYToString(r,'*',1));
310: f = ppMult(cf,r_0);
311: while (syz != POLYNULL) {
312: grd = syz->m->e[0].x;
313: i = syz->m->e[0].D;
314: c = syz->coeffp;
315: if (c->tag == POLY_COEFF) {
316: q = c->val.f;
317: }else{
318: q = POLYNULL;
319: }
320: fprintf(stderr,"[grd,idx]=[%d,%d], %s\n",grd,i,
321: POLYToString(((gg->polys[grd])->g)[i],'*',1));
322: /* f = ppAdd(f,ppMult(q,((gg->polys[grd])->g)[i])); */
323: syz = syz->next;
324: }
325: /* f = ppSub(f,r); */
326: }
327:
1.5 takayama 328:
329: POLY reduction_ecart(r,gset,needSyz,syzp)
330: POLY r;
331: struct gradedPolySet *gset;
332: int needSyz;
333: struct syz0 *syzp; /* set */
334: {
1.8 takayama 335: POLY rn;
1.12 takayama 336: if (TraceLift && needSyz) {
337: warningGradedSet("TraceLift cannot be used to get syzygy. TraceLift is turned off.\n");
338: TraceLift = 0;
339: }
1.8 takayama 340: if (TraceLift) {
341: if (EcartAutomaticHomogenization) {
342: if (TraceLift_ringmod == NULL) {
343: warningPoly("reduction_ecart: TraceLift_ringmod is not set.\n");
344: return reduction_ecart1(r,gset,needSyz,syzp);
345: }
1.12 takayama 346: rn = reduction_ecart1_mod(r,gset);
1.8 takayama 347: if (rn == POLYNULL) return rn;
348: else return reduction_ecart1(r,gset,needSyz,syzp);
349: }else{
350: return reduction_ecart0(r,gset,needSyz,syzp);
351: }
1.5 takayama 352: }else{
1.8 takayama 353: if (EcartAutomaticHomogenization) {
354: return reduction_ecart1(r,gset,needSyz,syzp);
355: }else{
356: return reduction_ecart0(r,gset,needSyz,syzp);
357: }
1.5 takayama 358: }
359: }
360:
1.1 takayama 361: /*
362: r and gset are assumed to be (0,1)-homogenized (h-homogenized)
1.5 takayama 363: Polynomials r and gset are assumed
1.3 takayama 364: to be double homogenized (h-homogenized and s(H)-homogenized)
1.1 takayama 365: */
1.5 takayama 366: static POLY reduction_ecart0(r,gset,needSyz,syzp)
1.1 takayama 367: POLY r;
368: struct gradedPolySet *gset;
369: int needSyz;
370: struct syz0 *syzp; /* set */
371: {
372: int reduced,reduced1,reduced2;
373: int grd;
374: struct polySet *set;
375: POLY cf,syz;
376: int i;
377: POLY cc,cg;
378: struct ecartReducer ells;
379: struct ecartPolyArray *gg;
380: POLY pp;
381: int ell;
1.5 takayama 382: POLY cf_o;
383: POLY syz_o;
384: POLY r_0;
1.15 takayama 385: int se;
1.14 takayama 386: struct coeff *cont;
1.1 takayama 387:
388: extern struct ring *CurrentRingp;
389: struct ring *rp;
390:
1.5 takayama 391: r_0 = r;
1.1 takayama 392: gg = NULL;
393: if (needSyz) {
394: if (r ISZERO) { rp = CurrentRingp; } else { rp = r->m->ringp; }
395: cf = cxx(1,0,0,rp);
396: syz = ZERO;
397: }
398:
399: if (r != POLYNULL) {
1.4 takayama 400: rp = r->m->ringp;
401: if (! rp->weightedHomogenization) {
402: errorKan1("%s\n","ecart.c: the given ring must be declared with [(weightedHomogenization) 1]");
403: }
1.1 takayama 404: }
405:
1.14 takayama 406: if (DoCancel && (r != POLYNULL)) shouldReduceContent(r,1);
407:
1.20 takayama 408: if (DebugReductionEcart&1) printf("reduction_ecart0--------------------------------------\n");
1.5 takayama 409: do {
1.18 takayama 410: if (DebugReductionRed&1) printf("r=%s\n",POLYToString(r,'*',1));
1.6 takayama 411: if (DebugReductionEcart & 1) printf("r=%s+...\n",POLYToString(head(r),'*',1));
1.5 takayama 412: ells = ecartFindReducer(r,gset,gg);
413: ell = ells.ell;
414: if (ell > 0) {
1.6 takayama 415: if (DebugReductionEcart & 2) printf("*");
1.5 takayama 416: if (needSyz) {
417: gg = ecartPutPolyInG(r,gg,cf,syz);
418: }else{
419: gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL);
420: }
421: }
422: if (ell >= 0) {
423: if (ells.first) {
424: pp = ((gset->polys[ells.grade])->gh)[ells.gseti];
425: }else{
1.6 takayama 426: if (DebugReductionEcart & 4) printf("#");
1.5 takayama 427: pp = (gg->pa)[ells.ggi];
428: }
1.16 takayama 429: if (ell > 0) r = ppMult(cxx(1,0,ell,rp),r); /* r = s^ell r */
1.5 takayama 430: r = (*reduction1)(r,pp,needSyz,&cc,&cg);
1.14 takayama 431:
432: if (DoCancel && (r != POLYNULL)) {
433: if (shouldReduceContent(r,0)) {
434: r = reduceContentOfPoly(r,&cont);
435: shouldReduceContent(r,1);
436: if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("CONT=%s ",coeffToString(cont));
437: }
438: }
439:
440:
1.5 takayama 441: if (needSyz) {
442: if (ells.first) {
443: if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp));
444: cf = ppMult(cc,cf);
445: syz = cpMult(toSyzCoeff(cc),syz);
446: syz = ppAdd(syz,toSyzPoly(cg,ells.grade,ells.gseti));
447: }else{
448: if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp));
449: cf_o = (gg->cf)[ells.ggi];
450: syz_o = (gg->syz)[ells.ggi];
451: cf = ppMult(cc,cf);
452: cf = ppAdd(cf,ppMult(cg,cf_o));
453: syz = cpMult(toSyzCoeff(cc),syz);
454: syz = ppAdd(syz,cpMult(toSyzCoeff(cg),syz_o));
1.6 takayama 455: /* Note. 2003.07.19 */
1.5 takayama 456: }
1.18 takayama 457: if (DebugReductionRed && needSyz) {
1.6 takayama 458: POLY tp;
459: tp = ecartCheckSyz0(cf,r_0,syz,gset,r);
460: if (tp != POLYNULL) {
461: fprintf(stderr,"reduction_ecart0(): sygyzy is broken. Return the Current values.\n");
1.16 takayama 462: fprintf(stderr,"tp=%s\n",POLYToString(tp,'*',1));
463: ecartCheckSyz0_printinfo(cf,r_0,syz,gset,r);
1.6 takayama 464: syzp->cf = cf;
465: syzp->syz = syz;
466: return r;
467: }
1.5 takayama 468: }
469: }
470: if (r ISZERO) goto ss;
1.15 takayama 471:
472: /* r = r/s^? Don't do this?? */
473: r = ecartDivideSv(r,&se);
474: if (needSyz && (se > 0)) {
475: POLY tt;
476: tt = cxx(1,0,-se,rp);
1.16 takayama 477: cf = ppMult(tt,cf);
1.15 takayama 478: syz = cpMult(toSyzCoeff(tt),syz);
479: }
480:
1.5 takayama 481: }
482: }while (ell >= 0);
483:
484: ss: ;
485: if (needSyz) {
1.16 takayama 486: syzp->cf = cf; /* cf is in the ring of r */
487: syzp->syz = syz; /* syz is in the syzRing of r */
1.14 takayama 488: }
489:
490: if (DoCancel && (r != POLYNULL)) {
491: if (r->m->ringp->p == 0) {
492: r = reduceContentOfPoly(r,&cont);
493: if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("cont=%s ",coeffToString(cont));
494: }
1.5 takayama 495: }
496:
497: return(r);
498: }
499:
500: /*
501: r and gset are assumed to be (0,1)-homogenized (h-homogenized)
502: r and gset are not assumed
503: to be double homogenized (h-homogenized and s(H)-homogenized)
504: They are automatically s(H)-homogenized, and s is set to 1
505: when this function returns.
506: */
507: static POLY reduction_ecart1(r,gset,needSyz,syzp)
508: POLY r;
509: struct gradedPolySet *gset;
510: int needSyz;
511: struct syz0 *syzp; /* set */
512: {
513: int reduced,reduced1,reduced2;
514: int grd;
515: struct polySet *set;
516: POLY cf,syz;
517: int i;
518: POLY cc,cg;
519: struct ecartReducer ells;
520: struct ecartPolyArray *gg;
521: POLY pp;
522: int ell;
1.12 takayama 523: POLY cf_o;
524: POLY syz_o;
525: POLY r_0;
1.18 takayama 526: POLY r0;
1.5 takayama 527: int se;
1.9 takayama 528: struct coeff *cont;
1.5 takayama 529:
530: extern struct ring *CurrentRingp;
531: struct ring *rp;
1.9 takayama 532: extern struct ring *SmallRingp;
1.5 takayama 533:
1.12 takayama 534: r_0 = r;
1.5 takayama 535: gg = NULL;
536: if (needSyz) {
537: if (r ISZERO) { rp = CurrentRingp; } else { rp = r->m->ringp; }
538: cf = cxx(1,0,0,rp);
539: syz = ZERO;
540: }
541:
542: if (r != POLYNULL) {
543: rp = r->m->ringp;
544: if (! rp->weightedHomogenization) {
545: errorKan1("%s\n","ecart.c: the given ring must be declared with [(weightedHomogenization) 1]");
546: }
1.3 takayama 547: }
1.5 takayama 548:
549: r = goHomogenize11(r,DegreeShifto_vec,DegreeShifto_size,-1,1);
1.1 takayama 550: /* 1 means homogenize only s */
1.10 takayama 551: if (DoCancel && (r != POLYNULL)) shouldReduceContent(r,1);
1.1 takayama 552:
1.20 takayama 553: if (DebugReductionEcart&1) printf("reduction_ecart1=======================================\n");
1.1 takayama 554: do {
1.18 takayama 555: if (DebugReductionRed & 1) printf("(ecart1(d)) r=%s\n",POLYToString(r,'*',1));
1.7 takayama 556: if (DebugReductionEcart & 1) printf("r=%s+,,,\n",POLYToString(head(r),'*',1));
557:
1.1 takayama 558: ells = ecartFindReducer(r,gset,gg);
559: ell = ells.ell;
560: if (ell > 0) {
1.7 takayama 561: if (DebugReductionEcart & 2) printf("%");
1.16 takayama 562: if (needSyz) {
563: gg = ecartPutPolyInG(r,gg,cf,syz);
564: }else{
565: gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL);
566: }
1.1 takayama 567: }
568: if (ell >= 0) {
569: if (ells.first) {
570: pp = ((gset->polys[ells.grade])->gh)[ells.gseti];
571: }else{
1.9 takayama 572: if (DebugReductionEcart & 4) {printf("+"); fflush(NULL);}
1.1 takayama 573: pp = (gg->pa)[ells.ggi];
574: }
1.16 takayama 575: if (ell > 0) r = ppMult(cxx(1,0,ell,rp),r); /* r = s^ell r */
1.18 takayama 576: r0 = r;
1.1 takayama 577: r = (*reduction1)(r,pp,needSyz,&cc,&cg);
1.18 takayama 578: if (DebugReductionEcart & 8) {
579: if (ell > 0) {printf("ell+ "); fflush(NULL);
580: }else {
581: printf("ell0 "); fflush(NULL);
582: if ((*mmLarger)(r,r0) >= 1) {
583: printf("error in reduction.");
584: printf(" r0=%s\n",POLYToString(r0,'*',1));
585: printf("==>r=%s\n",POLYToString(r,'*',1));
586: getchar(); getchar();
587: }
588: }
589: }
1.10 takayama 590:
1.12 takayama 591: if (DoCancel && (r != POLYNULL)) {
1.10 takayama 592: if (shouldReduceContent(r,0)) {
593: r = reduceContentOfPoly(r,&cont);
594: shouldReduceContent(r,1);
1.11 takayama 595: if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("CONT=%s ",coeffToString(cont));
1.10 takayama 596: }
597: }
598:
1.1 takayama 599: if (needSyz) {
600: if (ells.first) {
1.16 takayama 601: if (ell > 0) cc = ppMult(cc,cxx(1,0,ell,rp));
1.1 takayama 602: cf = ppMult(cc,cf);
603: syz = cpMult(toSyzCoeff(cc),syz);
1.16 takayama 604: syz = ppAdd(syz,toSyzPoly(cg,ells.grade,ells.gseti));
1.1 takayama 605: }else{
1.12 takayama 606: if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp));
607: cf_o = (gg->cf)[ells.ggi];
608: syz_o = (gg->syz)[ells.ggi];
609: cf = ppMult(cc,cf);
610: cf = ppAdd(cf,ppMult(cg,cf_o));
611: syz = cpMult(toSyzCoeff(cc),syz);
612: syz = ppAdd(syz,cpMult(toSyzCoeff(cg),syz_o));
613: /* Note. 2003.07.19 */
1.1 takayama 614: }
615: }
1.12 takayama 616:
1.18 takayama 617: if (DebugReductionRed && needSyz) {
1.12 takayama 618: POLY tp;
619: tp = ecartCheckSyz0(cf,r_0,syz,gset,r);
1.16 takayama 620: tp = goDeHomogenizeS(tp);
1.12 takayama 621: if (tp != POLYNULL) {
1.16 takayama 622: fprintf(stderr,"Error: reduction_ecart1(): sygyzy is broken. Return the Current values.\n");
1.12 takayama 623: fprintf(stderr,"%s\n",POLYToString(tp,'*',1));
624: syzp->cf = cf;
625: syzp->syz = syz;
626: return r;
627: }
628: }
629:
1.5 takayama 630: if (r ISZERO) goto ss1;
631: r = ecartDivideSv(r,&se); /* r = r/s^? */
1.15 takayama 632:
1.16 takayama 633: if (needSyz && (se > 0)) { /* It may not necessary because of dehomo*/
634: POLY tt;
635: /*printf("!1/H!"); fflush(NULL);*/ /* misc-2003/ecart/t1.sm1 foo4 */
636: tt = cxx(1,0,-se,rp);
637: cf = ppMult(tt,cf);
638: syz = cpMult(toSyzCoeff(tt),syz);
639: }
640:
641: /* For debug */
642: if (DebugReductionRed && needSyz) {
643: POLY tp;
644: tp = ecartCheckSyz0(cf,r_0,syz,gset,r);
645: tp = goDeHomogenizeS(tp);
646: if (tp != POLYNULL) {
647: fprintf(stderr,"Error: reduction_ecart1() after divide: sygyzy is broken. Return the Current values.\n");
648: fprintf(stderr,"%s\n",POLYToString(tp,'*',1));
649: syzp->cf = cf;
650: syzp->syz = syz;
651: return r;
652: }
653: }
1.15 takayama 654:
1.1 takayama 655: }
656: }while (ell >= 0);
657:
1.5 takayama 658: ss1: ;
1.1 takayama 659: if (needSyz) {
1.12 takayama 660: /* dehomogenize the syzygy. BUG, this may be inefficient. */
1.16 takayama 661: cf = goDeHomogenizeS(cf);
662: syz = goDeHomogenizeS(syz);
663: /*printf("cf=%s\n",POLYToString(cf,'*',1));
664: printf("syz=%s\n",POLYToString(syz,'*',1));*/
1.1 takayama 665: syzp->cf = cf; /* cf is in the CurrentRingp */
666: syzp->syz = syz; /* syz is in the SyzRingp */
667: }
1.8 takayama 668:
669: r = goDeHomogenizeS(r);
1.12 takayama 670: if (DoCancel && (r != POLYNULL)) {
1.10 takayama 671: if (r->m->ringp->p == 0) {
1.16 takayama 672: r = reduceContentOfPoly(r,&cont);
673: if (DebugReductionEcart || DebugReductionRed || DebugContentReduction) printf("cont=%s ",coeffToString(cont));
674: }
675: }
676:
677: /* For debug */
678: if (DebugReductionRed && needSyz) {
679: POLY tp;
680: tp = ecartCheckSyz0(cf,r_0,syz,gset,r);
681: tp = goDeHomogenizeS(tp);
682: if (tp != POLYNULL) {
683: fprintf(stderr,"Error: reduction_ecart1() last step: sygyzy is broken. Return the Current values.\n");
684: fprintf(stderr,"%s\n",POLYToString(tp,'*',1));
685: syzp->cf = cf;
686: syzp->syz = syz;
687: return r;
1.10 takayama 688: }
1.9 takayama 689: }
1.8 takayama 690:
691: return(r);
692: }
693:
694: /* Functions for trace lift */
695: static struct ecartReducer ecartFindReducer_mod(POLY r,
696: struct gradedPolySet *gset,
697: struct ecartPolyArray *epa)
698: {
699: int grd;
700: struct polySet *set;
701: int minGrade = 0;
702: int minGseti = 0;
703: int minGgi = 0;
704: int ell1 = LARGE;
705: int ell2 = LARGE;
706: int ell;
707: int i;
708: struct ecartReducer er;
709: /* Try to find a reducer in gset; */
710: grd = 0;
711: while (grd < gset->maxGrade) {
712: set = gset->polys[grd];
713: for (i=0; i<set->size; i++) {
714: if (set->gh[i] == POLYNULL) {
715: /* goHomogenize set->gh[i] */
716: if (EcartAutomaticHomogenization) {
717: set->gh[i] = goHomogenize11(set->g[i],DegreeShifto_vec,DegreeShifto_size,-1,1);
718: }else{
719: set->gh[i] = set->g[i];
720: }
721: }
722: if (TraceLift && (set->gmod[i] == POLYNULL)) {
723: set->gmod[i] = modulop(set->gh[i],TraceLift_ringmod);
724: }
725: if (TraceLift) {
726: ell = ecartGetEll(r,set->gmod[i]);
727: }else{
728: ell = ecartGetEll(r,set->gh[i]);
729: }
730: if ((ell>=0) && (ell < ell1)) {
731: ell1 = ell;
732: minGrade = grd; minGseti=i;
733: }
734: }
735: grd++;
736: }
737: if (epa != NULL) {
738: /* Try to find in the second group. */
739: for (i=0; i< epa->size; i++) {
1.21 ! takayama 740: ell = ecartGetEllPartial(r,(epa->pa)[i]);
1.8 takayama 741: if ((ell>=0) && (ell < ell2)) {
742: ell2 = ell;
743: minGgi = i;
744: }
745: }
746: }
747:
1.18 takayama 748: if ((DebugReductionRed & 1)|| (DebugReductionEcart&1)) {
1.8 takayama 749: printf("ecartFindReducer_mod(): ell1=%d, ell2=%d, minGrade=%d, minGseti=%d, minGgi=%d, p=%d\n",ell1,ell2,minGrade,minGseti,minGgi,TraceLift_ringmod->p);
750: }
751: if (ell1 <= ell2) {
752: if (ell1 == LARGE) {
753: er.ell = -1;
754: return er;
755: }else{
756: er.ell = ell1;
757: er.first = 1;
758: er.grade = minGrade;
759: er.gseti = minGseti;
760: return er;
761: }
762: }else{
763: er.ell = ell2;
764: er.first = 0;
765: er.ggi = minGgi;
766: return er;
767: }
768: }
769:
770: static POLY reduction_ecart1_mod(r,gset)
771: POLY r;
772: struct gradedPolySet *gset;
773: {
774: int reduced,reduced1,reduced2;
775: int grd;
776: struct polySet *set;
777: int i;
778: POLY cc,cg;
779: struct ecartReducer ells;
780: struct ecartPolyArray *gg;
781: POLY pp;
1.18 takayama 782: POLY r0;
1.8 takayama 783: int ell;
784: int se;
785:
786: extern struct ring *CurrentRingp;
787: struct ring *rp;
788:
789: gg = NULL;
790:
791: if (r != POLYNULL) {
792: rp = r->m->ringp;
793: if (! rp->weightedHomogenization) {
794: errorKan1("%s\n","ecart.c: the given ring must be declared with [(weightedHomogenization) 1]");
795: }
796: }
797:
798: r = goHomogenize11(r,DegreeShifto_vec,DegreeShifto_size,-1,1);
799: /* 1 means homogenize only s */
800: /*printf("r=%s (mod 0)\n",POLYToString(head(r),'*',1));
801: KshowRing(TraceLift_ringmod); **/
802:
803: r = modulop(r,TraceLift_ringmod);
1.17 takayama 804: if (r != POLYNULL) rp = r->m->ringp; /* reset rp */
1.8 takayama 805:
806: /* printf("r=%s (mod p)\n",POLYToString(head(r),'*',1)); **/
807:
808: if (DebugReductionEcart&1) printf("=====================================mod\n");
809: do {
1.18 takayama 810: if (DebugReductionRed & 1) printf("(ecart1_mod(d)) r=%s\n",POLYToString(r,'*',1));
1.8 takayama 811: if (DebugReductionEcart & 1) printf("r=%s+,,,\n",POLYToString(head(r),'*',1));
812:
813: ells = ecartFindReducer_mod(r,gset,gg);
814: ell = ells.ell;
815: if (ell > 0) {
816: if (DebugReductionEcart & 2) printf("%");
817: gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL);
818: }
819: if (ell >= 0) {
820: if (ells.first) {
821: pp = ((gset->polys[ells.grade])->gmod)[ells.gseti];
822: }else{
1.9 takayama 823: if (DebugReductionEcart & 4) {printf("+"); fflush(NULL);}
1.8 takayama 824: pp = (gg->pa)[ells.ggi];
825: }
1.16 takayama 826: if (ell > 0) r = ppMult(cxx(1,0,ell,rp),r); /* r = s^ell r */
1.18 takayama 827:
828: r0 = r;
1.8 takayama 829: r = (*reduction1)(r,pp,0,&cc,&cg);
1.18 takayama 830:
831: if (DebugReductionEcart & 8) {
832: if (ell > 0) {printf("ell+ "); fflush(NULL);
833: }else {
834: printf("ell0 "); fflush(NULL);
835: if ((*mmLarger)(r,r0) >= 1) {
836: printf("error in reduction.");
837: printf(" r0=%s\n",POLYToString(r0,'*',1));
838: printf("==>r=%s\n",POLYToString(r,'*',1));
839: getchar(); getchar();
840: }
841: }
842: }
843:
1.8 takayama 844: if (r ISZERO) goto ss1;
845: r = ecartDivideSv(r,&se); /* r = r/s^? */
846: }
847: }while (ell >= 0);
848:
849: ss1: ;
1.5 takayama 850:
1.1 takayama 851: r = goDeHomogenizeS(r);
1.5 takayama 852:
1.1 takayama 853: return(r);
1.10 takayama 854: }
1.20 takayama 855:
856: #ifdef EXPERIMENT
857: static struct ecartReducer ecartFindReducer(POLY r,struct gradedPolySet *gset,
858: struct ecartPolyArray *epa)
859: {
860: int grd;
861: struct polySet *set;
862: int minGrade = 0;
863: int minGseti = 0;
864: int minGgi = 0;
865: int ell1 = LARGE;
866: int ell2 = LARGE;
867: int ell;
868: int i;
869: struct ecartReducer er;
870: int ellmin ;
871: ell1 = ell2 = ellmin = LARGE;
872: /* Try to find a reducer in gset; */
873: grd = 0;
874: while (grd < gset->maxGrade) {
875: set = gset->polys[grd];
876: for (i=0; i<set->size; i++) {
877: if (set->gh[i] == POLYNULL) {
878: /* goHomogenize set->gh[i] */
879: if (EcartAutomaticHomogenization) {
880: set->gh[i] = goHomogenize11(set->g[i],DegreeShifto_vec,DegreeShifto_size,-1,1);
881: }else{
882: set->gh[i] = set->g[i];
883: }
884: }
885: ell = ecartGetEll(r,set->gh[i]);
886: if ((ell>=0) && (ell < ell1)) {
887: ell1 = ell;
888: minGrade = grd; minGseti=i;
889: }
890: }
891: grd++;
892: }
893: if (epa != NULL) {
894: /* Try to find in the second group. */
895: for (i=0; i< epa->size; i++) {
896: ell = ecartGetEll(r,(epa->pa)[i]);
897: if ((ell>=0) && (ell < ell2)) {
898: ell2 = ell;
899: minGgi = i;
900: }
901: }
902: }
903:
904: if ((DebugReductionRed&1) || (DebugReductionEcart&1)) {
905: printf("ecartFindReducer(): ell1=%d, ell2=%d, minGrade=%d, minGseti=%d, minGgi=%d\n",ell1,ell2,minGrade,minGseti,minGgi);
906: }
907:
908: #ifdef EXPERIMENTAL1
909: if (Sugar) { /* experimental */
910: if (ell1 <= ell2) {
911: if (ell1 == LARGE) {
912: er.ell = -1;
913: return er;
914: }else{
915: int new_s;
916: er.ell = ell1;
917: er.first = 1;
918: er.grade = minGrade;
919: er.gseti = minGseti;
920: /* reduce if and only if Ecart_sugarGrade does not increase. */
921: new_s = grade_gen(r)-grade_gen(gset->polys[minGrade]->gh[minGseti]);
922: if (new_s + minGrade <= Ecart_sugarGrade+1) {
923: return er;
924: }else{
925: printf("new_s=%d, minGrade=%d, sugarGrade=%d\n",new_s,minGrade,
926: Ecart_sugarGrade);
927: er.ell = -1;
928: return er;
929: }
930: }
931: }else{
932: er.ell = ell2;
933: er.first = 0;
934: er.ggi = minGgi;
935: return er;
936: }
937: }
938: #endif
939:
940: if ((ell1 == LARGE) && (ell2 == LARGE)) {
941: er.ell = -1;
942: return er;
943: }
944:
945: if (ell1 < ell2) {
946: ellmin = ell1;
947: }else{
948: ellmin = ell2;
949: }
950: {
951: #define M 100
952: int aminGrade[M];
953: int aminGseti[M];
954: int aminGgi[M];
955: int sp1,sp2;
956: int i;
957: int gmin;
958: int gtmp;
959: sp1 = sp2 = 0;
960:
961: if (ell1 == ellmin) {
962: grd = 0;
963: while (grd < gset->maxGrade) {
964: set = gset->polys[grd];
965: for (i=0; i<set->size; i++) {
966: ell = ecartGetEll(r,set->gh[i]);
967: if (ell == ellmin) {
968: aminGrade[sp1] = grd; aminGseti[sp1]=i;
969: sp1++;
970: if (sp1 >= M) {
971: fprintf(stderr,"aminGrade, buffer overflow. sp1 is set to 1.\n");
972: sp1 = 1;
973: }
974: }
975: }
976: grd++;
977: }
978: }
979: if (ell2 == ellmin) {
980: if (epa != NULL) {
981: for (i=0; i< epa->size; i++) {
982: ell = ecartGetEll(r,(epa->pa)[i]);
983: if (ell == ellmin) {
984: aminGgi[sp2] = i;
985: sp2++;
986: if (sp2 >= M) {
987: fprintf(stderr,"aminGgi, buffer overflow. sp2 is set to 1.\n");
988: sp2 = 1;
989: }
990: }
991: }
992: }
993: }
994: /* print summary */
995: printf("summary -----------------------\n");
996: for (i=0; i<sp1; i++) {
997: printf("i=%d: minGrade=%d, minGseti=%d,",i,aminGrade[i],aminGseti[i]);
998: printf("grade=%d\n",grade_gen(gset->polys[aminGrade[i]]->gh[aminGseti[i]]));
999: }
1000: gmin = LARGE;
1001: for (i=0; i<sp2; i++) {
1002: printf("i=%d: minGgi=%d,",i,aminGgi[i]);
1003: gtmp = grade_gen((epa->pa)[aminGgi[i]]);
1004: printf("grade=%d\n",gtmp);
1005: if (gtmp < gmin) {
1006: gmin = gtmp;
1007: minGgi = aminGgi[i];
1008: }
1009: }
1010:
1011: if (ell1 <= ell2) {
1012: er.ell = ell1;
1013: er.first = 1;
1014: er.grade = minGrade;
1015: er.gseti = minGseti;
1016: return er;
1017: }else{
1018: er.ell = ell2;
1019: er.first = 0;
1020: er.ggi = minGgi;
1021: return er;
1022: }
1023:
1024: }
1025:
1026: if (ell1 <= ell2) {
1027: er.ell = ell1;
1028: er.first = 1;
1029: er.grade = minGrade;
1030: er.gseti = minGseti;
1031: return er;
1032: }else{
1033: er.ell = ell2;
1034: er.first = 0;
1035: er.ggi = minGgi;
1036: return er;
1037: }
1038: }
1039: #endif
1.10 takayama 1040:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>