Annotation of OpenXM/src/kan96xx/Kan/ecart.c, Revision 1.6
1.6 ! takayama 1: /* $OpenXM: OpenXM/src/kan96xx/Kan/ecart.c,v 1.5 2003/07/19 06:03:57 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);
24: static int ecartCheckPoly(POLY f); /* check if it does not contain s-variables */
25: static int ecartCheckEnv(); /* check if the environment is OK for ecart div*/
1.5 takayama 26: static struct ecartPolyArray *ecartPutPolyInG(POLY g,struct ecartPolyArray *eparray,POLY cf, POLY syz);
1.1 takayama 27: static int ecartGetEll(POLY r,POLY g);
1.5 takayama 28: static POLY ecartDivideSv(POLY r,int *d);
29: /* No automatic homogenization and s is used as a standart var. */
30: static POLY reduction_ecart0(POLY r,struct gradedPolySet *gset,
31: int needSyz, struct syz0 *syzp);
32: /* Automatic homogenization and s->1 */
33: static POLY reduction_ecart1(POLY r,struct gradedPolySet *gset,
34: int needSyz, struct syz0 *syzp);
35: static POLY ecartCheckSyz0(POLY cf,POLY r_0,POLY syz,
36: struct gradedPolySet *gg,POLY r);
1.1 takayama 37:
38: extern int DebugReductionRed;
1.6 ! takayama 39: int DebugReductionEcart = 0;
1.1 takayama 40:
41: /* This is used for goHomogenization */
42: extern int DegreeShifto_size;
43: extern int *DegreeShifto_vec;
44:
1.3 takayama 45: /* It is used reduction_ecart() and ecartFindReducer()
46: to determine if we homogenize in this function */
47: extern int EcartAutomaticHomogenization;
48:
1.1 takayama 49: #define LARGE 0x7fffffff
50:
51:
1.5 takayama 52: static POLY ecartDivideSv(POLY r,int *d) {
1.1 takayama 53: POLY f;
54: int k;
1.5 takayama 55: *d = 0;
1.1 takayama 56: if (r == POLYNULL) return r;
57: f = r; k = LARGE;
58: while (r != POLYNULL) {
59: if (r->m->e[0].x < k) {
60: k = r->m->e[0].x;
61: }
62: r = r->next;
63: }
64: if (k > 0) {
1.5 takayama 65: *d = k;
1.1 takayama 66: f = r;
67: while (r != POLYNULL) {
68: r->m->e[0].x -= k;
69: r = r->next;
70: }
71: }
72: return f;
73: }
74:
75: static int ecartGetEll(POLY f,POLY g) {
76: int n,i,p;
77: MONOMIAL tf;
78: MONOMIAL tg;
79:
80: if (f ISZERO) return(-1);
81: if (g ISZERO) return(-1);
82:
83: checkRingIsR(f,g);
84:
85: if (!isSameComponent_x(f,g)) return(-1);
86: tf = f->m; tg = g->m; n = tf->ringp->n;
87: for (i=1; i<n; i++) {
88: if (tf->e[i].x < tg->e[i].x) return(-1);
89: if (tf->e[i].D < tg->e[i].D) return(-1);
90: }
91: if (tf->e[0].D < tg->e[0].D) return(-1); /* h */
1.2 takayama 92: p = tf->e[0].x - tg->e[0].x; /* H, s */
1.1 takayama 93: if (p >=0 ) return 0;
94: else return(-p);
95: }
96:
97:
98: #define EP_SIZE 10
1.5 takayama 99: static struct ecartPolyArray *ecartPutPolyInG(POLY g,struct ecartPolyArray *eparray,POLY cf,POLY syz)
1.1 takayama 100: {
101: struct ecartPolyArray *a;
102: POLY *ep;
103: int i;
104: if (eparray == (struct ecartPolyArray *)NULL) {
105: a = (struct ecartPolyArray *) sGC_malloc(sizeof(struct ecartPolyArray));
106: outofmemory(a);
107: ep = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
108: outofmemory(ep);
1.5 takayama 109: a->cf = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
110: outofmemory(a->cf);
111: a->syz = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
112: outofmemory(a->syz);
1.1 takayama 113: a->limit = EP_SIZE;
114: a->size = 0;
115: a->pa = ep;
116: eparray = a;
117: }
118: if (eparray->size >= eparray->limit) {
119: a = (struct ecartPolyArray *) sGC_malloc(sizeof(struct ecartPolyArray));
120: outofmemory(a);
121: ep = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
122: outofmemory(ep);
1.5 takayama 123: a->cf = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
124: outofmemory(a->cf);
125: a->syz = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
126: outofmemory(a->syz);
1.1 takayama 127: a->limit = (eparray->limit)*2;
128: a->size = eparray->size;
129: a->pa = ep;
130: for (i=0; i<eparray->size; i++) {
131: (a->pa)[i] = (eparray->pa)[i];
1.5 takayama 132: (a->cf)[i] = (eparray->cf)[i];
133: (a->syz)[i] = (eparray->syz)[i];
1.1 takayama 134: }
135: eparray = a;
136: }
137: i = eparray->size;
138: (eparray->pa)[i] = g;
1.5 takayama 139: (eparray->cf)[i] = cf;
140: (eparray->syz)[i] = syz;
1.1 takayama 141: (eparray->size)++;
142: return eparray;
143: }
144:
145: static struct ecartReducer ecartFindReducer(POLY r,struct gradedPolySet *gset,
146: struct ecartPolyArray *epa)
147: {
148: int grd;
149: struct polySet *set;
150: int minGrade = 0;
151: int minGseti = 0;
152: int minGgi = 0;
153: int ell1 = LARGE;
154: int ell2 = LARGE;
155: int ell;
156: int i;
157: struct ecartReducer er;
158: /* Try to find a reducer in gset; */
159: grd = 0;
160: while (grd < gset->maxGrade) {
161: set = gset->polys[grd];
162: for (i=0; i<set->size; i++) {
1.2 takayama 163: if (set->gh[i] == POLYNULL) {
164: /* goHomogenize set->gh[i] */
1.4 takayama 165: if (EcartAutomaticHomogenization) {
166: set->gh[i] = goHomogenize11(set->g[i],DegreeShifto_vec,DegreeShifto_size,-1,1);
167: }else{
168: set->gh[i] = set->g[i];
169: }
1.2 takayama 170: }
171: ell = ecartGetEll(r,set->gh[i]);
1.1 takayama 172: if ((ell>=0) && (ell < ell1)) {
173: ell1 = ell;
174: minGrade = grd; minGseti=i;
175: }
176: }
177: grd++;
178: }
179: if (epa != NULL) {
180: /* Try to find in the second group. */
181: for (i=0; i< epa->size; i++) {
182: ell = ecartGetEll(r,(epa->pa)[i]);
183: if ((ell>=0) && (ell < ell2)) {
184: ell2 = ell;
185: minGgi = i;
186: }
187: }
188: }
189:
1.6 ! takayama 190: if (DebugReductionRed || (DebugReductionEcart&1)) {
1.1 takayama 191: printf("ecartFindReducer(): ell1=%d, ell2=%d, minGrade=%d, minGseti=%d, minGgi=%d\n",ell1,ell2,minGrade,minGseti,minGgi);
192: }
193: if (ell1 <= ell2) {
194: if (ell1 == LARGE) {
195: er.ell = -1;
196: return er;
197: }else{
198: er.ell = ell1;
199: er.first = 1;
200: er.grade = minGrade;
201: er.gseti = minGseti;
202: return er;
203: }
204: }else{
205: er.ell = ell2;
206: er.first = 0;
207: er.ggi = minGgi;
208: return er;
209: }
210: }
211:
1.5 takayama 212: static POLY ecartCheckSyz0(POLY cf,POLY r_0,POLY syz,
213: struct gradedPolySet *gg,POLY r)
214: {
215: POLY f;
216: int grd,i;
217: POLY q;
218: struct coeff *c;
219: f = ppMult(cf,r_0);
220: while (syz != POLYNULL) {
221: grd = syz->m->e[0].x;
222: i = syz->m->e[0].D;
223: c = syz->coeffp;
224: if (c->tag == POLY_COEFF) {
225: q = c->val.f;
226: }else{
227: q = POLYNULL;
228: }
229: f = ppAdd(f,ppMult(q,((gg->polys[grd])->g)[i]));
230: /* not gh, works only for _ecart0 */
231: syz = syz->next;
232: }
233: f = ppSub(f,r);
234: return f;
235: }
236:
237:
238: POLY reduction_ecart(r,gset,needSyz,syzp)
239: POLY r;
240: struct gradedPolySet *gset;
241: int needSyz;
242: struct syz0 *syzp; /* set */
243: {
244: if (EcartAutomaticHomogenization) {
245: return reduction_ecart1(r,gset,needSyz,syzp);
246: }else{
247: return reduction_ecart0(r,gset,needSyz,syzp);
248: }
249: }
250:
1.1 takayama 251: /*
252: r and gset are assumed to be (0,1)-homogenized (h-homogenized)
1.5 takayama 253: Polynomials r and gset are assumed
1.3 takayama 254: to be double homogenized (h-homogenized and s(H)-homogenized)
1.1 takayama 255: */
1.5 takayama 256: static POLY reduction_ecart0(r,gset,needSyz,syzp)
1.1 takayama 257: POLY r;
258: struct gradedPolySet *gset;
259: int needSyz;
260: struct syz0 *syzp; /* set */
261: {
262: int reduced,reduced1,reduced2;
263: int grd;
264: struct polySet *set;
265: POLY cf,syz;
266: int i;
267: POLY cc,cg;
268: struct ecartReducer ells;
269: struct ecartPolyArray *gg;
270: POLY pp;
271: int ell;
1.5 takayama 272: POLY cf_o;
273: POLY syz_o;
274: POLY r_0;
1.1 takayama 275:
276: extern struct ring *CurrentRingp;
277: struct ring *rp;
278:
1.5 takayama 279: r_0 = r;
1.1 takayama 280: gg = NULL;
281: if (needSyz) {
282: if (r ISZERO) { rp = CurrentRingp; } else { rp = r->m->ringp; }
283: cf = cxx(1,0,0,rp);
284: syz = ZERO;
285: }
286:
287: if (r != POLYNULL) {
1.4 takayama 288: rp = r->m->ringp;
289: if (! rp->weightedHomogenization) {
290: errorKan1("%s\n","ecart.c: the given ring must be declared with [(weightedHomogenization) 1]");
291: }
1.1 takayama 292: }
293:
1.6 ! takayama 294: if (DebugReductionEcart&1) printf("--------------------------------------\n");
1.5 takayama 295: do {
296: if (DebugReductionRed) printf("r=%s\n",POLYToString(r,'*',1));
1.6 ! takayama 297: if (DebugReductionEcart & 1) printf("r=%s+...\n",POLYToString(head(r),'*',1));
1.5 takayama 298: ells = ecartFindReducer(r,gset,gg);
299: ell = ells.ell;
300: if (ell > 0) {
1.6 ! takayama 301: if (DebugReductionEcart & 2) printf("*");
1.5 takayama 302: if (needSyz) {
303: gg = ecartPutPolyInG(r,gg,cf,syz);
304: }else{
305: gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL);
306: }
307: }
308: if (ell >= 0) {
309: if (ells.first) {
310: pp = ((gset->polys[ells.grade])->gh)[ells.gseti];
311: }else{
1.6 ! takayama 312: if (DebugReductionEcart & 4) printf("#");
1.5 takayama 313: pp = (gg->pa)[ells.ggi];
314: }
315: if (ell > 0) r = mpMult(cxx(1,0,ell,rp),r); /* r = s^ell r */
316: r = (*reduction1)(r,pp,needSyz,&cc,&cg);
317: if (needSyz) {
318: if (ells.first) {
319: if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp));
320: cf = ppMult(cc,cf);
321: syz = cpMult(toSyzCoeff(cc),syz);
322: syz = ppAdd(syz,toSyzPoly(cg,ells.grade,ells.gseti));
323: }else{
324: if (ell >0) cc = ppMult(cc,cxx(1,0,ell,rp));
325: cf_o = (gg->cf)[ells.ggi];
326: syz_o = (gg->syz)[ells.ggi];
327: cf = ppMult(cc,cf);
328: cf = ppAdd(cf,ppMult(cg,cf_o));
329: syz = cpMult(toSyzCoeff(cc),syz);
330: syz = ppAdd(syz,cpMult(toSyzCoeff(cg),syz_o));
1.6 ! takayama 331: /* Note. 2003.07.19 */
1.5 takayama 332: }
333: if (DebugReductionRed) {
1.6 ! takayama 334: POLY tp;
! 335: tp = ecartCheckSyz0(cf,r_0,syz,gset,r);
! 336: if (tp != POLYNULL) {
! 337: fprintf(stderr,"reduction_ecart0(): sygyzy is broken. Return the Current values.\n");
! 338: fprintf(stderr,"%s\n",POLYToString(tp,'*',1));
! 339: syzp->cf = cf;
! 340: syzp->syz = syz;
! 341: return r;
! 342: }
1.5 takayama 343: }
344: }
345: if (r ISZERO) goto ss;
346: /*r = ecartDivideSv(r,&se); r = r/s^? Don't do this. */
347: }
348: }while (ell >= 0);
349:
350: ss: ;
351: if (needSyz) {
352: syzp->cf = cf; /* cf is in the CurrentRingp */
353: syzp->syz = syz; /* syz is in the SyzRingp */
354: }
355:
356: return(r);
357: }
358:
359: /*
360: r and gset are assumed to be (0,1)-homogenized (h-homogenized)
361: r and gset are not assumed
362: to be double homogenized (h-homogenized and s(H)-homogenized)
363: They are automatically s(H)-homogenized, and s is set to 1
364: when this function returns.
365: */
366: static POLY reduction_ecart1(r,gset,needSyz,syzp)
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;
382: int se;
383:
384: extern struct ring *CurrentRingp;
385: struct ring *rp;
386:
387: gg = NULL;
388: if (needSyz) {
389: if (r ISZERO) { rp = CurrentRingp; } else { rp = r->m->ringp; }
390: cf = cxx(1,0,0,rp);
391: syz = ZERO;
392: }
393:
394: if (r != POLYNULL) {
395: rp = r->m->ringp;
396: if (! rp->weightedHomogenization) {
397: errorKan1("%s\n","ecart.c: the given ring must be declared with [(weightedHomogenization) 1]");
398: }
1.3 takayama 399: }
1.5 takayama 400:
401: r = goHomogenize11(r,DegreeShifto_vec,DegreeShifto_size,-1,1);
1.1 takayama 402: /* 1 means homogenize only s */
403:
404: do {
405: if (DebugReductionRed) printf("r=%s\n",POLYToString(r,'*',1));
406: ells = ecartFindReducer(r,gset,gg);
407: ell = ells.ell;
408: if (ell > 0) {
1.5 takayama 409: gg = ecartPutPolyInG(r,gg,POLYNULL,POLYNULL);
1.1 takayama 410: }
411: if (ell >= 0) {
412: if (ells.first) {
413: pp = ((gset->polys[ells.grade])->gh)[ells.gseti];
414: }else{
415: pp = (gg->pa)[ells.ggi];
416: }
417: if (ell > 0) r = mpMult(cxx(1,0,ell,rp),r); /* r = s^ell r */
418: r = (*reduction1)(r,pp,needSyz,&cc,&cg);
419: if (needSyz) {
420: if (ells.first) {
421: cf = ppMult(cc,cf);
422: syz = cpMult(toSyzCoeff(cc),syz);
423: syz = ppAddv(syz,toSyzPoly(cg,ells.grade,ells.gseti));
424: }else{
1.6 ! takayama 425: fprintf(stderr,"It has not yet implemented.\n");
! 426: exit(10);
1.1 takayama 427: /* BUG: not yet */
428: }
429: }
1.5 takayama 430: if (r ISZERO) goto ss1;
431: r = ecartDivideSv(r,&se); /* r = r/s^? */
1.1 takayama 432: }
433: }while (ell >= 0);
434:
1.5 takayama 435: ss1: ;
1.1 takayama 436: if (needSyz) {
437: syzp->cf = cf; /* cf is in the CurrentRingp */
438: syzp->syz = syz; /* syz is in the SyzRingp */
439: /* BUG: dehomogenize the syzygy */
1.6 ! takayama 440: fprintf(stderr,"It has not yet implemented.\n");
! 441: exit(10);
1.1 takayama 442: }
1.5 takayama 443:
1.1 takayama 444: r = goDeHomogenizeS(r);
1.5 takayama 445:
1.1 takayama 446: return(r);
447: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>