Annotation of OpenXM/src/kan96xx/Kan/gb.c, Revision 1.1.1.1
1.1 maekawa 1: #include <stdio.h>
2: #include "datatype.h"
3: #include "extern2.h"
4: #include "matrix.h"
5: #include "gradedset.h"
6:
7: #define INITGRADE 4
8: #define INITSIZE 2
9:
10: int KanGBmessage = 1;
11:
12: static int Debug = 0;
13: static int sugarGrade(struct pair *inode,struct gradedPolySet *grG,
14: POLY gt,int gtSugarGrade);
15: extern int Sugar;
16: extern int Statistics;
17: extern int Criterion1;
18: extern int UseCriterion1;
19: extern int UseCriterion2B;
20: extern int Spairs;
21: extern int Criterion2B, Criterion2F, Criterion2M;
22:
23:
24: struct gradedPairs *updatePairs(grD,gt,gtGrade,t,grG)
25: struct gradedPairs *grD; /* set of pairs */
26: POLY gt; /* new polynomial */
27: int gtGrade;
28: int t;
29: struct gradedPolySet *grG; /* (f,gt), f \in grG, should be added to grD. */
30: {
31: int gmax,newGrade;
32: struct pair *node,*new,*inode,*jnode;
33: int i,k;
34: struct polySet *g;
35: POLY lcmp;
36: POLY it,jt;
37: extern int Verbose;
38:
39:
40: if (Verbose) {
41: printf("\nupdatePairs(): ");
42: }
43:
44: gmax = grG->maxGrade;
45: node = new = newPair((struct pair *)NULL);
46: /* (f,gt) are stored in new as a list. */
47: for (k=0; k<gmax; k++) {
48: g = grG->polys[k];
49: for (i=0; i<g->size; i++) {
50: if (g->del[i] == 0 && (gtGrade != k || t != i)) {
51: lcmp = (*lcm)(g->g[i],gt);
52: if (lcmp == ZERO) {
53:
54: }else {
55: new->del = 0;
56: new->next = newPair(new);
57: new = new->next;
58: new->lcm = lcmp;
59: new->ig = k; new->ii = i; /* g->g[i] */
60: new->jg = gtGrade; new->ji = t; /* gt */
61: new->grade = -1; /* one do not need to set grade here. */
62: }
63: }
64: }
65: }
66:
67:
68: /* Check the criterion */
69: /* See Gebauer and Mora, 1988, Journal of Symbolic Computation
70: 279. We must not use the test T(i)T(j) = T(i,j) because it
71: does not hold in the ring of diff op.
72: */
73:
74: inode = node->next;
75: while (inode != (struct pair *)NULL) {
76: jnode = inode->next;
77: while (jnode != (struct pair *)NULL) {
78: it = inode->lcm;
79: jt = jnode->lcm;
80: if ((*mmLarger)(it,jt) == 2) {
81: /* F(j,t), i<j */
82: jnode->del = 1;
83: if (Verbose) printf("F[%d,%d]([%d,%d],[%d,%d]) ",
84: inode->ig,inode->ii,
85: jnode->ig,jnode->ii,
86: gtGrade,t);
87: Criterion2F++;
88: }else {
89: /* g[i] > g[j] ?, M(i,t) */
90: if ((*isReducible)(it,jt)) {
91: inode->del = 1;
92: if (Verbose) printf("M[%d,%d]([%d,%d],[%d,%d]) ",
93: jnode->ig,jnode->ii,
94: inode->ig,inode->ii,
95: gtGrade,t);
96: Criterion2M++;
97: }
98: /* M(j,t) */
99: if ((*isReducible)(jt,it)) {
100: jnode->del = 1;
101: if (Verbose) printf("M[%d,%d]([%d,%d],[%d,%d]) ",
102: inode->ig,inode->ii,
103: jnode->ig,jnode->ii,
104: gtGrade,t);
105: Criterion2M++;
106: }
107: }
108: jnode = jnode->next;
109: }
110: inode = inode->next;
111: }
112:
113: if (UseCriterion1) {
114: inode = node->next;
115: while (inode != NULL) {
116: if (inode->del == 0) {
117: if (criterion1(gt,grG->polys[inode->ig]->g[inode->ii],inode->lcm)) {
118: inode->del = 1;
119: if (Verbose) printf("Criterion1([%d,%d],[%d,%d]) ",
120: inode->ig,inode->ii,
121: gtGrade,t);
122: Criterion1++;
123: }
124: }
125: inode = inode->next;
126: }
127: }
128:
129: if (UseCriterion2B) {
130: Criterion2B += deletePairByCriterion2B(grD,gt,grG);
131: }
132:
133: /* Merge to grD */
134: inode = node->next;
135: if (Debug) outputNode(inode);
136: while (inode != (struct pair *)NULL) {
137: if (inode->del == 0) {
138: if (Sugar) {
139: inode->grade = sugarGrade(inode,grG,gt,gtGrade);
140: }else{
141: inode->grade = (*grade)(inode->lcm);
142: }
143: if (grD->lim <= inode->grade) {
144: grD = enlargeGradedPairs(2*(inode->grade)+1,grD);
145: }
146: insertPair(inode,grD->pairs[inode->grade]);
147: grD->maxGrade = max(grD->maxGrade,inode->grade+1); /* don't forget */
148: }
149: inode = inode->next;
150: }
151: if (Debug) printf("OK.\n");
152: if (Verbose) {
153: printf(" Remaining pairs are %d. maxGrade=%d\n",countPairs(grD),grD->maxGrade);
154: }
155: return(grD);
156: }
157:
158: struct gradedPolySet *groebner_gen(f,needBack,needSyz,grP,countDown,forceReduction)
159: struct arrayOfPOLY *f;
160: int needBack;
161: int needSyz;
162: struct pair **grP; /* if (needSyz), it is set. */
163: int countDown;
164: int forceReduction;
165: {
166: int r;
167: struct gradedPolySet *g;
168: struct gradedPairs *d;
169: int i;
170: int grade,indx;
171: POLY gt;
172: struct pair *top;
173: int ig,ii,jg,ji;
174: POLY gi,gj;
175: struct spValue h;
176: struct syz0 syz;
177: int pgrade = 0;
178: POLY rd;
179: POLY syzPoly;
180: POLY syzCf;
181: struct syz0 *syzp;
182: int serial;
183: struct pair *listP;
184: extern int Verbose;
185: extern int ReduceLowerTerms;
186: extern int StopDegree;
187: extern struct ring *CurrentRingp;
188: extern char *F_mpMult;
189: struct ring *rp;
190:
191: r = f->n;
192: if (r<=0) return((struct gradedPolySet *)NULL);
193: if (UseCriterion1) {
194: if (needSyz) {
195: warningGradedSet("You cannot use both the options needSyz and UseCriterion1.\nInput [(UseCriterion1) 0] system_variable to turn off the use of the Criterion 1.");
196: }
197: if (strcmp(F_mpMult,"poly")!=0) {
198: warningGradedSet("You can use the option UseCriterion1 only for the computation in the ring of polynomials.\nInput [(UseCriterion1) 0] system_variable to turn off the use of the Criterion 1.");
199: }
200: }
201:
202: Spairs = Criterion1 = Criterion2B = Criterion2F = Criterion2M = 0;
203:
204: g = newGradedPolySet(INITGRADE);
205: d = newGradedPairs(INITGRADE*2);
206: for (i=0; i<g->lim; i++) {
207: g->polys[i] = newPolySet(INITSIZE);
208: }
209:
210: for (i=0; i<r; i++) {
211: gt = getArrayOfPOLY(f,i);
212: if (gt ISZERO) { rp = CurrentRingp; } else { rp = gt->m->ringp; }
213: grade = -1; whereInG(g,gt,&grade,&indx,Sugar);
214: if (KanGBmessage == 2) {
215: printf("init=%s, ",POLYToString(head(gt),'*',1));
216: printf("(gr,indx)=(%d,%d).\n",grade,indx);
217: }
218: d = updatePairs(d,gt,grade,indx,g);
219: serial = i;
220: if (!needBack) {
221: g = putPolyInG(g,gt,grade,indx,
222: (struct syz0 *)NULL,1,serial);
223: }else {
224: syzp = newSyz0();
225: syzp->cf = cxx(1,0,0,rp);
226: syzp->syz = toSyzPoly(cxx(1,0,0,rp),grade,indx);
227: g = putPolyInG(g,gt,grade,indx,syzp,1,serial);
228: }
229:
230: markRedundant0(g,grade,indx);
231: if (Debug) {
232: outputGradedPairs(d); outputGradedPolySet(g,needSyz);
233: }
234: }
235: if (needSyz) {
236: *grP = newPair((struct pair *)NULL);
237: listP = *grP;
238: }
239:
240: while ((top = getPair(d)) != (struct pair *)NULL) {
241: ig = top->ig; ii = top->ii; /* [ig,ii] */
242: jg = top->jg; ji = top->ji; /* [jg,ji] */
243: gi = g->polys[ig]->g[ii];
244: gj = g->polys[jg]->g[ji];
245:
246: Spairs++;
247: h = (*sp)(gi,gj);
248: rd = ppAddv(ppMult(h.a,gi),ppMult(h.b,gj));
249: if (!Sugar || forceReduction) {
250: rd = (*reduction)(rd,g,needBack,&syz);
251: }else{
252: rd = reduction_sugar(rd,g,needBack,&syz,top->grade);
253: }
254: syzPoly = syz.syz;
255: syzCf = syz.cf;
256:
257: if (KanGBmessage) {
258: if (KanGBmessage == 2) {
259: printf("sp([%d,%d],[%d,%d]), ",ig,ii,jg,ji);
260: if (rd ISZERO) {
261: printf(" 0 \n");
262: } else{
263: printf(" terms=%d, grade=%d, ",pLength(rd),top->grade);
264: printf(" init=%s, ",POLYToString(head(rd),'*',1));
265: }
266: }else{
267: if (pgrade != top->grade) {
268: pgrade = top->grade;
269: printf(" %d",pgrade);
270: fflush(stdout);
271: }
272: if (rd ISZERO) {
273: printf("o"); fflush(stdout);
274: }else{
275: printf("."); fflush(stdout);
276: }
277: }
278: }
279:
280: if (!(rd ISZERO)) {
281: if (needBack || needSyz) {
282: syzp = newSyz0();
283: syzp->cf = syzCf; /* no meaning */
284: syzp->syz = ppAdd(toSyzPoly(h.a,ig,ii),toSyzPoly(h.b,jg,ji));
285: syzp->syz = cpMult(toSyzCoeff(syzCf),syzp->syz);
286: syzp->syz = ppAdd(syzp->syz,syzPoly);
287: }
288:
289: if (ReduceLowerTerms && !(Sugar)) {
290: rd = (*reductionCdr)(rd,g,needBack,&syz);
291: if (needBack || needSyz) {
292: /* syzp->cf = ppMult(syz.cf,syzp->cf); no meaning */
293: syzp->syz = ppAdd(syz.syz,
294: cpMult(toSyzCoeff(syz.cf),syzp->syz));
295: }
296: }
297:
298: if(Sugar && (!forceReduction)){grade=top->grade;}else{grade=-1;}whereInG(g,rd,&grade,&indx,Sugar);
299: if (KanGBmessage == 2) {
300: printf("(gr,indx)=(%d,%d).\n",grade,indx);
301: /*
302: printf("sp(%s,%s)-->%s\n",POLYToString(gi,' ',1),
303: POLYToString(gj,' ',1),
304: POLYToString(rd,' ',1));
305: */
306: }
307:
308: d = updatePairs(d,rd,grade,indx,g);
309: g = putPolyInG(g,rd,grade,indx,syzp,0,-1);
310: if (Sugar) { markRedundant0(g,grade,indx);}
311: else {markRedundant(g,rd,grade,indx,Sugar);}
312:
313: if (KanGBmessage && (StopDegree < pgrade)) {
314: printf("Computation of the Groebner basis is suspended bacause of StopDegree < computing grade.\n");
315: printf("Note that the result is NOT groebner basis.\n");
316: break;
317: }
318: if (countDown) {
319: if (eliminated(rd) == 1) {
320: --countDown;
321: printf("x"); fflush(stdout);
322: if (countDown == 0) {
323: printf("\nThe computation of the Groebner basis is suspended because of countDown==0.\n");
324: printf("Note that the result is NOT groebner basis.\n");
325: break;
326: }
327: }
328: }
329: if (Debug) {
330: outputGradedPairs(d); outputGradedPolySet(g,needSyz);
331: }
332: }else{
333: if (needSyz) {
334: top->syz = ppAdd(toSyzPoly(h.a,ig,ii),toSyzPoly(h.b,jg,ji));
335: top->syz = cpMult(toSyzCoeff(syzCf),top->syz);
336: top->syz = ppAdd(top->syz,syzPoly);
337: listP->next = top; top->prev = listP; listP = listP->next;
338: }
339: }
340: }
341:
342: if (KanGBmessage == 2) {
343: outputGradedPolySet(g,needSyz);
344: }
345: if (KanGBmessage) {
346: if (Sugar) {
347: printf("\nCompleted (GB with sugar).\n");
348: }else{
349: printf("\nCompleted.\n");
350: }
351: }
352:
353: if (Statistics) {
354: printf("\nThe number of s-pairs is %d.\n",Spairs);
355: printf("Criterion1 is applied %d times.\n",Criterion1);
356: printf("Criterions M,F and B are applied M=%d, F=%d, B=%d times.\n",Criterion2M,Criterion2F,Criterion2B);
357: Spairs = Criterion1 = Criterion2M = Criterion2F = Criterion2B = 0;
358: }
359:
360:
361: return(g);
362: }
363:
364:
365: static int sugarGrade(struct pair *inode,struct gradedPolySet *grG,
366: POLY gt,int gtSugarGrade)
367: {
368: int a,b,ans;
369: int debug = 0;
370: a = grade_gen(inode->lcm)-grade_gen(grG->polys[inode->ig]->g[inode->ii]);
371: b = grade_gen(inode->lcm)-grade_gen(gt);
372: /* inode = lcm(f_i, gt) = p f_i = q gt modulo lower order terms.
373: a = tdeg(p), b = tdeg(gt);
374: */
375: if (debug) {
376: printf("Sugar grade of sp([%d,%d],[%d,%d]) ",inode->ig,inode->ii,
377: inode->jg,inode->ji);
378: printf("is max(%d+%d,%d+%d)\n",a,inode->ig,b,gtSugarGrade);
379: }
380: a = a+(inode->ig); /* new sugar degree of p f_i. */
381: b = b + gtSugarGrade; /* new sugar degree of q gt */
382: return( a > b ? a : b);
383: }
384:
385:
386:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to gb.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM / src / kan96xx / Kan