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>