Annotation of OpenXM/src/kan96xx/Kan/gbGM.c, Revision 1.6
1.6 ! ohara 1: /* $OpenXM: OpenXM/src/kan96xx/Kan/gbGM.c,v 1.5 2005/06/16 06:54:55 takayama Exp $ */
1.1 maekawa 2: /* gbGM.c GM=Gebauer and Moller
1.3 takayama 3: */
1.1 maekawa 4:
5: #include <stdio.h>
1.6 ! ohara 6: #include <stdlib.h>
1.1 maekawa 7: #include "datatype.h"
8: #include "extern2.h"
9: #include "matrix.h"
10: #include "gradedset.h"
11:
12:
13: /*#define DEBUG*/
14: /* Take statistics of mmLarger_matrix in orderT.c */
15: /* #define STATISTICS */
16:
17: struct polySet_gm {
18: POLY *g;
19: int lim;
20: int size;
21: int *del;
22: };
23:
24: struct pair_gm {
25: POLY lcm;
26: int i;
27: int j;
28: };
29:
30: struct pairSet {
31: struct pair_gm *p;
32: int lim;
33: int size;
34: int *del;
35: };
36:
37: /* prototypes */
38: struct polySet_gm newPolySet_gm(int n);
39: struct pairSet newPairSet(int n);
40: int pairSetSize(struct pairSet d);
41: struct pairSet pairSetJoin(struct pairSet a,struct pairSet b);
42: struct polySet_gm enlargePolySet_gm(struct polySet_gm g);
43: struct pairSet deletePair_gm(struct pairSet d,int index);
44: int minPair_gm(struct pairSet d);
45: struct pairSet updatePair_gm(struct pairSet d,int t,
1.3 takayama 46: struct polySet_gm g,POLY gt);
1.1 maekawa 47: struct polySet_gm markRedundant_gm(struct polySet_gm g,int j);
48: /*
1.3 takayama 49: struct gradedPolySet *groebner_gm(struct arrayOfPOLY f,
50: int needBack,
51: int needSyz, struct pair **grP);
1.1 maekawa 52: */
53: void outputPairSet(struct pairSet d);
54: void outputPolySet_gm(struct polySet_gm g);
55: void errorGBGM(char *s);
56: POLY mReductionBySet(POLY f, POLY *s,int size);
57: POLY mReductionCdrBySet(POLY f, POLY *s,int size);
58: int criterion1(POLY f,POLY g,POLY lc);
59:
60: extern int UseCriterion1;
61: extern int Verbose;
62:
63: #define GLIMIT 100
64: int Criterion1,Criterion2; /* for taking a statistics. */
65: int Spairs;
66: extern int Criterion2B, Criterion2F, Criterion2M;
67: extern int Statistics;
68:
69: #ifdef STATISTICS
70: int CountE;
71: int CountI[2*N0];
72: #endif
73:
74: struct polySet_gm newPolySet_gm(n)
1.3 takayama 75: int n;
1.1 maekawa 76: {
77: struct polySet_gm g;
78: int i;
79: g.g = (POLY *) sGC_malloc(sizeof(POLY)*(n+1));
80: g.del = (int *) sGC_malloc(sizeof(int)*(n+1));
81: if ((g.g == (POLY *)NULL) ||
82: (g.del == (int *)NULL))
83: errorGBGM("No more memory.");
84: g.lim = n;
85: g.size = 0;
86: for (i=0; i<n; i++) g.del[i] = 0;
87: return(g);
88: }
89:
90: struct pairSet newPairSet(n)
1.3 takayama 91: int n;
1.1 maekawa 92: {
93: struct pairSet g;
94: int i;
95: g.p = (struct pair_gm *) sGC_malloc(sizeof(struct pair_gm)*(n+1));
96: g.del = (int *) sGC_malloc(sizeof(int)*(n+1));
97: if ((g.p == (struct pair_gm *)NULL) ||
98: (g.del == (int *)NULL))
99: errorGBGM("No more memory.");
100: g.lim = n;
101: g.size = 0;
102: for (i=0; i<n; i++) g.del[i] = 0;
103: return(g);
104: }
105:
106: int pairSetSize(d)
1.3 takayama 107: struct pairSet d;
1.1 maekawa 108: {
109: int s,i;
110: s = 0;
111: for (i=0; i< d.size; i++) {
112: if (d.del[i] == 0) s++;
113: }
114: return(s);
115: }
116:
117: struct pairSet pairSetJoin(a,b)
1.3 takayama 118: struct pairSet a,b;
1.1 maekawa 119: {
120: int m,n,k,i;
121: struct pairSet ans;
122: m = pairSetSize(a);
123: n = pairSetSize(b);
124: ans = newPairSet(m+n);
125: k = 0;
126: for (i=0 ; i<a.size; i++) {
127: if (a.del[i] == 0) {
128: ans.p[k] = a.p[i]; ans.del[k++] = 0;
129: }
130: }
131: for (i=0; i<b.size; i++) {
132: if (b.del[i] == 0) {
133: ans.p[k] = b.p[i]; ans.del[k++] = 0;
134: }
135: }
136: if (k != (m+n)) errorGBGM("Internal error in pairSetJoin().");
137: ans.size = m+n;
138: return(ans);
139: }
140:
141: struct polySet_gm enlargePolySet_gm(g)
1.3 takayama 142: struct polySet_gm g;
1.1 maekawa 143: {
144: int i;
145: struct polySet_gm ans;
146: if (g.lim < g.size+3) {
147: ans = newPolySet_gm((g.lim)*2);
148: for (i=0; i<g.size; i++) {
149: ans.g[i] = g.g[i];
150: ans.del[i] = g.del[i];
151: }
152: ans.size = g.size;
153: }else{ ans = g; }
154: return(ans);
155: }
156:
157: struct pairSet deletePair_gm(d,index)
1.3 takayama 158: struct pairSet d;
159: int index;
160: /* delete d[index] */
1.1 maekawa 161: {
162: int i;
163: d.del[index] = 0;
164: for (i= index; i<d.size-1; i++) {
165: d.p[i] = d.p[i+1];
166: }
167: d.size--;
168: return(d);
169: }
170:
171: int minPair_gm(d)
1.3 takayama 172: struct pairSet d;
1.1 maekawa 173: {
174: POLY min;
175: int index,i;
176: min = d.p[0].lcm;
177: index = 0;
178: for (i=0; i < d.size; i++) {
179: if ( (*mmLarger)(min, d.p[i].lcm)) {
180: index = i;
181: min = d.p[i].lcm;
182: }
183: }
184: return(index);
185: }
186:
187: struct pairSet updatePair_gm(d,t,g,gt)
1.3 takayama 188: struct pairSet d;
189: int t;
190: struct polySet_gm g;
191: POLY gt;
1.1 maekawa 192: {
193: int i,j,k;
194: struct pairSet new;
195: POLY it;
196: POLY jt;
197:
198: if (Verbose >= 2) {
199: printf("\nupdatepair_gm(): ");
200: printf("g.size=%d ",g.size);
201: printf("d.size=%d; ",d.size);
202: }
203:
204: /* make list (1,t), (2,t), .... ====> new ( D1 ) */
205: new = newPairSet(g.size); new.size = 0;
206: for (i=0; i<g.size; i++) {
207: if (g.del[i] == 0) {
208: new.p[i].i = i;
209: new.p[i].j = t;
210: new.p[i].lcm = (*lcm)(g.g[i],gt);
211: }else{
212: new.p[i].i = i;
213: new.p[i].j = t;
214: new.p[i].lcm = ZERO; /* Compute it when it is necessary. */
215: }
216: new.del[i] = g.del[i];
217: new.size++;
218: }
219: #ifdef DEBUG
220: /*printf("\nnew is ...");
1.3 takayama 221: outputPairSet(new);*/
1.1 maekawa 222: #endif
223:
224: /* Candel in D=d all (i,j) such that B_t(i,j) */
225:
226: for (k=0; k<d.size; k++) {
227: if ((d.del[k] == 0) && ((*isReducible)(d.p[k].lcm,gt)) ) {
228: /* check T(it) != T(i,j) != T(j,t) */
229: i = d.p[k].i; j = d.p[k].j;
230: if ((new.del[i] == 1) || (new.del[j] == 1)) {
1.3 takayama 231: /* fprintf(stderr,"Warning in updatepair_gm(): i=%d, j=%d; rewriting new.\n",i,j); */
232: if (new.del[i] == 1) {
233: new.p[i].lcm = (*lcm)(g.g[i],gt);
234: }
235: if (new.del[j] == 1) {
236: new.p[j].lcm = (*lcm)(g.g[j],gt);
237: }
1.1 maekawa 238: }
239: if (((*mmLarger)((new.p[i].lcm),(new.p[j].lcm)) != 2) &&
1.3 takayama 240: ((*mmLarger)((new.p[i].lcm),(d.p[k].lcm)) != 2) &&
241: ((*mmLarger)((new.p[j].lcm),(d.p[k].lcm)) != 2)) {
242: /* delete T(i,j) in d */
243: d.del[k] = 1;
244: if (Verbose >= 2) printf("B%d(%d,%d) ",t,i,j);
245: Criterion2B++;
1.1 maekawa 246: }
247: }
248: }
249:
250: /* Then, d is D' */
251: /* Next, check M(i,t) and F(i,t) */
252: for (i=0; i<g.size; i++) {
253: for (j=i+1; j<g.size; j++) {
254: if ((g.del[i] == 0) && (g.del[j] == 0)) {
1.3 takayama 255: it = new.p[i].lcm;
256: jt = new.p[j].lcm;
257: switch ( (*mmLarger)(it,jt)) {
258: case 2:
259: /* F(j,t), i<j */
260: new.del[j] = 1;
261: if (Verbose >= 2) printf("F%d(%d,%d) ",i,j,t);
262: Criterion2F++;
263: break;
264: case 1:
265: /* g[i] > g[j], M(i,t) */
266: if ((*isReducible)(it,jt)) {
267: new.del[i] = 1;
268: if (Verbose >=2) printf("M%d(%d,%d) ",j,i,t);
269: Criterion2M++;
270: }
271: break;
272: case 0: /* M(j,t) */
273: if ((*isReducible)(jt,it)) {
274: new.del[j] = 1;
275: if (Verbose >=2) printf("M%d(%d,%d) ",i,j,t);
276: Criterion2M++;
277: }
278: break;
1.1 maekawa 279:
1.3 takayama 280: }
1.1 maekawa 281: }
282: }
283: }
284:
285: /* Finally, check T(i) T(t) = T(i,t) */
286: if (UseCriterion1) {
287: for (i=0; i<g.size; i++) {
288: if (new.del[i] == 0) {
1.3 takayama 289: if (criterion1(g.g[i],gt,new.p[i].lcm)) {
290: new.del[i] = 1;
291: if (Verbose >=2) printf("1(%d,%d) ",i,t);
292: Criterion1++;
293: }
1.1 maekawa 294: }
295: }
296: }
297:
298: /* d <-- D' \cup new */
299: d = pairSetJoin(d,new);
300: if (Verbose >=2) printf("new d.size = %d.\n",d.size);
301: #ifdef DEBUG
302: outputPairSet(d);
303: #endif
304: return(d);
305: }
306:
307:
308: struct polySet_gm markRedundant_gm(g,j)
1.3 takayama 309: struct polySet_gm g;
310: int j;
311: /* compare only with g[j] */
1.1 maekawa 312: {
313: int i;
314: for (i=0; i<g.size; i++) {
315: if ((g.del[i] == 0) && (i != j)) {
316: if (g.del[j] == 0) {
1.3 takayama 317: switch((*mmLarger)((g.g[i]),(g.g[j]))) {
318: case 2:
319: g.del[j] = 1;
320: break;
321: case 1:
322: if ((*isReducible)((g.g[i]),(g.g[j]))) {
323: g.del[i] = 1;
324: }
325: break;
326: case 0:
327: if ((*isReducible)((g.g[j]),(g.g[i]))) {
328: g.del[j] = 1;
329: }
330: break;
331: }
1.1 maekawa 332: }
333: }
334: }
335: #ifdef DEBUG
336: /*printf("\nend of markRedundant_gm...");
1.3 takayama 337: outputPolySet_gm(g);*/
1.1 maekawa 338: #endif
339: return(enlargePolySet_gm(g));
340: }
341:
342:
343:
1.5 takayama 344: struct gradedPolySet *groebner_gm(f,needBack,needSyz,grP,countDown,forceReduction,reduceOnly,gbCheck)
1.3 takayama 345: struct arrayOfPOLY *f;
346: int needBack;
347: int needSyz;
348: struct pair **grP;
349: int countDown;
350: int forceReduction;
1.4 takayama 351: int reduceOnly;
1.5 takayama 352: int gbCheck;
1.1 maekawa 353: {
354: int r;
355: struct pair_gm top;
356: int tindex;
357: POLY h;
358: int i,j;
359: int ss;
360: struct polySet_gm g;
361: struct pairSet d;
362: int t;
363: struct spValue sv;
364: struct gradedPolySet *ans;
365: int grade,indx;
366: extern int ReduceLowerTerms;
367:
368: if (needBack || needSyz) {
369: fprintf(stderr,"Warning: groebner_gm() does not compute the backward transformation and syzygies.\n");
1.4 takayama 370: }
1.5 takayama 371: if (reduceOnly || gbCheck) {
372: fprintf(stderr,"Warning: groebner_gm() does not implement reduceOnly/gbCheck.\n");
1.1 maekawa 373: }
374:
375: #ifdef STATISTICS
376: CountE = 0;
377: for (i=0; i<2*N; i++) CountI[i]=0;
378: #endif
379: /* taking the statistics. */
380: Criterion1 = Criterion2B = Criterion2M = Criterion2F = 0;
381: Spairs = 0;
382:
383: r = f->n;
384: g = newPolySet_gm(r + GLIMIT);
385: d = newPairSet(1);
386:
387: g.size = 0;
388: g.g[g.size] = getArrayOfPOLY(f,g.size);
389: g.del[g.size] = 0; g.size ++;
390:
391: for (t=1; t<r; t++) {
392: d = updatePair_gm(d,t,g,getArrayOfPOLY(f,t));
393: g.g[t] = getArrayOfPOLY(f,t); g.del[t] = 0; g.size ++ ;
394: g = markRedundant_gm(g,t);
395: }
396:
397: #ifdef DEBUG
398: outputPolySet_gm(g);
399: #endif
400:
401: while (d.size > 0) {
402: tindex = minPair_gm(d);
403: top = d.p[tindex];
404: if (Verbose) {
405: printf("\nRemaining pairs = %d : ",d.size);
406: }
407: i = top.i;
408: j = top.j;
409: if (Verbose >=2) printf("sp(%d,%d) ",i,j);
410: Spairs++;
411: sv = (*sp)(g.g[i],g.g[j]);
412: h = ppAddv(ppMult(sv.a,g.g[i]),ppMult(sv.b,g.g[j]));
413: h = mReductionBySet(h,g.g,g.size);
414: if (Verbose >=2) printf("->%s ",POLYToString(h,' ',0));
415:
416: if (!Verbose) {
417: printf(" <%d>",d.size);
418: fflush(stdout);
419: }
420:
421: if (!(h ISZERO)) {
422: if (ReduceLowerTerms) {
1.3 takayama 423: h = mReductionCdrBySet(h,g.g,g.size);
1.1 maekawa 424: }
425: d = updatePair_gm(d,r,g,h);
426: g.g[r] = h; g.del[r] = 0;
427: r++; g.size++;
428: g = markRedundant_gm(g,r-1);
429: }
430: d = deletePair_gm(d,tindex);
431: #ifdef DEBUG
432: outputPolySet_gm(g);
433: #endif
434: }
435:
436:
437: if (Verbose || Statistics) {
438: printf("\n\nThe number of s-pairs is %d.\n",Spairs);
439: printf("Criterion1 is applied %d times.\n",Criterion1);
440: printf("Criterions M,F and B are applied M=%d, F=%d, B=%d times.\n",Criterion2M,Criterion2F,Criterion2B);
441: Criterion1 = Criterion2M = Criterion2F = Criterion2B = 0;
442: Spairs = 0;
443: }
444: #ifdef STATISTICS
445: printf("\n\nEqual : %d\n",CountE);
446: for (i=0; i<2*N; i++) {
447: printf("%d times of loop : %d\n",i,CountI[i]);
448: }
449: #endif
450:
451: ans = newGradedPolySet(1);
1.5 takayama 452: ans->gb = 1;
1.1 maekawa 453: for (i=0; i<ans->lim; i++) {
454: ans->polys[i] = newPolySet(1);
455: }
456:
457: for (i=0; i<g.size; i++) {
458: if (g.del[i] == 0) {
459: grade=-1;whereInG(ans,g.g[i],&grade,&indx,0); /* we do not use sugar. */
460: ans = putPolyInG(ans,g.g[i],grade,indx,(struct syz0 *)NULL,0,-1);
461: }
462: }
463:
464: printf("\n");
465: return(ans);
466: }
467:
468:
469: void outputPairSet(d)
1.3 takayama 470: struct pairSet d;
1.1 maekawa 471: { int i;
1.3 takayama 472: printf("\nOutput struct pairSet. ");
473: printf(".size=%d, .lim=%d :\n",d.size,d.lim);
474: for (i=0; i<d.size; i++) {
475: printf("%d[%d] i=%d,j=%d, lcm=%s ",i,d.del[i],d.p[i].i,d.p[i].j,
476: POLYToString(d.p[i].lcm,' ',0));
477: }
478: printf("\n"); fflush(stdout);
1.1 maekawa 479: }
480:
481: void outputPolySet_gm(g)
1.3 takayama 482: struct polySet_gm g;
1.1 maekawa 483: { int i;
1.3 takayama 484: printf("\nOutput struct polySet_gm. ");
485: printf(".size=%d, .lim=%d :\n",g.size,g.lim);
486: for (i=0; i<g.size; i++) {
487: printf("%d[%d] %s\n",i,g.del[i],POLYToString(g.g[i],' ',0));
488: }
489: printf("\n"); fflush(stdout);
1.1 maekawa 490: }
491:
492: int criterion1(f0,g0,lc0)
1.3 takayama 493: POLY f0;
494: POLY g0;
495: POLY lc0;
1.1 maekawa 496: {
497: /* This is used only for commutative case. */
498: register int i;
499: MONOMIAL f,g,lc;
500: int n;
501: f = f0->m; g = g0->m; lc = lc0->m;
502: if (f->ringp != g->ringp) return(0);
503: if (f->ringp != lc->ringp) return(0);
504: n = f->ringp->n;
505: for (i = 0; i<n ; i++) {
506: if (( f->e[i].x + g->e[i].x) != (lc->e[i].x)) return(0);
507: }
508: for (i = 0; i<n ; i++) {
509: if (( f->e[i].D + g->e[i].D) != (lc->e[i].D)) return(0);
510: }
511: return(1);
512: }
513:
514:
515: POLY mReductionBySet(f,s,size)
1.3 takayama 516: POLY f;
517: POLY *s;
518: int size;
1.1 maekawa 519: {
520: int reduced1;
521: int i;
522: POLY cc,cg;
523: do {
524: reduced1 = 0;
525: for (i=0; i<size; i++) {
526: if (f ISZERO) goto ss;
527: if ((*isReducible)(f,s[i])) {
1.3 takayama 528: f = (*reduction1)(f,s[i],0,&cc,&cg);
529: reduced1 = 1;
1.1 maekawa 530: }
531: }
532: } while (reduced1 != 0);
533: ss: ;
534: return(f);
535: }
536:
537: POLY mReductionCdrBySet(f,s,size)
1.3 takayama 538: POLY f;
539: POLY *s;
540: int size;
1.1 maekawa 541: {
542: int reduced1;
543: int i;
544: POLY cc,cg;
545: POLY fs;
546: do {
547: reduced1 = 0;
548: for (i=0; i<size; i++) {
549: if (f ISZERO) goto ss;
550: if ((fs=(*isCdrReducible)(f,s[i])) != ZERO) {
1.3 takayama 551: f = (*reduction1Cdr)(f,fs,s[i],0,&cc,&cg);
552: reduced1 = 1;
1.1 maekawa 553: }
554: }
555: } while (reduced1 != 0);
556: ss: ;
557: return(f);
558: }
559:
560:
561: void errorGBGM(s)
1.3 takayama 562: char *s;
1.1 maekawa 563: {
564: fprintf(stderr,"Fatal error in yaGbasis.c: %s \n",s);
565: exit(10);
566: }
567:
568:
569:
570:
571:
572:
573:
574:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to gbGM.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM / src / kan96xx / Kan