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