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

Annotation of OpenXM/src/kan96xx/Kan/replace.c, Revision 1.1.1.1

1.1       maekawa     1: #include <stdio.h>
                      2: #include "datatype.h"
                      3: #include "extern2.h"
                      4:
                      5: static int badLRule(POLY set[],int num);
                      6:
                      7: POLY mReplace(mm,lSideX,rSideX,sizex,lSideD,rSideD,sized,commutative)
                      8: POLY mm;
                      9: int lSideX[];
                     10: POLY rSideX[];  /* Rule: a=lSideX[i], x_a ---> rSideX[i] */
                     11: int sizex;
                     12: int lSideD[];   /* Rule: b=lSideD[i], D_b ---> rSideD[i] */
                     13: POLY rSideD[];
                     14: int sized;
                     15: int commutative;
                     16: {
                     17: /* The function should be tuned by using a table. */
                     18:   POLY rp;
                     19:   POLY fac;
                     20:   int i;
                     21:   int j;
                     22:   int flag;
                     23:   MONOMIAL mp;
                     24:   int n;
                     25:
                     26:   mp = mm->m;
                     27:   rp = newCell(coeffCopy(mm->coeffp),newMonomial(mp->ringp));
                     28:   n = mp->ringp->n;
                     29:
                     30:   /* Multiply backwards. It's faster. */
                     31:   for (i=n-1; i>=0; i--) {
                     32:     if (mp->e[i].D != 0) {
                     33:       flag = 1;
                     34:       for (j=0; j<sized; j++) {
                     35:        if (lSideD[j] == i) {
                     36:          if (commutative) fac = pPower_poly(rSideD[j],mp->e[i].D);
                     37:          else fac = pPower(rSideD[j],mp->e[i].D);
                     38:          /* if negative, this does not work  well!! */
                     39:          flag = 0;
                     40:          break;
                     41:        }
                     42:       }
                     43:       if (flag) fac = cdd(1,i,mp->e[i].D,mp->ringp);
                     44:       if (commutative) rp = ppMult_poly(fac,rp);
                     45:       else rp = ppMult(fac,rp);
                     46:     }
                     47:   }
                     48:
                     49:   for (i=n-1; i>=0; i--) {
                     50:     if (mp->e[i].x != 0) {
                     51:       flag = 1;
                     52:       for (j=0; j<sizex; j++) {
                     53:        if (lSideX[j] == i) {
                     54:          if (commutative) fac = pPower_poly(rSideX[j],mp->e[i].x);
                     55:          else fac = pPower(rSideX[j],mp->e[i].x);
                     56:          /* if negative, this does not work well!! */
                     57:          flag = 0;
                     58:          break;
                     59:        }
                     60:       }
                     61:       if (flag) fac = cxx(1,i,mp->e[i].x,mp->ringp);
                     62:       if (commutative) rp = ppMult_poly(fac,rp);
                     63:       else rp = ppMult(fac,rp);
                     64:     }
                     65:   }
                     66:
                     67:   return(rp);
                     68: }
                     69:
                     70: /*
                     71:   lRule[i] ---> rRule[i]
                     72: */
                     73: POLY replace(f,lRule,rRule,num)
                     74: POLY f;
                     75: POLY lRule[];  /* lRule[i] must be x0 or ... or D{N-1} */
                     76: POLY rRule[];
                     77: int num;
                     78: {
                     79:   POLY rSideX[N0];
                     80:   POLY rSideD[N0];
                     81:   int  lSideX[N0];
                     82:   int  lSideD[N0];
                     83:   int i,j,size;
                     84:   int sizex,sized;
                     85:   POLY rp;
                     86:   int n;
                     87:
                     88:   /***printf("f=%s\n",POLYToString(f,'*',0));
                     89:   for (i=0; i<num; i++) {
                     90:     printf("%s --> %s\n",POLYToString(lRule[i],'*',0),
                     91:           POLYToString(rRule[i],'*',0));
                     92:   }
                     93:   printf("\n"); ***/
                     94:
                     95:   if (f ISZERO) return(ZERO);
                     96:   sizex = sized = 0;
                     97:   if (badLRule(lRule,num)) {
                     98:     warningPoly("replace(): lRule must be a variable.");
                     99:     return(ZERO);
                    100:   }
                    101:   n = f->m->ringp->n;
                    102:   /* x-part */
                    103:   for (i=0; i<n; i++) {
                    104:     for (j=0; j<num; j++) {
                    105:       if (lRule[j]->m->e[i].x == 1) {
                    106:        lSideX[sizex] = i;
                    107:        rSideX[sizex] = rRule[j];
                    108:        sizex++;
                    109:        if (sizex >= n) {
                    110:          warningPoly(" replace(): too many x-rules . ");
                    111:          sizex--;
                    112:        }
                    113:        break;
                    114:       }
                    115:     }
                    116:   }
                    117:   /* D-part */
                    118:   for (i=0; i<n; i++) {
                    119:     for (j=0; j<num; j++) {
                    120:       if (lRule[j]->m->e[i].D == 1) {
                    121:        lSideD[sized] = i;
                    122:        rSideD[sized] = rRule[j];
                    123:        sized++;
                    124:        if (sized >= n) {
                    125:          warningPoly(" replacen(): too many D-rules . ");
                    126:          sized--;
                    127:        }
                    128:        break;
                    129:       }
                    130:     }
                    131:   }
                    132:
                    133:   rp = ZERO;
                    134:   if (f ISZERO) return(ZERO);
                    135:   while(f != POLYNULL) {
                    136:     rp = ppAddv(rp,mReplace(f,lSideX,rSideX,sizex,lSideD,rSideD,sized,0));
                    137:     f = f->next;
                    138:   }
                    139:   return(rp);
                    140: }
                    141:
                    142:
                    143: /* For the dirty trick of mpMult_difference */
                    144: POLY replace_poly(f,lRule,rRule,num)
                    145: POLY f;
                    146: POLY lRule[];  /* lRule[i] must be x0 or ... or D{N-1} */
                    147: POLY rRule[];
                    148: int num;
                    149: {
                    150:   POLY rSideX[N0];
                    151:   POLY rSideD[N0];
                    152:   int  lSideX[N0];
                    153:   int  lSideD[N0];
                    154:   int i,j,size;
                    155:   int sizex,sized;
                    156:   POLY rp;
                    157:   int n;
                    158:
                    159:   /***printf("f=%s\n",POLYToString(f,'*',0));
                    160:   for (i=0; i<num; i++) {
                    161:     printf("%s --> %s\n",POLYToString(lRule[i],'*',0),
                    162:           POLYToString(rRule[i],'*',0));
                    163:   }
                    164:   printf("\n"); ***/
                    165:
                    166:   if (f ISZERO) return(ZERO);
                    167:   sizex = sized = 0;
                    168:   if (badLRule(lRule,num)) {
                    169:     warningPoly("replace(): lRule must be a variable.");
                    170:     return(ZERO);
                    171:   }
                    172:   n = f->m->ringp->n;
                    173:   /* x-part */
                    174:   for (i=0; i<n; i++) {
                    175:     for (j=0; j<num; j++) {
                    176:       if (lRule[j]->m->e[i].x == 1) {
                    177:        lSideX[sizex] = i;
                    178:        rSideX[sizex] = rRule[j];
                    179:        sizex++;
                    180:        if (sizex >= n) {
                    181:          warningPoly(" replace(): too many x-rules . ");
                    182:          sizex--;
                    183:        }
                    184:        break;
                    185:       }
                    186:     }
                    187:   }
                    188:   /* D-part */
                    189:   for (i=0; i<n; i++) {
                    190:     for (j=0; j<num; j++) {
                    191:       if (lRule[j]->m->e[i].D == 1) {
                    192:        lSideD[sized] = i;
                    193:        rSideD[sized] = rRule[j];
                    194:        sized++;
                    195:        if (sized >= n) {
                    196:          warningPoly(" replacen(): too many D-rules . ");
                    197:          sized--;
                    198:        }
                    199:        break;
                    200:       }
                    201:     }
                    202:   }
                    203:
                    204:   rp = ZERO;
                    205:   if (f ISZERO) return(ZERO);
                    206:   while(f != POLYNULL) {
                    207:     rp = ppAddv(rp,mReplace(f,lSideX,rSideX,sizex,lSideD,rSideD,sized,1));
                    208:     f = f->next;
                    209:   }
                    210:   return(rp);
                    211: }
                    212:
                    213:
                    214: static int badLRule(set,num)
                    215: POLY set[];
                    216: int num;
                    217: { int i;
                    218:   for (i=0; i<num; i++) {
                    219:     if (set[0] ISZERO) {
                    220:       return(1);
                    221:     }
                    222:   }
                    223:   return(0);
                    224: }
                    225:

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