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

Annotation of OpenXM/src/kan96xx/Kan/ecart.c, Revision 1.2

1.2     ! takayama    1: /* $OpenXM: OpenXM/src/kan96xx/Kan/ecart.c,v 1.1 2003/07/17 09:10:54 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;
                     19: };
                     20:
                     21: static struct ecartReducer ecartFindReducer(POLY r,struct gradedPolySet *gset,struct ecartPolyArray *epa);
                     22: static int ecartCheckPoly(POLY f);  /* check if it does not contain s-variables */
                     23: static int ecartCheckEnv();         /* check if the environment is OK for ecart div*/
                     24: static struct ecartPolyArray *ecartPutPolyInG(POLY g,struct ecartPolyArray *eparray);
                     25: static int ecartGetEll(POLY r,POLY g);
                     26: static POLY ecartDivideSv(POLY r);
                     27:
                     28: extern int DebugReductionRed;
                     29:
                     30: /* This is used for goHomogenization */
                     31: extern int DegreeShifto_size;
                     32: extern int *DegreeShifto_vec;
                     33:
                     34: #define LARGE 0x7fffffff
                     35:
                     36:
                     37: static POLY ecartDivideSv(POLY r) {
                     38:   POLY f;
                     39:   int k;
                     40:   if (r == POLYNULL) return r;
                     41:   f = r; k = LARGE;
                     42:   while (r != POLYNULL) {
                     43:     if (r->m->e[0].x < k) {
                     44:       k = r->m->e[0].x;
                     45:     }
                     46:     r = r->next;
                     47:   }
                     48:   if (k > 0) {
                     49:     f = r;
                     50:     while (r != POLYNULL) {
                     51:       r->m->e[0].x -= k;
                     52:       r = r->next;
                     53:     }
                     54:   }
                     55:   return f;
                     56: }
                     57:
                     58: static int ecartGetEll(POLY f,POLY g) {
                     59:   int n,i,p;
                     60:   MONOMIAL tf;
                     61:   MONOMIAL tg;
                     62:
                     63:   if (f ISZERO) return(-1);
                     64:   if (g ISZERO) return(-1);
                     65:
                     66:   checkRingIsR(f,g);
                     67:
                     68:   if (!isSameComponent_x(f,g)) return(-1);
                     69:   tf = f->m; tg = g->m; n = tf->ringp->n;
                     70:   for (i=1; i<n; i++) {
                     71:     if (tf->e[i].x < tg->e[i].x) return(-1);
                     72:     if (tf->e[i].D < tg->e[i].D) return(-1);
                     73:   }
                     74:   if (tf->e[0].D < tg->e[0].D) return(-1);  /*  h  */
1.2     ! takayama   75:   p = tf->e[0].x - tg->e[0].x;  /* H,  s */
1.1       takayama   76:   if (p >=0 ) return 0;
                     77:   else return(-p);
                     78: }
                     79:
                     80:
                     81: #define EP_SIZE 10
                     82: static struct ecartPolyArray *ecartPutPolyInG(POLY g,struct ecartPolyArray *eparray)
                     83: {
                     84:   struct ecartPolyArray *a;
                     85:   POLY *ep;
                     86:   int i;
                     87:   if (eparray == (struct ecartPolyArray *)NULL) {
                     88:     a = (struct ecartPolyArray *) sGC_malloc(sizeof(struct ecartPolyArray));
                     89:     outofmemory(a);
                     90:     ep = (POLY *) sGC_malloc(sizeof(POLY)*EP_SIZE);
                     91:     outofmemory(ep);
                     92:     a->limit = EP_SIZE;
                     93:     a->size = 0;
                     94:     a->pa = ep;
                     95:     eparray = a;
                     96:   }
                     97:   if (eparray->size >= eparray->limit) {
                     98:     a = (struct ecartPolyArray *) sGC_malloc(sizeof(struct ecartPolyArray));
                     99:     outofmemory(a);
                    100:     ep = (POLY *) sGC_malloc(sizeof(POLY)*((eparray->limit)*2));
                    101:     outofmemory(ep);
                    102:     a->limit = (eparray->limit)*2;
                    103:     a->size = eparray->size;
                    104:     a->pa = ep;
                    105:     for (i=0; i<eparray->size; i++) {
                    106:       (a->pa)[i] = (eparray->pa)[i];
                    107:     }
                    108:     eparray  = a;
                    109:   }
                    110:   i = eparray->size;
                    111:   (eparray->pa)[i] = g;
                    112:   (eparray->size)++;
                    113:   return eparray;
                    114: }
                    115:
                    116: static struct ecartReducer ecartFindReducer(POLY r,struct gradedPolySet *gset,
                    117:                                             struct ecartPolyArray *epa)
                    118: {
                    119:   int grd;
                    120:   struct polySet *set;
                    121:   int minGrade = 0;
                    122:   int minGseti = 0;
                    123:   int minGgi   = 0;
                    124:   int ell1 = LARGE;
                    125:   int ell2 = LARGE;
                    126:   int ell;
                    127:   int i;
                    128:   struct ecartReducer er;
                    129:   /* Try to find a reducer in gset; */
                    130:   grd = 0;
                    131:   while (grd < gset->maxGrade) {
                    132:     set = gset->polys[grd];
                    133:     for (i=0; i<set->size; i++) {
1.2     ! takayama  134:       if (set->gh[i] == POLYNULL) {
        !           135:         /* goHomogenize set->gh[i] */
        !           136:         set->gh[i] = goHomogenize11(set->g[i],DegreeShifto_vec,DegreeShifto_size,-1,1);
        !           137:       }
        !           138:       ell = ecartGetEll(r,set->gh[i]);
1.1       takayama  139:       if ((ell>=0) && (ell < ell1)) {
                    140:         ell1 = ell;
                    141:         minGrade = grd; minGseti=i;
                    142:       }
                    143:     }
                    144:     grd++;
                    145:   }
                    146:   if (epa != NULL) {
                    147:     /* Try to find in the second group. */
                    148:     for (i=0; i< epa->size; i++) {
                    149:       ell = ecartGetEll(r,(epa->pa)[i]);
                    150:       if ((ell>=0) && (ell < ell2)) {
                    151:         ell2 = ell;
                    152:         minGgi = i;
                    153:       }
                    154:     }
                    155:   }
                    156:
                    157:   if (DebugReductionRed) {
                    158:     printf("ecartFindReducer(): ell1=%d, ell2=%d, minGrade=%d, minGseti=%d, minGgi=%d\n",ell1,ell2,minGrade,minGseti,minGgi);
                    159:   }
                    160:   if (ell1 <= ell2) {
                    161:     if (ell1 == LARGE) {
                    162:       er.ell = -1;
                    163:       return er;
                    164:     }else{
                    165:       er.ell = ell1;
                    166:       er.first = 1;
                    167:       er.grade = minGrade;
                    168:       er.gseti = minGseti;
                    169:       return er;
                    170:     }
                    171:   }else{
                    172:       er.ell = ell2;
                    173:       er.first = 0;
                    174:       er.ggi = minGgi;
                    175:       return er;
                    176:   }
                    177: }
                    178:
                    179: /*
                    180:   r and gset are assumed to be (0,1)-homogenized (h-homogenized)
                    181:  */
                    182: POLY reduction_ecart(r,gset,needSyz,syzp)
                    183:      POLY r;
                    184:      struct gradedPolySet *gset;
                    185:      int needSyz;
                    186:      struct syz0 *syzp; /* set */
                    187: {
                    188:   int reduced,reduced1,reduced2;
                    189:   int grd;
                    190:   struct polySet *set;
                    191:   POLY cf,syz;
                    192:   int i;
                    193:   POLY cc,cg;
                    194:   struct ecartReducer ells;
                    195:   struct ecartPolyArray *gg;
                    196:   POLY pp;
                    197:   int ell;
                    198:
                    199:   extern struct ring *CurrentRingp;
                    200:   struct ring *rp;
                    201:
                    202:   gg = NULL;
                    203:   if (needSyz) {
                    204:     if (r ISZERO) { rp = CurrentRingp; } else { rp = r->m->ringp; }
                    205:     cf = cxx(1,0,0,rp);
                    206:     syz = ZERO;
                    207:   }
                    208:
                    209:   if (r != POLYNULL) {
1.2     ! takayama  210:        rp = r->m->ringp;
        !           211:        if (! rp->weightedHomogenization) {
1.1       takayama  212:          errorKan1("%s\n","ecart.c: the given ring must be declared with [(weightedHomogenization) 1]");
                    213:        }
                    214:   }
                    215:
                    216:   r = goHomogenize11(r,DegreeShifto_vec,DegreeShifto_size,-1,1);
                    217:   /* 1 means homogenize only s */
                    218:
                    219:   do {
                    220:     if (DebugReductionRed) printf("r=%s\n",POLYToString(r,'*',1));
                    221:     ells = ecartFindReducer(r,gset,gg);
                    222:     ell = ells.ell;
                    223:     if (ell > 0) {
                    224:       gg = ecartPutPolyInG(r,gg);
                    225:     }
                    226:     if (ell >= 0) {
                    227:       if (ells.first) {
                    228:         pp = ((gset->polys[ells.grade])->gh)[ells.gseti];
                    229:       }else{
                    230:         pp = (gg->pa)[ells.ggi];
                    231:       }
                    232:       if (ell > 0) r = mpMult(cxx(1,0,ell,rp),r); /* r = s^ell r */
                    233:       r = (*reduction1)(r,pp,needSyz,&cc,&cg);
                    234:       if (needSyz) {
                    235:         if (ells.first) {
                    236:           cf = ppMult(cc,cf);
                    237:           syz = cpMult(toSyzCoeff(cc),syz);
                    238:           syz = ppAddv(syz,toSyzPoly(cg,ells.grade,ells.gseti));
                    239:         }else{
                    240:           /* BUG: not yet */
                    241:         }
                    242:       }
                    243:       if (r ISZERO) goto ss;
                    244:       r = ecartDivideSv(r); /* r = r/s^? */
                    245:     }
                    246:   }while (ell >= 0);
                    247:
                    248:  ss: ;
                    249:   if (needSyz) {
                    250:     syzp->cf = cf;   /* cf is in the CurrentRingp */
                    251:     syzp->syz = syz; /* syz is in the SyzRingp */
                    252:     /* BUG: dehomogenize the syzygy */
                    253:   }
1.2     ! takayama  254:   /*
1.1       takayama  255:   r = goDeHomogenizeS(r);
1.2     ! takayama  256:   */
1.1       takayama  257:   return(r);
                    258: }

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>