Annotation of OpenXM/src/k097/lib/restriction/complex.k, Revision 1.2
1.2 ! takayama 1: /* $OpenXM: OpenXM/src/k097/lib/restriction/complex.k,v 1.1 2001/01/04 12:29:33 takayama Exp $ */
! 2: /* Document of this module is at k097/Doc/complex.texi */
1.1 takayama 3:
4: load["lib/restriction/restriction.k"];;
5:
6: def Res_solv(m,d,rng) {
7: local r,rr,ans,ac;
8: ac = Length(Arglist);
9: r = GetRing(Poly("1")); /* Save the current ring. */
10: if (ac < 3) {
11: rng = null;
12: rr = GetRing(m);
13: if (Tag(rr) == 0) rr = GetRing(d);
14: if (Tag(rr) != 0) SetRing(rr);
15: }else{
16: SetRing(rng);
17: }
18: m=DC(m,"polynomial"); d = DC(d,"polynomial");
19:
20: if (IsRing(rng)) {
21: sm1(" [ m d rng] res-solv /ans set ");
22: }else{
23: sm1(" [ m d] res-solv /ans set ");
24: }
25:
26: SetRing(r);
27: return(ans);
28: }
29:
30: /* m : D^p ---> D^q/jj, u m = d mod jj */
31: def Res_solv2(m,d,jj,rng) {
32: local r,rr,ans,ac,pp,qq,kk;
33:
34: ac = Length(Arglist);
35: r = GetRing(Poly("1")); /* Save the current ring. */
36: if (ac < 4) {
37: rng = null;
38: rng = GetRing(m);
39: if (Tag(rng) == 0) rng = GetRing(d);
40: }
41:
42: pp = Length(m);
43: if (!IsArray(m[0])) {
44: sm1(" m { [ 2 1 roll ] } map /m set ");
45: }
46: qq = Length(m[0]);
47: if (Length(jj) > 0) {
48: if (!IsArray(jj[0])) {
49: sm1(" jj { [ 2 1 roll ] } map /jj set ");
50: }
51: if (qq != Length(jj[0])) {
52: Error(" Matrix size mismatch in m and jj of Kernel2(m,jj,r).");
53: }
54: }
55: m = Join(m,jj);
56:
57: ans = Res_solv(m,d,rng);
58: /* Println(ans); */
59: SetRing(r);
60: return([Firstn(ans[0],pp),ans[1]]);
61: }
62: /* Res_solv2([x,y],[x^2+y^2],[x]):*/
63:
64: def Res_solv_h(m,d,rng) {
65: local r,rr,ans,ac;
66: ac = Length(Arglist);
67: r = GetRing(Poly("1")); /* Save the current ring. */
68: if (ac < 3) {
69: rng = null;
70: rr = GetRing(m);
71: if (Tag(rr) == 0) rr = GetRing(d);
72: if (Tag(rr) != 0) SetRing(rr);
73: }else{
74: SetRing(rng);
75: }
76: m=DC(m,"polynomial"); d = DC(d,"polynomial");
77:
78: if (IsRing(rng)) {
79: sm1(" [ m d rng] res-solv-h /ans set ");
80: }else{
81: sm1(" [ m d] res-solv-h /ans set ");
82: }
83:
84: SetRing(r);
85: return(ans);
86: }
87:
88: /* m : D^p ---> D^q/jj, u m = d mod jj */
89: def Res_solv2_h(m,d,jj,rng) {
90: local r,rr,ans,ac,pp,qq,kk;
91:
92: ac = Length(Arglist);
93: r = GetRing(Poly("1")); /* Save the current ring. */
94: if (ac < 4) {
95: rng = null;
96: rng = GetRing(m);
97: if (Tag(rng) == 0) rng = GetRing(d);
98: }
99:
100: pp = Length(m);
101: if (!IsArray(m[0])) {
102: sm1(" m { [ 2 1 roll ] } map /m set ");
103: }
104: qq = Length(m[0]);
105: if (Length(jj) > 0) {
106: if (!IsArray(jj[0])) {
107: sm1(" jj { [ 2 1 roll ] } map /jj set ");
108: }
109: if (qq != Length(jj[0])) {
110: Error(" Matrix size mismatch in m and jj of Kernel2(m,jj,r).");
111: }
112: }
113: m = Join(m,jj);
114:
115: ans = Res_solv_h(m,d,rng);
116: /* Println(ans); */
117: SetRing(r);
118: return([Firstn(ans[0],pp),ans[1]]);
119: }
120: /* Res_solv2_h([x,y],[x^2+y^2],[x]):*/
121:
122: def Getxvars() {
123: local v,n,i,ans,ans2;
124: sm1(" getvNames /v set ");
125: sm1(" [(NN)] system_variable (universalNumber) dc /n set ");
126: ans = [];
127: for (i=1; i<n; i++) {
128: ans = Append(ans,v[i]);
129: }
130: sm1(" ans from_records /ans2 set ");
131: return([ans,ans2]);
132: }
133:
134: /* This works only for D cf. Getxvars(). */
135: def Intersection(i1,i2,rng) {
136: local r,rr,ans,n,tt,vv,ac;
137: ac = Length(Arglist);
138: r = GetRing(Poly("1")); /* Save the current ring */
139: if (ac < 3) {
140: rr = GetRing(i1);
141: if (Tag(rr) == 0) rr = GetRing(i2);
142: if (Tag(rr) != 0) SetRing(rr);
143: }else{
144: SetRing(rng);
145: }
146: /*
147: i1=DC(i1,"polynomial"); i2 = DC(i2,"polynomial");
148: */
149: i1 = ReParse(i1); i2 = ReParse(i2);
150:
151: vv = Getxvars();
152: vv = vv[1];
153: if (Length(i1) == 0) {
154: ans = i2;
155: }else if (!IsArray(i1[0])) {
156: sm1(" [i1 i2 vv] intersection /ans set ");
157: }else {
158: n = Length(i1[0]);
159: sm1(" i1 fromVectors /i1 set ");
160: sm1(" i2 fromVectors /i2 set ");
161: sm1(" [i1 i2 vv] intersection /tt set ");
162: sm1(" [n (integer) dc tt] toVectors /ans set ");
163: }
164:
165: SetRing(r);
166: return(ans);
167: }
168:
169: def Firstn(a,n) {
170: local r,i,j,ans;
171: if (Length(a) == 0) {
172: return([ ]);
173: }
174: if (!IsArray(a[0])) {
175: r = NewArray(n);
176: for (i=0; i<n; i++) {
177: r[i] = a[i];
178: }
179: return(r);
180: }else{
181: ans = [ ];
182: for (j=0; j < Length(a); j++) {
183: r = NewArray(n);
184: for (i=0; i<n; i++) {
185: r[i] = a[j,i];
186: }
187: ans = Append(ans,r);
188: }
189: return(ans);
190: }
191: }
192:
193: /* Kernel is also defined in lib/minimal/minimal.k */
194: def Kernel(f,v) {
195: local ans;
196: /* v : string or ring */
197: if (Length(Arglist) < 2) {
198: sm1(" [f] syz /ans set ");
199: }else{
200: sm1(" [f v] syz /ans set ");
201: }
202: return(ans);
203: }
204:
205: def Kernel_h(f,v) {
206: local ans;
207: /* v : string or ring */
208: if (Length(Arglist) < 2) {
209: sm1(" [f] syz_h /ans set ");
210: }else{
211: sm1(" [f v] syz_h /ans set ");
212: }
213: return(ans);
214: }
215:
216: /* Kernel of (D^p --- m ---> D^q/jj) */
217: def Kernel2(m,jj,r)
218: {
219: local crng,ac,pp,qq,kk;
220: ac = Length(Arglist);
221: crng = GetRing(Poly("1"));
222: if (ac < 3) {
223: r = GetRing(m);
224: }
225: pp = Length(m);
226: if (!IsArray(m[0])) {
227: sm1(" m { [ 2 1 roll ] } map /m set ");
228: }
229: qq = Length(m[0]);
230: if (Length(jj) > 0) {
231: if (!IsArray(jj[0])) {
232: sm1(" jj { [ 2 1 roll ] } map /jj set ");
233: }
234: if (qq != Length(jj[0])) {
235: Error(" Matrix size mismatch in m and jj of Kernel2(m,jj,r).");
236: }
237: }
238: m = Join(m,jj);
239: kk = Kernel(m,r);
240: SetRing(crng);
241: return(Firstn(kk[0],pp));
242: }
243:
244: /* Kernel of (D^p --- m ---> D^q/jj) */
245: def Kernel2_h(m,jj,r)
246: {
247: local crng,ac,pp,qq,kk;
248: ac = Length(Arglist);
249: crng = GetRing(Poly("1"));
250: if (ac < 3) {
251: r = GetRing(m);
252: }
253: pp = Length(m);
254: if (!IsArray(m[0])) {
255: sm1(" m { [ 2 1 roll ] } map /m set ");
256: }
257: qq = Length(m[0]);
258: if (Length(jj) > 0) {
259: if (!IsArray(jj[0])) {
260: sm1(" jj { [ 2 1 roll ] } map /jj set ");
261: }
262: if (qq != Length(jj[0])) {
263: Error(" Matrix size mismatch in m and jj of Kernel2(m,jj,r).");
264: }
265: }
266: m = Join(m,jj);
267: kk = Kernel_h(m,r);
268: SetRing(crng);
269: return(Firstn(kk[0],pp));
270: }
271:
272: /* From lib/minimal/minimal.k */
273: def Gb(m,rng) {
274: local r,rr,ans,ac;
275: ac = Length(Arglist);
276: r = GetRing(Poly("1")); /* Save the current ring. */
277: if (ac < 2) {
278: rng = null;
279: rr = GetRing(m);
280: if (Tag(rr) != 0) SetRing(rr);
281: }else{
282: rr = rng;
283: SetRing(rr);
284: }
285: /* m=DC(m,"polynomial"); */
286: m = ReParse(m);
287: sm1(" [m rr] gb /ans set ");
288: SetRing(r);
289: return(ans);
290: }
291:
292: def Gb_h(m,rng) {
293: local r,rr,ans,ac;
294: ac = Length(Arglist);
295: r = GetRing(Poly("1")); /* Save the current ring. */
296: if (ac < 2) {
297: rng = null;
298: rr = GetRing(m);
299: if (Tag(rr) != 0) SetRing(rr);
300: }else{
301: rr = rng;
302: SetRing(rr);
303: }
304: /* m=DC(m,"polynomial"); */
305: m = ReParse(m);
306: sm1(" [m rr] gb_h /ans set ");
307: SetRing(r);
308: return(ans);
309: }
310:
311: def Res_shiftMatrix(m,v,rng) {
312: local n,ans,r,ac,i,j,b1,b2;
313: sm1(" 40 (string) dc /b1 set ");
314: sm1(" 41 (string) dc /b2 set ");
315: ac = Length(Arglist);
316: r = GetRing(Poly("1")); /* Save the current ring. */
317: if (ac < 3) {
318: }else{
319: SetRing(rng);
320: }
321: n = Length(m);
322: ans = NewVector(n);
323: for (i=0; i<n; i++) {
324: ans[i] = NewVector(n);
325: for (j=0; j<n; j++) {
326: ans[i,j] = Poly("0");
327: }
328: ans[i,i] = Poly(AddString([
329: DC(v,"string"),"^",b1,DC(m[i],"string"),b2]));
330: }
331: SetRing(r);
332: return(ans);
333: }
334:
1.2 ! takayama 335: def ChangeRing(f) {
! 336: local r;
! 337: r = GetRing(f);
! 338: if (Tag(r) == 14) {
! 339: SetRing(r);
! 340: return(true);
1.1 takayama 341: }else{
1.2 ! takayama 342: return(false);
1.1 takayama 343: }
344: }
345:
346:
347: def test2() {
348: RingD("x,y,z");
349: /* Step 1. J */
350: Println(" ---------- J --------------");
351: mm = [[1],
352: [x*Dx],
353: [y*Dy],
354: [z*Dz]];
355: b1 = Res_solv(mm,[Dz]);
356: Println(b1);
357: b2 = Res_solv(mm,[Dy]);
358: Println(b2);
359: /* Step 2. K */
360: Println(" --------- K -------------");
361: mm2 = [[Dz],
362: [Dy],
363: [x*Dx],
364: [y*Dy],
365: [z*Dz]];
366: k1 = Kernel(mm2);
367: Pmat(Firstn(k1[0],2));
368: aa=Kernel2([[Dz],[Dy]],[x*Dx,y*Dy,z*Dz]);
369: Pmat(aa);
370: }
371:
372: def test3() {
373: RingD("x,y,z");
374: mm = [[-Dz],[Dy],[x*Dx],[0]];
375: kk = Kernel(mm);
376: Pmat(kk[0]);
377: rrr= RingD("x,y,z,t",[["t",1,"x",-1,"y",-1,"z",-1,
378: "Dx",1,"Dy",1,"Dz",1]]);
379: kk0 = Reparse(kk[0]);
380: gg = Gb(kk0*Res_shiftMatrix([1,1,0,2],t));
381: Pmat(gg[0]); Pmat(gg[1]);
382: gg2= Substitute(Substitute(gg[0],t,1),h,1);
383: Pmat(gg2);
384: Println("----------------------------");
385: mm2 = [[0,0],
386: [0,0],
387: [0,0],
388: [-Dy,-Dz]];
389: jj2 = [[x*Dx,0],
390: [y*Dy,0],
391: [z,0],
392: [0,x*Dx],
393: [0,y],
394: [0,z*Dz]];
395: kk2 = Kernel2(mm2,jj2);
396: Pmat(kk2);
397: Println("-----------------------");
398: ii = Intersection(gg2,kk2);
399: Pmat(ii);
400: SetRing(rrr);
401: ii = Reparse(ii);
402: gg3 = Gb(ii*Res_shiftMatrix([1,1,0,2],t));
403: Pmat(gg3[0]);
404: gg4= Substitute(Substitute(gg3[0],t,1),h,1);
405: Println("---- page 20, observation 4 -----");
406: Pmat(gg4);
407: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>