1 | // version string automatically expanded by CVS |
---|
2 | version="$Id: powers.lib,v 1.3 2008-06-26 13:17:43 motsak Exp $"; |
---|
3 | category="Linear algebra"; |
---|
4 | |
---|
5 | info=" |
---|
6 | LIBRARY: powers.lib Symmetric and Exterior k-powers |
---|
7 | AUTHOR: Oleksandr Motsak, email: {U@D}, where U={motsak}, D={mathematik.uni-kl.de}. |
---|
8 | |
---|
9 | KEYWORDS: symmetric power, exterior power |
---|
10 | |
---|
11 | PROCEDURES: |
---|
12 | symmetricBasis(int, int) return a basis of a k-th symmetric power |
---|
13 | symmetricPower(module,int) return k-th symmetric power of a matrix |
---|
14 | |
---|
15 | exteriorBasis(int, int) return a basis of a k-th exterior power |
---|
16 | exteriorPower(module,int) return k-th exterior power of a matrix |
---|
17 | |
---|
18 | "; |
---|
19 | |
---|
20 | |
---|
21 | LIB "general.lib"; // for sort |
---|
22 | LIB "nctools.lib"; // for SuperCommutative |
---|
23 | |
---|
24 | // thanks to Oleksandr Iena for his persistence ;-) |
---|
25 | |
---|
26 | proc symmetricBasis(int n, int k) |
---|
27 | "USAGE: symmetricBasis(n, k); n int, k int |
---|
28 | RETURN: ring: the basis of the k-th symmetric power of a n-dim vector space |
---|
29 | NOTE: The output consists of a polynomial ring in n variables, together with the ideal called I. |
---|
30 | The ideal I is the basis of the k-th symmetric power of a n-dim vector space (ordered lexicographically). |
---|
31 | The char. of the returned ring is the same as of current basis ring or zero. |
---|
32 | SEE ALSO: exteriorBasis |
---|
33 | KEYWORDS: symmetric basis |
---|
34 | EXAMPLE: example symmetricBasis; shows an example" |
---|
35 | { |
---|
36 | int p = 0; |
---|
37 | if( nameof(basering) != "basering" ) |
---|
38 | { |
---|
39 | p = char(basering); |
---|
40 | } |
---|
41 | |
---|
42 | ring T = (p), (e(1..n)), dp; |
---|
43 | ideal I = simplify( NF(maxideal(k), std(0)), 1 + 2 + 8 ); |
---|
44 | int N = ncols(I); |
---|
45 | I = sort(I)[1]; // lex |
---|
46 | I = I[N..1]; |
---|
47 | |
---|
48 | export I; |
---|
49 | return(T); |
---|
50 | } |
---|
51 | example |
---|
52 | { "EXAMPLE:"; echo = 2; |
---|
53 | |
---|
54 | def r = symmetricBasis(2, 3); setring r; r; |
---|
55 | I; // basis of 3rd sym. power of a 2-dim v.s. |
---|
56 | kill r; |
---|
57 | |
---|
58 | ring r = (0),(a, b, c, d), dp; r; |
---|
59 | def g = symmetricBasis(3, 2); setring g; g; I; |
---|
60 | |
---|
61 | kill g, r; |
---|
62 | |
---|
63 | ring r = (32003),(a, b, c, d), dp; r; |
---|
64 | def g = symmetricBasis(4, 2); setring g; g; I; |
---|
65 | } |
---|
66 | |
---|
67 | proc exteriorBasis(int n, int k) |
---|
68 | "USAGE: exteriorBasis(n, k); n int, k int |
---|
69 | RETURN: ring: the basis of the k-th exterior power of a n-dim vector space |
---|
70 | NOTE: The output consists of a polynomial ring in n variables, together with the ideal called I. |
---|
71 | The ideal I is the basis of the k-th exterior power of a n-dim vector space (ordered lexicographically). |
---|
72 | The char. of the returned ring is the same as of current basis ring or zero. |
---|
73 | |
---|
74 | SEE ALSO: symmetricBasis |
---|
75 | KEYWORDS: exterior basis |
---|
76 | EXAMPLE: example exteriorBasis; shows an example" |
---|
77 | { |
---|
78 | int p = 0; |
---|
79 | if( nameof(basering) != "basering" ) |
---|
80 | { |
---|
81 | p = char(basering); |
---|
82 | } |
---|
83 | |
---|
84 | ring S = (p), (e(1..n)), dp; |
---|
85 | def T = SuperCommutative(); setring T; |
---|
86 | ideal I = simplify( NF(maxideal(k), std(0)), 1 + 2 + 8 ); |
---|
87 | int N = ncols(I); |
---|
88 | I = sort(I)[1]; |
---|
89 | I = I[N..1 ]; |
---|
90 | |
---|
91 | export I; |
---|
92 | return(T); |
---|
93 | } |
---|
94 | example |
---|
95 | { "EXAMPLE:"; echo = 2; |
---|
96 | |
---|
97 | def r = exteriorBasis(2, 3); setring r; r; |
---|
98 | "Basis: ", I; // basis of 3rd sym. power of a 2-dim v.s. |
---|
99 | kill r; |
---|
100 | |
---|
101 | ring r = (0),(a, b, c, d), dp; r; |
---|
102 | def g = exteriorBasis(3, 2); setring g; g; "Basis: ", I; |
---|
103 | |
---|
104 | kill g, r; |
---|
105 | |
---|
106 | ring r = (32003),(a, b, c, d), dp; r; |
---|
107 | def g = exteriorBasis(4, 2); setring g; g; "Basis: ", I; |
---|
108 | } |
---|
109 | |
---|
110 | |
---|
111 | |
---|
112 | proc symmetricPower(module A, int k) |
---|
113 | "USAGE: symmetricPower(A, k); A module, k int |
---|
114 | RETURN: module: the k-th symmetric power of A |
---|
115 | NOTE: the chosen bases (ordered lexicographically) and |
---|
116 | temporary data may be shown if printlevel is big enough |
---|
117 | SEE ALSO: exteriorPower |
---|
118 | KEYWORDS: symmetric power |
---|
119 | EXAMPLE: example symmetricPower; shows an example" |
---|
120 | { |
---|
121 | int p = printlevel - voice + 2; |
---|
122 | def save = basering; |
---|
123 | |
---|
124 | int M = nrows(A); |
---|
125 | int N = ncols(A); |
---|
126 | |
---|
127 | int i, j; |
---|
128 | |
---|
129 | ///////////////////////////////////////////////////////////////// |
---|
130 | def Tn = symmetricBasis(N, k); setring Tn; |
---|
131 | ideal Ink = I; |
---|
132 | dbprint(p-3, "Temporary Source Ring", basering); |
---|
133 | dbprint(p, "S^k(Source Basis)", Ink); |
---|
134 | |
---|
135 | ///////////////////////////////////////////////////////////////// |
---|
136 | def Tm = symmetricBasis(M, k); setring Tm; |
---|
137 | ideal Imk = I; |
---|
138 | dbprint(p-3, "Temporary Image Ring", basering); |
---|
139 | dbprint(p, "S^k(Image Basis)", Imk); |
---|
140 | |
---|
141 | ///////////////////////////////////////////////////////////////// |
---|
142 | def Rm = save + Tm; |
---|
143 | setring Rm; |
---|
144 | dbprint(p-2, "Temporary Working Ring", Rm); |
---|
145 | |
---|
146 | module A = imap(save, A); |
---|
147 | |
---|
148 | ideal B; poly t; |
---|
149 | |
---|
150 | for( i = N; i > 0; i-- ) |
---|
151 | { |
---|
152 | t = 0; |
---|
153 | for( j = M; j > 0; j-- ) |
---|
154 | { |
---|
155 | t = t + A[i][j] * var( nvars(save) + j); // tensor product!!! |
---|
156 | } |
---|
157 | |
---|
158 | B[i] = t; |
---|
159 | } |
---|
160 | |
---|
161 | dbprint(p-1, "Matrix of simgle images", B); |
---|
162 | |
---|
163 | map TMap = Tn, B; ideal C = TMap(Ink); // apply S^k A to Source basis vectors... (Ink) |
---|
164 | C = NF(C, std(0)); |
---|
165 | |
---|
166 | dbprint(p-1, "Image Matrix: ", C); |
---|
167 | |
---|
168 | // and write it in Image basis (Imk) |
---|
169 | |
---|
170 | ideal Imk = imap(Tm, Imk); |
---|
171 | |
---|
172 | module D; poly lm; vector tt; |
---|
173 | |
---|
174 | for( i = ncols(C); i > 0; i-- ) |
---|
175 | { |
---|
176 | t = C[i]; |
---|
177 | tt = 0; |
---|
178 | |
---|
179 | while( t != 0 ) |
---|
180 | { |
---|
181 | lm = leadmonom(t); |
---|
182 | // lm; |
---|
183 | for( j = ncols(Imk); j > 0; j-- ) |
---|
184 | { |
---|
185 | if( lm / Imk[j] != 0 ) |
---|
186 | { |
---|
187 | tt = tt + (lead(t) / Imk[j]) * gen(j); |
---|
188 | break; |
---|
189 | } |
---|
190 | } |
---|
191 | t = t - lead(t); |
---|
192 | } |
---|
193 | |
---|
194 | D[i] = tt; |
---|
195 | // tt; |
---|
196 | } |
---|
197 | |
---|
198 | // D; print(D); |
---|
199 | |
---|
200 | |
---|
201 | setring save; |
---|
202 | |
---|
203 | // basering; |
---|
204 | |
---|
205 | module AA = imap(Rm, D); |
---|
206 | |
---|
207 | // nrows(AA) = Nk; // ???? |
---|
208 | // ncols(AA) = Mk; |
---|
209 | |
---|
210 | // Nk, Mk; |
---|
211 | |
---|
212 | // AA[Mk] = AA[Mk] * 1; |
---|
213 | |
---|
214 | return(AA); |
---|
215 | } |
---|
216 | example |
---|
217 | { "EXAMPLE:"; echo = 2; |
---|
218 | int save = printlevel; printlevel = 1; |
---|
219 | |
---|
220 | ring r = (0),(a, b, c, d), dp; r; |
---|
221 | |
---|
222 | module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1); |
---|
223 | |
---|
224 | print(symmetricPower(A, 2)); |
---|
225 | print(symmetricPower(A, 3)); |
---|
226 | |
---|
227 | module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2); |
---|
228 | |
---|
229 | print(symmetricPower(B, 2)); |
---|
230 | print(symmetricPower(B, 3)); |
---|
231 | |
---|
232 | kill r; |
---|
233 | ring r = (0),(a, b, c, d), dp;def g = SuperCommutative();setring g; |
---|
234 | |
---|
235 | module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1); |
---|
236 | |
---|
237 | print(symmetricPower(A, 2)); |
---|
238 | print(symmetricPower(A, 3)); |
---|
239 | |
---|
240 | module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2); |
---|
241 | |
---|
242 | print(symmetricPower(B, 2)); |
---|
243 | print(symmetricPower(B, 3)); |
---|
244 | |
---|
245 | kill r; |
---|
246 | |
---|
247 | printlevel = save; |
---|
248 | } |
---|
249 | |
---|
250 | |
---|
251 | |
---|
252 | proc exteriorPower(module A, int k) |
---|
253 | "USAGE: exteriorPower(A, k); A module, k int |
---|
254 | RETURN: module: the k-th exterior power of A |
---|
255 | NOTE: the chosen bases and temporary data |
---|
256 | may be shown if printlevel is big enough. |
---|
257 | Last rows may be invisible if zero. |
---|
258 | SEE ALSO: symmetricPower |
---|
259 | KEYWORDS: exterior power |
---|
260 | EXAMPLE: example exteriorPower; shows an example" |
---|
261 | { |
---|
262 | int p = printlevel - voice + 2; |
---|
263 | def save = basering; |
---|
264 | |
---|
265 | int M = nrows(A); |
---|
266 | int N = ncols(A); |
---|
267 | |
---|
268 | int i, j; |
---|
269 | |
---|
270 | ///////////////////////////////////////////////////////////////// |
---|
271 | def Tn = exteriorBasis(N, k); setring Tn; |
---|
272 | ideal Ink = I; |
---|
273 | dbprint(p-3, "Temporary Source Ring", basering); |
---|
274 | dbprint(p, "E^k(Source Basis)", Ink); |
---|
275 | |
---|
276 | ///////////////////////////////////////////////////////////////// |
---|
277 | def Tm = exteriorBasis(M, k); setring Tm; |
---|
278 | ideal Imk = I; |
---|
279 | dbprint(p-3, "Temporary Image Ring", basering); |
---|
280 | dbprint(p, "E^k(Image Basis)", Imk); |
---|
281 | |
---|
282 | |
---|
283 | def Rm = save + Tm; |
---|
284 | setring Rm; |
---|
285 | dbprint(p-2, "Temporary Working Ring", Rm); |
---|
286 | |
---|
287 | module A = imap(save, A); |
---|
288 | |
---|
289 | ideal B; poly t; |
---|
290 | |
---|
291 | for( i = N; i > 0; i-- ) |
---|
292 | { |
---|
293 | t = 0; |
---|
294 | for( j = M; j > 0; j-- ) |
---|
295 | { |
---|
296 | t = t + A[i][j] * var( nvars(save) + j); // tensor product!!! |
---|
297 | } |
---|
298 | |
---|
299 | B[i] = t; |
---|
300 | } |
---|
301 | |
---|
302 | dbprint(p-1, "Matrix of simgle images", B); |
---|
303 | |
---|
304 | map TMap = Tn, B; ideal C = TMap(Ink); // apply S^k A to Source basis vectors... (Ink) |
---|
305 | C = NF(C, std(0)); |
---|
306 | |
---|
307 | dbprint(p-1, "Image Matrix: ", C); |
---|
308 | |
---|
309 | // and write it in Image basis (Imk) |
---|
310 | |
---|
311 | ideal Imk = imap(Tm, Imk); |
---|
312 | |
---|
313 | module D; poly lm; vector tt; |
---|
314 | |
---|
315 | for( i = ncols(C); i > 0; i-- ) |
---|
316 | { |
---|
317 | t = C[i]; |
---|
318 | tt = 0; |
---|
319 | |
---|
320 | while( t != 0 ) |
---|
321 | { |
---|
322 | lm = leadmonom(t); |
---|
323 | // lm; |
---|
324 | for( j = ncols(Imk); j > 0; j-- ) |
---|
325 | { |
---|
326 | if( lm / Imk[j] != 0 ) |
---|
327 | { |
---|
328 | tt = tt + (lead(t) / Imk[j]) * gen(j); |
---|
329 | break; |
---|
330 | } |
---|
331 | } |
---|
332 | t = t - lead(t); |
---|
333 | } |
---|
334 | |
---|
335 | D[i] = tt; |
---|
336 | // tt; |
---|
337 | } |
---|
338 | |
---|
339 | // D; print(D); |
---|
340 | |
---|
341 | |
---|
342 | setring save; |
---|
343 | |
---|
344 | // basering; |
---|
345 | |
---|
346 | module AA = imap(Rm, D); |
---|
347 | |
---|
348 | // nrows(AA) = Nk; // ???? |
---|
349 | // ncols(AA) = Mk; |
---|
350 | |
---|
351 | // Nk, Mk; |
---|
352 | |
---|
353 | // AA[Mk] = AA[Mk] * 1; |
---|
354 | |
---|
355 | return(AA); |
---|
356 | } |
---|
357 | example |
---|
358 | { "EXAMPLE:"; echo = 2; |
---|
359 | int save = printlevel; printlevel = 1; |
---|
360 | |
---|
361 | ring r = (0),(a, b, c, d), dp; r; |
---|
362 | |
---|
363 | module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1); |
---|
364 | |
---|
365 | print(exteriorPower(A, 2)); |
---|
366 | print(exteriorPower(A, 3)); |
---|
367 | |
---|
368 | module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2); |
---|
369 | |
---|
370 | print(exteriorPower(B, 2)); |
---|
371 | print(exteriorPower(B, 3)); |
---|
372 | |
---|
373 | kill r; |
---|
374 | ring r = (0),(a, b, c, d), dp;def g = SuperCommutative();setring g; |
---|
375 | |
---|
376 | module A = a*gen(1), b * gen(1), c*gen(1), d * gen(1); |
---|
377 | |
---|
378 | print(exteriorPower(A, 2)); |
---|
379 | print(exteriorPower(A, 3)); |
---|
380 | |
---|
381 | module B = a*gen(1) + c* gen(2), b * gen(1) + d * gen(2); |
---|
382 | |
---|
383 | print(exteriorPower(B, 2)); |
---|
384 | print(exteriorPower(B, 3)); |
---|
385 | |
---|
386 | kill r; |
---|
387 | |
---|
388 | printlevel = save; |
---|
389 | } |
---|