Annotation of OpenXM/src/kan96xx/Kan/gb.c, Revision 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>