[BACK]Return to gb.c CVS log [TXT][DIR] Up to [local] / OpenXM / src / kan96xx / Kan

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>