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