[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.2

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

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