1 /***********************************************************************/
2 /***********************************************************************/
3 /* */
4 /* Programm zur Berechnung einer symmetrieangepassten Zerlegung */
5 /* eines durch eine Matrix gegebenen Operators M. */
6 /* Dieser Operator habe die symmetrien einer endlichen */
7 /* Permutationsgruppe G. Es wird die symmetrieangepasste Basis */
8 /* in der Matrix B berechnet und der Basiswechsel durchgefuehrt. */
9 /* Nach dem Programmlauf steht in M ein Vektor von Matrizen, */
10 /* die die Diagonalbloecke widerspiegeln. Im Vektor vf stehen deren */
11 /* Vielfachheiten. */
12 /* */
13 /* Written by: Ralf Hager September 1991 */
14 /* */
15 /***********************************************************************/
16 /***********************************************************************/
17
18 #include"def.h"
19 #include"macro.h"
20
21
22 #define SAB_get_index index_vector /* AK 010393 */
23 #define war_schon_da(a,b) ((index_vector(a,b) == -1L ) ? 1L : 0L)
24
25 INT glm_get_BV();
26 INT bilde_htupel();
27 INT werte_Polynom_aus();
28 INT bestimme_D();
29 INT glm_sab();
30 INT _homtest();
31 INT bestimme_zufallsmatrizen();
32 INT glm_B_W();
33 INT reell_gram_schmidt();
34 INT B_W();
35 INT vec_mat_mult();
36 static INT _bdg();
37 static INT _xdg();
38 static INT glmn_B_W();
39 static INT get_poly_self();
40 static INT copy_V_in_B();
41 static INT get_vectors();
42 static INT get_BV();
43 static INT add_to_PROJ();
44 static INT mhochexpo();
45 static INT TP();
46 static INT get_nr_of_tupel();
47 static INT lls_vgl();
48 static INT ti_eq_tk();
49 static INT get_axial_distance();
50 static INT ti_eq_transpo_tk();
51 static INT get_main_diag();
52 static INT D_calc();
53 static INT D_row_calc();
54 static INT transponiere();
55 static INT get_T();
56 static INT get_Q_perm();
57 static INT get_H_perm();
58 static INT cut_col_row();
59 static INT ziffern_existieren();
60 static INT Pcut_col_row();
61 static INT createP();
62 static INT lex_vgl();
63 static INT sortieren();
64 static INT zuweisen();
65 static INT _del();
66 static INT _ins();
67 static INT _kt();
68 static INT get_col();
69 static INT get_row();
70 static INT Pziffern_existieren();
71 static INT alpha_beta_tabs();
72 static INT get_submatrix();
73 static INT _sab();
74 /***********************************************************************/
75 /* */
76 /* Routine: sab */
77 /* Diese Routine ruft zunaechst _sab auf, wo die Basis B berechnet */
78 /* wird. Im Unterprogramm B_W wird dann der Basiswechsel durchgefuehrt.*/
79 /* */
80 /***********************************************************************/
81 /* RH 011091 */
82
83 #ifdef SABTRUE
sab(D_i,D,basis,M,sn_dim)84 INT sab(D_i,D,basis,M,sn_dim) OP D_i,D,basis,M,sn_dim;
85 {
86 OP sym_dim;
87 INT erg = OK;
88 CTO(VECTOR,"sab(2)",D);
89
90 sym_dim = callocobject();
91
92 erg += _sab(D_i,D,basis,sym_dim,sn_dim);
93 erg += B_W(basis,M,sym_dim,sn_dim);
94
95 erg += freeall(sym_dim);
96 ENDR("sab");
97 }
98 #endif /* SABTRUE */
99
100
101 /***********************************************************************/
102 /* */
103 /* Routine: _sab */
104 /* Die symmetrieangepasste Basis wird berechnet und in B abgespeichert.*/
105 /* */
106 /***********************************************************************/
107 /* RH 011091 */
108
109 #ifdef SABTRUE
_sab(D_i,D,basis,sym_dim,sn_dim)110 static INT _sab(D_i,D,basis,sym_dim,sn_dim) OP D_i, D, basis, sym_dim, sn_dim;
111 {
112 INT i;
113 INT j;
114 INT k;
115 INT l;
116 INT s;
117 INT merk = 0L;
118
119 OP count = callocobject();
120 OP G = callocobject();
121 OP N = callocobject();
122 OP f_D = callocobject();
123 OP irr_dim = callocobject();
124 OP faktor = callocobject();
125 OP dim = callocobject();
126 OP PROJ = callocobject();
127 OP V;
128 OP X;
129
130 m_i_i(0L,count);
131 m_i_i(0L,dim);
132 m_i_i(S_V_LI(D),G);
133 m_i_i(S_P_LI(S_V_I(D,0L)),f_D);
134 m_i_i(S_V_LI(D_i),N);
135 m_lh_nm(f_D,f_D,basis);
136 m_l_nv(N,sym_dim);
137 m_l_nv(N,sn_dim);
138 m_lh_nm(f_D,f_D,PROJ);
139
140 for(i=0L;i<S_I_I(N);++i)
141 {
142 for(l=0;l<S_I_I(f_D);++l)
143 for(s=0;s<S_I_I(f_D);++s)
144 m_i_i(0L,S_M_IJ(PROJ,l,s));
145
146 m_i_i(S_V_LI(S_V_I(S_V_I(D_i,i),0L)),irr_dim);
147 for(j=0L;j<S_I_I(G);++j)
148 {
149 copy(S_V_I(S_V_I(S_V_I(D_i,i),j),0L),faktor);
150 if(!nullp(faktor))
151 add_to_PROJ(PROJ,faktor,S_V_I(D,j));
152 }
153 copy(count,dim);
154 get_BV(PROJ,basis,count);
155 sub(count,dim,dim);
156 merk+= S_I_I(dim)*S_I_I(irr_dim);
157 if(S_I_I(dim) != 0L)
158 {
159 V = callocobject();
160 m_i_i(S_I_I(dim),S_V_I(sym_dim,i));
161 m_i_i(S_I_I(irr_dim),S_V_I(sn_dim,i));
162 get_vectors(basis,count,dim,V);
163 copy_V_in_B(V,basis,count,dim);
164 for(k=1L;k<S_I_I(irr_dim);++k)
165 {
166 for(l=0;l<S_I_I(f_D);++l)
167 for(s=0;s<S_I_I(f_D);++s)
168 m_i_i(0L,S_M_IJ(PROJ,l,s));
169 for(j=0L;j<S_I_I(G);++j)
170 {
171 copy(S_V_I(S_V_I(S_V_I(D_i,i),j),k),
172 faktor);
173 if(!nullp(faktor))
174 {
175 add_to_PROJ(PROJ,
176 faktor,
177 S_V_I(D,j));
178 }
179 }
180 copy(G,faktor);
181 div(irr_dim,faktor,faktor);
182 mult_apply(faktor,PROJ); /* AK 020692 statt mult */
183
184 X = callocobject();
185 mult(PROJ,V,X);
186 m_i_i(S_I_I(count)+S_I_I(dim),count);
187 copy_V_in_B(X,basis,count,dim);
188 freeall(X);
189 }
190 freeall(V);
191 }
192 }
193 freeall(count);
194 freeall(dim);
195 freeall(G);
196 freeall(N);
197 freeall(f_D);
198 freeall(irr_dim);
199 freeall(faktor);
200 freeall(PROJ);
201 return OK;
202 }
203 #endif /* SABTRUE */
204
205 /***********************************************************************/
206 /* */
207 /* Routine: add_to_PROJ */
208 /* Zur Projektormatrix PROJ wird d_1i(pi^-1)D(pi) dazuaddiert. */
209 /* */
210 /***********************************************************************/
211 /* RH 011091 */
212
add_to_PROJ(PROJ,d,D)213 static INT add_to_PROJ(PROJ,d,D) OP PROJ, d, D;
214 {
215 INT i;
216
217 for(i=0L;i<S_P_LI(D);++i)
218 add_apply( d, S_M_IJ(PROJ,S_P_II(D,i)-1L,i));
219 return OK;
220 }
221
222
223 /***********************************************************************/
224 /* */
225 /* Routine: get_BV */
226 /* Aus der Projektormatrix PROJ werden die linear unabhaengigen */
227 /* Vektoren durch ein Gaussverfahren herausgesucht. */
228 /* In count wird deren Anzahl abgespeichert. Sie werden in B */
229 /* uebertragen. */
230 /* */
231 /***********************************************************************/
232 /* RH 011091 */
233
get_BV(PROJ,basis,count)234 static INT get_BV(PROJ,basis,count) OP PROJ; OP basis; OP count;
235 {
236 INT i;
237 INT j;
238 INT k;
239 INT l;
240 INT s;
241 INT f_D;
242
243 OP Proj = callocobject();
244 OP val = callocobject();
245 OP faktor = callocobject();
246
247 copy(PROJ,Proj);
248 f_D = S_M_LI(basis);
249 k = 0L;
250 for(i=0L;i<f_D;++i)
251 {
252 for(j=k;j<f_D;++j)
253 {
254 if(!nullp(S_M_IJ(Proj,j,i)))
255 {
256 for(l=0L;l<f_D;++l)
257 copy(S_M_IJ(PROJ,l,i),
258 S_M_IJ(basis,l,S_I_I(count)));
259 m_i_i(S_I_I(count)+1L,count);
260 if(j != k)
261 for(l=0L;l<f_D;++l)
262 swap(S_M_IJ(Proj,k,l),
263 S_M_IJ(Proj,j,l));
264 for(l=k+1L;l<f_D;++l)
265 if(!nullp(S_M_IJ(Proj,l,i)))
266 {
267 div(S_M_IJ(Proj,l,i),
268 S_M_IJ(Proj,k,i),
269 faktor);
270 for(s=i;s<f_D;++s)
271 {
272 mult(faktor,
273 S_M_IJ(Proj,k,s),
274 val);
275 sub(S_M_IJ(Proj,l,s),
276 val,
277 S_M_IJ(Proj,l,s));
278 }
279 }
280 }
281 }
282 k++;
283 }
284
285 freeall(Proj);
286 freeall(val);
287 freeall(faktor);
288 return OK;
289 }
290
291
292 /***********************************************************************/
293 /* */
294 /* Routine: get_vectors */
295 /* Holt aus basis dim Vektoren und speichert sie in V ab. */
296 /* */
297 /***********************************************************************/
298 /* RH 011091 */
299
300 #ifdef MATRIXTRUE
get_vectors(basis,count,dim,V)301 static INT get_vectors(basis,count,dim,V) OP basis, count, dim, V;
302 {
303 INT i;
304 INT j;
305
306 m_lh_nm(dim,S_M_L(basis),V);
307
308 for(i=0L;i<S_M_LI(basis);++i)
309 for(j=0L;j<S_I_I(dim);++j)
310 copy(S_M_IJ(basis,i,S_I_I(count)-S_I_I(dim)+j),
311 S_M_IJ(V,i,j));
312
313 return OK;
314 }
315 #endif /* MATRIXTRUE */
316
317 /***********************************************************************/
318 /* */
319 /* Routine: copy_V_in_B */
320 /* V wird auf basis kopiert an die Stelle count -dim bis count-1L */
321 /* */
322 /***********************************************************************/
323 /* RH 011091 */
324
copy_V_in_B(V,basis,count,dim)325 static INT copy_V_in_B(V,basis,count,dim) OP V,basis,count,dim;
326 {
327 INT i;
328 INT j;
329
330 for(i=0L;i<S_M_HI(V);++i)
331 for(j=0L;j<S_M_LI(V);++j)
332 copy(S_M_IJ(V,i,j),
333 S_M_IJ(basis,i,S_I_I(count)-S_I_I(dim)+j));
334 return OK;
335 }
336
337 /***********************************************************************/
338 /* */
339 /* Routine: gram_schmidt */
340 /* Nach dem Gram-Schmidtschen-Orthonormalisierungsverfahren wird eine */
341 /* Matrix, bestehend aus l.u. Vektoren orthonormalisiert. */
342 /* */
343 /***********************************************************************/
344 /* RH 011091 */
345
346 #ifdef MATRIXTRUE
gram_schmidt(a)347 INT gram_schmidt(a) OP a;
348 {
349 INT i;
350 INT j;
351 INT k;
352
353 OP r = callocobject();
354 OP h = callocobject();
355 OP n = callocobject();
356 OP p = callocobject();
357 OP s = callocobject();
358
359 m_lh_nm(S_M_H(a),S_M_L(a),r);
360 m_i_i(S_M_HI(a),n);
361 m_i_i(S_M_LI(a),p);
362
363 for(k=0L;k<S_I_I(p);++k)
364 {
365 for(j=0L;j<k;++j)
366 {
367 m_i_i(0L,S_M_IJ(r,j,k));
368 for(i=0L;i<S_I_I(n);++i)
369 {
370 mult(S_M_IJ(a,i,j),S_M_IJ(a,i,k),h);
371 add_apply(h,S_M_IJ(r,j,k));
372 }
373 for(i=0L;i<S_I_I(n);++i)
374 {
375 mult(S_M_IJ(r,j,k),S_M_IJ(a,i,j),h);
376 sub(S_M_IJ(a,i,k),h,S_M_IJ(a,i,k));
377 }
378 }
379 m_i_i(0L,s);
380 for(i=0L;i<S_I_I(n);++i)
381 {
382 complex_conjugate(S_M_IJ(a,i,k),h);
383 mult_apply(S_M_IJ(a,i,k),h); /* AK 020692 statt mult */
384 add_apply(h,s); /* AK 010692 statt add */
385 }
386 squareroot(s,S_M_IJ(r,k,k));
387 for(i=0L;i<S_I_I(n);++i)
388 div(S_M_IJ(a,i,k),S_M_IJ(r,k,k),S_M_IJ(a,i,k));
389 }
390
391 freeall(r);
392 freeall(h);
393 freeall(n);
394 freeall(p);
395 freeall(s);
396 return OK;
397 }
398 #endif /* MATRIXTRUE */
399
400
401
402 /***********************************************************************/
403 /* */
404 /* Routine: group_gen */
405 /* Ausgehend von einem Erzeugendensystem S von Permutationen */
406 /* wir die Gruppe G erzeugt und in D abgespeichert. In SMat sind die */
407 /* darstellenden irreduziblen Matrizen abgespeichert. Mit ihnen */
408 /* werden die irreduziblen Matrixdarstellungen saemtlicher */
409 /* Gruppenelemente erzeugt und in Di abgespeichert. */
410 /* Der Algorithmus ist eine vereinfachte Version des Dimino-Algorithmus*/
411 /* */
412 /***********************************************************************/
413 /* RH 011091 */
414
group_gen(S,SMat,D,Di)415 INT group_gen(S,SMat,D,Di) OP S,SMat,D,Di;
416 {
417 INT i;
418 INT j;
419 INT k;
420
421 OP irr_dim = callocobject();
422 OP hperm = callocobject();
423 OP hvek = callocobject();
424 OP besucht = callocobject();
425
426 m_il_v(S_V_LI(S)+1L,D);
427 m_il_nv(S_P_LI(S_V_I(S,0L)),S_V_I(D,0L));
428 first_permutation(S_P_L(S_V_I(S,0L)),S_V_I(D,0L));
429 for(i=1L;i<S_V_LI(D);++i)
430 copy(S_V_I(S,i-1L),S_V_I(D,i));
431 m_l_v(S_V_L(SMat),Di);
432 for(i=0L;i<S_V_LI(SMat);++i)
433 {
434 m_il_v(S_V_LI(S_V_I(SMat,i))+1L,S_V_I(Di,i));
435 m_i_i(S_M_LI(S_V_I(S_V_I(SMat,i),0L)),irr_dim);
436 m_l_nv(irr_dim,S_V_I(S_V_I(Di,i),0L));
437 m_i_i(1L,S_V_I(S_V_I(S_V_I(Di,i),0L),0L));
438 for(j=1L;j<S_V_LI(S_V_I(Di,i));++j)
439 {
440 m_i_i(S_M_LI(S_V_I(S_V_I(SMat,i),0L)),irr_dim);
441 m_l_v(irr_dim,S_V_I(S_V_I(Di,i),j));
442 for(k=0L;k<S_I_I(irr_dim);++k)
443 copy(S_M_IJ(S_V_I(S_V_I(SMat,i),j-1L),0L,k),
444 S_V_I(S_V_I(S_V_I(Di,i),j),k));
445 }
446 }
447
448 for(i=0L;i<S_V_LI(D);++i)
449 {
450 for(j=0L;j<S_V_LI(S);++j)
451 {
452 mult(S_V_I(D,i),S_V_I(S,j),hperm);
453 if(war_schon_da(hperm,D) != 0L)
454 {
455 inc(D);
456 copy(hperm,S_V_I(D,S_V_LI(D)-1L));
457 for(k=0L;k<S_V_LI(Di);++k)
458 {
459 inc(S_V_I(Di,k));
460 m_l_nv(S_M_L(S_V_I(S_V_I(SMat,k),0L)),
461 hvek);
462 vec_mat_mult(S_V_I(S_V_I(Di,k),i),
463 S_V_I(S_V_I(SMat,k),j),
464 hvek);
465 copy(hvek,
466 S_V_I(S_V_I(Di,k),
467 S_V_LI(S_V_I(Di,k))-1L));
468 }
469 }
470 }
471 }
472
473 m_l_nv(S_V_L(D),besucht);
474 for(i=0L;i<S_V_LI(D);++i)
475 {
476 m_i_i(1L,S_V_I(besucht,i));
477 invers(S_V_I(D,i),hperm);
478 k = SAB_get_index(hperm,D);
479 if(k>=0L && k != i && S_V_II(besucht,k) == 0L)
480 {
481 m_i_i(1L,S_V_I(besucht,k));
482 for(j=0L;j<S_V_LI(Di);++j)
483 swap(S_V_I(S_V_I(Di,j),i),
484 S_V_I(S_V_I(Di,j),k));
485 }
486 }
487
488 freeall(irr_dim);
489 freeall(hperm);
490 freeall(hvek);
491 freeall(besucht);
492 return OK;
493 }
494
495 #ifdef UNDEF
496 /***********************************************************************/
497 /* */
498 /* Routine: SAB_get_index */
499 /* Hilfsroutine bei der Konstruktion der Gruppe. */
500 /* */
501 /***********************************************************************/
502 /* RH 011091 */
503
SAB_get_index(hperm,D)504 static INT SAB_get_index(hperm,D) OP hperm; OP D;
505 {
506 INT i;
507
508 for(i=0L;i< S_V_LI(D);++i)
509 if(comp(hperm,S_V_I(D,i)) == 0) return(i);
510 return(-1);
511 }
512 #endif /* UNDEF */
513
514 /***********************************************************************/
515 /* */
516 /* Routine: vec_mat_mult */
517 /* Hilfsroutine bei der Konstruktion der Gruppe. */
518 /* */
519 /***********************************************************************/
520 /* RH 011091 */
521
vec_mat_mult(vec,mat,hvek)522 INT vec_mat_mult(vec,mat,hvek) OP vec,mat,hvek;
523 {
524 INT i;
525 INT j;
526
527 OP z = callocobject();
528
529 for(i=0L;i<S_V_LI(vec);++i)
530 {
531 for(j=0L;j<S_V_LI(vec);++j)
532 {
533 mult(s_v_i(vec,j),S_M_IJ(mat,j,i),z);
534 add_apply(z,S_V_I(hvek,i)); /* AK 010692 statt add */
535 }
536 }
537
538 freeall(z);
539 return OK;
540 }
541
542 #ifdef UNDEF
543 /***********************************************************************/
544 /* */
545 /* Routine: war_schon_da */
546 /* Hilfsroutine bei der Konstruktion der Gruppe. */
547 /* */
548 /***********************************************************************/
549 /* RH 011091 */
550
war_schon_da(hperm,D)551 war_schon_da(hperm,D) OP hperm,D;
552 {
553 INT i;
554
555 for(i=0L;i< S_V_LI(D);++i)
556 if(comp(hperm,S_V_I(D,i)) == 0) return(0L);
557 return(1);
558 }
559 #endif
560
561 /***********************************************************************/
562 /* */
563 /* Routine: B_W */
564 /* Der Basiswechsel von M wird durchgefuerht. Dazu wird zunaechst die */
565 /* Inverse Matrix basis^-1 von basis bestimmt. Anschliessend wird die */
566 /* Matrixmultiplikation durchgefuehrt, unter Beruecksichtigung */
567 /* der Vielfachheiten der vorkommenden Darstellungen. */
568 /* */
569 /***********************************************************************/
570 /* RH 011091 */
571 #ifdef SABTRUE
B_W(basis,M,sym_dim,sn_dim)572 INT B_W(basis,M,sym_dim,sn_dim) OP M,basis,sym_dim,sn_dim;
573 {
574 INT i,j,k,l,s;
575
576 OP B_1 = callocobject();
577 OP erg = callocobject();
578 OP mneu = callocobject();
579 OP count = callocobject();
580 OP Mneu = callocobject();
581 OP r = callocobject();
582 OP vf = callocobject();
583
584 invers(basis,B_1);
585
586 m_i_i(0L,count);
587 for(s=0L;s<S_V_LI(sym_dim);++s)
588 if(S_V_II(sym_dim,s) != 0L)
589 inc(count);
590 m_l_v(count,vf);
591 m_i_i(0L,count);
592 for(s=0L;s<S_V_LI(sym_dim);++s)
593 if(S_V_II(sym_dim,s) != 0L)
594 {
595 copy(S_V_I(sn_dim,s),S_V_I(vf,S_I_I(count)));
596 inc(count);
597 }
598 m_l_v(count,Mneu);
599 m_i_i(0L,count);
600 for(i=0L;i<S_V_LI(sym_dim);++i)
601 if(S_V_II(sym_dim,i) != 0L)
602 {
603 m_lh_m(S_V_I(sym_dim,i),S_V_I(sym_dim,i),S_V_I(Mneu,S_I_I(count)));
604 inc(count);
605 }
606 m_i_i(0L,count);
607 m_i_i(0L,r);
608 for(s=0L;s<S_V_LI(sym_dim);++s)
609 {
610 if(S_V_II(sym_dim,s) != 0L)
611 {
612 for(i=S_I_I(count);
613 i<S_I_I(count)+S_V_II(sym_dim,s);
614 ++i)
615 {
616 for(j=S_I_I(count);
617 j<S_I_I(count)+S_V_II(sym_dim,s);
618 ++j)
619 {
620 m_i_i(0L,mneu);
621 for(k=0;k<S_M_LI(basis);++k)
622 {
623 for(l=0L;l<S_M_LI(basis);++l)
624 {
625 if(!nullp(S_M_IJ(B_1,i,k)) &&
626 !nullp(S_M_IJ(basis,l,j)) &&
627 !nullp(S_M_IJ(M,k,l)))
628 {
629
630 mult(S_M_IJ(B_1,i,k),
631 S_M_IJ(basis,l,j),
632 erg);
633 mult(erg,
634 S_M_IJ(M,k,l),
635 erg);
636 add_apply(erg,mneu);
637 }
638 }
639 }
640 copy(mneu,S_M_IJ(S_V_I(Mneu,S_I_I(r)),i-S_I_I(count),j-S_I_I(count)));
641 }
642 }
643 inc(r);
644 }
645 mult(S_V_I(sym_dim,s),S_V_I(sn_dim,s),erg);
646 add_apply(erg,count);
647 }
648 copy(Mneu,M);
649 copy(vf,sn_dim);
650
651
652 freeall(B_1);
653 freeall(count);
654 freeall(erg);
655 freeall(mneu);
656 freeall(Mneu);
657 freeall(r);
658 freeall(vf);
659 return OK;
660 }
661 #endif /* SABTRUE */
662
663 /***********************************************************************/
664 /* */
665 /* Routine: sab_input */
666 /* Von der standardeingabe werden folgende Strukturen eingelesen: */
667 /* -------------------------------------------------------------- */
668 /* Anzahl der Erzeuger von G | orderS */
669 /* Erzeuger von G | S */
670 /* Anzahl irred. Darstellungen | anz_irr */
671 /* Darstellende Matrizen fuer diese | SMat */
672 /* Zu zerlegender Operator M | M */
673 /* -------------------------------------------------------------- */
674 /* */
675 /***********************************************************************/
676 /* RH 011091 */
677
sab_input(S,SMat,M)678 INT sab_input(S,SMat,M) OP S; OP SMat; OP M;
679 {
680 INT i;
681 INT j;
682
683 OP orderS = callocobject();
684 OP anz_irr = callocobject();
685
686 scan(INTEGER,orderS);
687 m_l_v(orderS,S);
688 for(i=0L;i<S_I_I(orderS);++i)
689 scan(PERMUTATION,S_V_I(S,i));
690 scan(INTEGER,anz_irr);
691
692 m_l_v(anz_irr,SMat);
693 for(i=0L;i<S_I_I(anz_irr);++i)
694 {
695 m_l_v(orderS,S_V_I(SMat,i));
696 for(j=0;j<S_I_I(orderS);++j)
697 scan(MATRIX,S_V_I(S_V_I(SMat,i),j));
698 }
699 scan(MATRIX,M);
700
701 freeall(orderS);
702 freeall(anz_irr);
703 return OK;
704 }
705
706 /***********************************************************************/
707 /***********************************************************************/
708 /* */
709 /* Programm zur Berechnung der natuerlichen Matrixdarstellung von S_n */
710 /* */
711 /* Written by: Ralf Hager September 1991 */
712 /* */
713 /***********************************************************************/
714 /***********************************************************************/
715
716
717 /***********************************************************************/
718 /* */
719 /* Routine: bdg */
720 /* Gibt zu gegebener Partition part und Permutation perm die */
721 /* zugehoerige irreduzible Matrixdarstellung in D zurueck. */
722 /* */
723 /***********************************************************************/
724 /* RH 011091 */
725
726 #ifdef DGTRUE
bdg(part,perm,D)727 INT bdg(part,perm,D) OP part; OP perm; OP D;
728 {
729 INT i,j;
730 OP t,P,_D;
731 INT erg = OK;
732
733 CTO(PARTITION,"bdg(1)",part);
734 CTO(PERMUTATION,"bdg(2)",perm);
735 CE3(part,perm,D,bdg);
736
737 t = CALLOCOBJECT();
738 P = CALLOCOBJECT();
739 _D = CALLOCOBJECT();
740
741 _bdg(part,perm,_D,t,P,-1L);
742 /* Matrix _D wird als Vektor von Vektoren berechnet.
743 Es ist auch der Aufruf mit var ungleich -1L erlaubt,
744 dann wird die var-te Zeile der darstellenden Matrix
745 berechnet. */
746
747 erg += m_ilih_m(S_V_LI(_D),S_V_LI(_D),D);
748 /* Der Vektor von Vektoren wird auf eine Matrix kopiert. */
749
750 for(i=0L;i<S_V_LI(_D);++i)
751 for(j=0L;j<S_M_LI(_D);++j)
752 COPY(S_V_I(S_V_I(_D,i),j),
753 S_M_IJ(D,i,j));
754
755 FREEALL3(t,P,_D);
756 ENDR("bdg");
757 }
758
759
760
761 /***********************************************************************/
762 /* */
763 /* Routine: _bdg */
764 /* Gibt zu gegebener Partition part und Permutation perm die */
765 /* zugehoerige irreduzible Matrixdarstellung in D zurueck. */
766 /* (D ist aber Vektor von Vektoren.) */
767 /* */
768 /***********************************************************************/
769 /* RH 011091 */
770
_bdg(part,perm,D,t,P,var)771 static INT _bdg(part,perm,D,t,P,var) OP part,perm,D,t,P; INT var;
772 {
773 INT i;
774 INT j;
775 INT lex_vgl();
776
777 OP p_inv = callocobject(); /* inverse Permutatione perm */
778 OP tpart = callocobject(); /* transponierte Partition von part */
779 OP beta = callocobject(); /* Inhalt der Tableaux */
780 OP Tt = callocobject(); /* Tableau mit transponiertem Umriss*/
781
782 if(emptyp(t))
783 {
784 m_il_v(S_P_LI(perm),beta);
785 for(i=0L;i<S_V_LI(beta);++i)
786 m_i_i(1L,S_V_I(beta,i));
787 /* Konstruktion der standardtableaux */
788 alpha_beta_tabs(part,beta,t);
789 /* Lex. aufst. Sortieren der Tabs */
790 sortieren(0L,S_V_LI(t)-1L,t,lex_vgl);
791 /* Struktur der Horizontal-
792 permutationen wird gebildet */
793 createP(P,S_P_L(perm),part,t);
794 }
795 conjugate(part,tpart);
796 m_il_v(S_PA_LI(tpart),Tt);
797 for(i=S_PA_LI(tpart)-1L;i>=0L;--i)
798 m_il_v(S_PA_II(tpart,i),
799 S_V_I(Tt,S_PA_LI(tpart)-1L-i));
800
801 invers(perm,p_inv); invers(perm,perm);
802 m_il_v(S_V_LI(t),D);
803 if(var == -1L) /* Gesamte Matrix berechnen */
804 {
805 for(i=0;i<S_V_LI(t);++i)
806 {
807 m_il_v(S_V_LI(t),S_V_I(D,i));
808 for(j=0L;j<S_V_LI(D);++j)
809 m_i_i(0L,S_V_I(S_V_I(D,i),j));
810 /* Die i-te Zeile der darstellenden Matrix wird
811 bestimmt */
812 D_calc(t,S_V_I(D,i),P,perm,p_inv,1L,i,i,i,Tt);
813 }
814 }
815 else /* NUR die var-te Zeile berechnen */
816 {
817 if(var < S_V_LI(D) && var >= 0L)
818 {
819 for(j=0L;j<S_V_LI(D);++j)m_i_i(0L,S_V_I(D,j));
820 D_calc(t,D,P,perm,p_inv,1L,var,var,var,Tt);
821 }
822 else return error("wrong linenr. !!!\n");
823 }
824 invers(perm,perm);
825
826 freeall(p_inv);
827 freeall(tpart);
828 freeall(Tt);
829 freeall(beta);
830 return OK;
831 }
832 #endif /* DGTRUE */
833 /***********************************************************************/
834 /* */
835 /* Routine: alpha_beta_tabs */
836 /* Programm zur Berechnung der Standardtableaux mit vorgegebenem */
837 /* Umriss alpha und Inhalt beta. */
838 /* */
839 /***********************************************************************/
840 /* RH 011091 */
841
842 #ifdef TABLEAUXTRUE
alpha_beta_tabs(alpha,beta,t)843 static INT alpha_beta_tabs(alpha,beta,t) OP alpha,beta,t;
844 {
845 OP n = callocobject();
846 OP h = callocobject();
847 OP ou = callocobject();
848 OP um = callocobject();
849 OP tou = callocobject();
850 OP ziel = callocobject();
851 OP dim = callocobject();
852 OP tab = callocobject();
853
854 INT i;
855 INT j;
856 INT ct = 0L;
857
858 /*
859 kostka_zahl(alpha,beta,dim);
860 */
861 kostka_number(beta,alpha,dim);
862 m_l_v(dim,t);
863 conjugate(alpha,h);
864
865
866 m_il_v(S_PA_LI(alpha)+1L,ou);
867 m_il_v(S_PA_LI(h)+1L,tou);
868 m_il_v(S_PA_LI(h)+1L,ziel);
869
870
871 for(i=1L;i<=S_PA_LI(alpha);++i)
872 m_i_i(S_PA_II(alpha,S_PA_LI(alpha)-i),S_V_I(ou,i));
873
874 for(i=1L;i<=S_PA_LI(h);++i)
875 m_i_i(S_PA_II(h,S_PA_LI(h)-i),S_V_I(tou,i));
876
877 m_il_v(S_V_LI(tou),tab);
878 for(i=1L;i<S_V_LI(tou);++i)
879 m_il_v(S_V_II(tou,i)+1L,S_V_I(tab,i));
880
881 for(i=1L;i<S_V_LI(tou);++i)
882 for(j=0L;j<=S_V_II(tou,i);++j)
883 m_i_i(0L,S_V_I(S_V_I(tab,i),j));
884
885 for(i=0L;i<S_I_I(dim);++i)
886 {
887 m_il_v(S_V_LI(ou)-1L,S_V_I(t,i));
888 for(j=0L;j<S_V_LI(ou)-1L;++j)
889 m_il_v(S_V_II(ou,j+1L),S_V_I(S_V_I(t,i),j));
890 }
891 copy(tou,um);
892 for(i=1L;i<S_V_LI(tou);++i)
893 m_i_i(S_V_II(um,1L)+i-1L-S_V_II(um,i),S_V_I(ziel,i));
894 m_i_i(-1L,S_V_I(um,0L));
895 for(i=0L;i<S_V_LI(um)-1L;++i)
896 m_i_i(S_V_II(um,1L)+i,S_V_I(um,i+1L));
897
898 m_i_i(0L,n);
899 m_i_i(0L,S_V_I(ziel,0L));
900 m_i_i(0L,S_V_I(ou,0L));
901 m_i_i(0L,S_V_I(tou,0L));
902
903 for(i=0L;i<S_V_LI(tou);++i)
904 add_apply(S_V_I(tou,i),n); /* AK 010692 statt add */
905
906 /* Rekursive Konstruktionsroutine */
907
908 _kt(n,t,tab,um,ziel,ou,tou,0L,0L,1L,&ct,
909 beta,S_V_II(beta,0L));
910
911
912 freeall(ou);
913 freeall(um);
914 freeall(tou);
915 freeall(ziel);
916 freeall(h);
917 freeall(n);
918 freeall(dim);
919 freeall(tab);
920 return OK;
921 }
922
923 /***********************************************************************/
924 /* */
925 /* Routine: _kt */
926 /* Rekursive Konstruktionsroutine zur Berechnung der */
927 /* Standardtableaux mit vorgegebenem */
928 /* Umriss alpha und Inhalt beta. */
929 /* Zum Algorithmus siehe [KO07], (S.8) */
930 /* */
931 /***********************************************************************/
932 /* RH 011091 */
_kt(n,t,tab,um,ziel,ou,tou,k,i,st,ct,beta,nr)933 static INT _kt(n,t,tab,um,ziel,ou,tou,k,i,st,ct,beta,nr)
934 OP n; OP t; OP tab; OP um; OP ziel;
935 OP ou; OP tou; INT k; INT i; INT st;
936 INT *ct; OP beta; INT nr;
937 {
938 INT l;
939
940 if(i==nr)
941 {
942 if(st==S_V_LI(beta))
943 {
944 zuweisen(t,tab,ou,*ct);
945 (*ct)++;
946 }
947 else _kt(n,t,tab,um,ziel,ou,tou,0L,0L,st+1L,ct,
948 beta,S_V_II(beta,st));
949 }
950 else
951 for(l=k+1L;l<=S_V_II(ou,1L);++l)
952 {
953 if(((S_V_II(um,l)-1L) > (S_V_II(um,l-1L)))&&
954 (S_V_II(um,l) > S_V_II(ziel,l)))
955 {
956 m_i_i(S_V_II(um,l)-1L,S_V_I(um,l));
957 _ins(S_V_I(tab,l),st);
958 _kt(n,t,tab,um,ziel,ou,tou,l,i+1L,st,ct,
959 beta,nr);
960 _del(S_V_I(tab,l));
961 m_i_i(S_V_II(um,l)+1L,S_V_I(um,l));
962 }
963 }
964 return OK;
965 }
966
967 /***********************************************************************/
968 /* */
969 /* Routine: _ins */
970 /* Hilfsroutine bei der Tableauberechnung */
971 /* */
972 /***********************************************************************/
973 /* RH 011091 */
974
_ins(line,z)975 static INT _ins(line,z) OP line; INT z;
976 {
977 INT i;
978
979 for(i=1L;i<S_V_LI(line);++i)
980 if(S_V_II(line,i)==0L)
981 {
982 m_i_i(z,S_V_I(line,i));
983 break;
984 }
985 return OK;
986 }
987 /***********************************************************************/
988 /* */
989 /* Routine: _del */
990 /* Hilfsroutine bei der Tableauberechnung */
991 /* */
992 /***********************************************************************/
993 /* RH 011091 */
994
_del(line)995 static INT _del(line) OP line;
996 {
997 INT i;
998
999 for(i=S_V_LI(line)-1L;i>=1L;--i)
1000 if(S_V_II(line,i)!=0L)
1001 {
1002 m_i_i(0L,S_V_I(line,i));
1003 break;
1004 }
1005 return OK;
1006 }
1007 /***********************************************************************/
1008 /* */
1009 /* Routine: _zuweisen */
1010 /* Hilfsroutine bei der Tableauberechnung */
1011 /* */
1012 /***********************************************************************/
1013 /* RH 011091 */
1014
zuweisen(t,tab,ou,ind)1015 static INT zuweisen(t,tab,ou,ind)
1016 OP t; OP tab; OP ou; INT ind;
1017 {
1018 INT i;
1019 INT j;
1020
1021 for(i=0;i<S_V_LI(ou)-1L;++i)
1022 for(j=0;j<S_V_II(ou,i+1L);++j)
1023 m_i_i(S_V_II(S_V_I(tab,j+1L),i+1L),
1024 S_V_I(S_V_I(S_V_I(t,ind),i),j));
1025 return OK;
1026 }
1027 #endif /* TABLEAUXTRUE */
1028
1029 /***********************************************************************/
1030 /* */
1031 /* Routine: _sortieren */
1032 /* Sortiert einen Vektor t von Strukturen bei gegebener */
1033 /* Vergleichsfunktion _vgl nach dem Algorithmus QUICKSORT */
1034 /* */
1035 /***********************************************************************/
1036 /* RH 011091 */
1037
sortieren(anf,end,t,_vgl)1038 static INT sortieren(anf,end,t,_vgl)
1039 INT anf; INT end; OP t; INT (*_vgl)();
1040 {
1041 INT laenge;
1042 INT ind;
1043 INT i;
1044 INT l;
1045 INT r;
1046 OP mitte = callocobject();
1047
1048
1049 laenge = end - anf + 1L;
1050
1051 if(laenge >= 2L)
1052 {
1053 copy(S_V_I(t,anf),mitte);
1054 ind = -1L;
1055 i = anf + 1L;
1056
1057 while((i<=end)&&(ind == -1L))
1058 {
1059 if((*_vgl)(mitte,S_V_I(t,i)) == -1L) ind = i;
1060 else
1061 {
1062 if((*_vgl)(mitte,S_V_I(t,i)) == 1L) ind = anf;
1063 else ++i;
1064 }
1065 if(ind > -1L)
1066 {
1067 l = anf;
1068 r = end;
1069 copy(S_V_I(t,ind),mitte);
1070 do
1071 {
1072 swap(S_V_I(t,l),S_V_I(t,r));
1073 while((*_vgl)(S_V_I(t,l),mitte) == -1L)
1074 ++l;
1075 while((*_vgl)(S_V_I(t,r),mitte) >= 0L)
1076 --r;
1077 }
1078 while(l <= r);
1079 sortieren(anf,l-1L,t,_vgl);
1080 sortieren(l,end,t,_vgl);
1081 }
1082 }
1083 }
1084 freeall(mitte); /* AK 010692 */
1085 return OK;
1086 }
1087
1088 /***********************************************************************/
1089 /* */
1090 /* Routine: lex_vgl */
1091 /* lexikographisch aufsteigender Vergleich zweier Operanden */
1092 /* */
1093 /***********************************************************************/
1094 /* RH 011091 */
lex_vgl(a,b)1095 static INT lex_vgl(a,b) OP a; OP b;
1096 {
1097 INT i;
1098 INT j;
1099
1100 for(i=0L;i<S_V_LI(a);++i)
1101 {
1102 for(j=0L;j<S_V_LI(S_V_I(a,i));++j)
1103 if(S_V_II(S_V_I(a,i),j) != S_V_II(S_V_I(b,i),j))
1104 {
1105 if(S_V_II(S_V_I(a,i),j) > S_V_II(S_V_I(b,i),j))
1106 {
1107 return(1L);
1108 }
1109 if(S_V_II(S_V_I(a,i),j) < S_V_II(S_V_I(b,i),j))
1110 {
1111 return(-1L);
1112 }
1113 }
1114 }
1115 return(0L);
1116 }
1117
1118 /***********************************************************************/
1119 /* */
1120 /* Routine: createP */
1121 /* Erstellen der Stuktur P aus Horizontalpermutationen */
1122 /* nach dem depth-first-search-Verfahren */
1123 /* */
1124 /***********************************************************************/
1125 /* RH 011091 */
createP(P,n,part,t)1126 static INT createP(P,n,part,t)
1127 /* Die Struktur P, bestehend aus Horizontalpermutationen
1128 wird gebildet. */
1129 OP P; OP n; OP part; OP t;
1130 {
1131 OP HP;
1132 OP neu;
1133 OP v = callocobject();
1134 OP tpart = callocobject();
1135
1136 INT i;
1137 INT j;
1138
1139 m_il_v(S_V_LI(t),P);
1140 m_il_v(20L,v);
1141 conjugate(part,tpart);
1142
1143 for(i=0L;i<S_V_LI(t);++i)
1144 {
1145 for(j=i+1L;j<S_V_LI(t);++j)
1146 {
1147 if( Pziffern_existieren(
1148 S_V_I(t,i),S_V_I(t,j),v,tpart) == 0L)
1149 {
1150 HP = S_V_I(P,i);
1151 if (emptyp(HP)) m_il_v(1L,HP);
1152 else inc(HP);
1153 neu = S_V_I(HP,S_V_LI(HP)-1L);
1154 m_il_v(2L,neu);
1155 m_i_i(j,S_V_I(neu,0L));
1156 b_ks_p(VECTOR,callocobject(),S_V_I(neu,1L));
1157 m_il_v(S_I_I(n),S_P_S(S_V_I(neu,1L)));
1158 get_H_perm(S_V_I(t,i),S_V_I(t,j),S_V_I(neu,1L));
1159 }
1160 }
1161 }
1162 freeall(v);
1163 freeall(tpart);
1164 return OK;
1165 }
1166
1167 /***********************************************************************/
1168 /* */
1169 /* Routine: Pziffern_existieren */
1170 /* Prueft, ob in zwei Tableaux zwei Ziffern gemeinsam in einer */
1171 /* Zeile und Spalte vorkommen */
1172 /* */
1173 /***********************************************************************/
1174 /* RH 011091 */
Pziffern_existieren(t_i,t_j,v,tpart)1175 static INT Pziffern_existieren(t_i,t_j,v,tpart) OP t_i,t_j,v,tpart;
1176 {
1177 static INT i;
1178 static INT j;
1179 static INT k;
1180 static INT tl;
1181 static OP zt_i,zv;
1182
1183 for(i=0L,tl=S_PA_LI(tpart)-1L;i<S_V_LI(S_V_I(t_i,0L));++i,tl--)
1184 {
1185 if(i < S_V_LI(S_V_I(t_i,1L)))
1186 {
1187 for(j=0L;j<S_V_LI(t_i);++j)
1188 {
1189 for(k=0,zv=S_V_S(v),zt_i=S_V_S(t_i);
1190 k<S_PA_II(tpart,tl);
1191 ++k,zt_i++,zv++)
1192 m_i_i(S_V_II(zt_i,i),zv);
1193 if(Pcut_col_row(
1194 v,S_PA_II(tpart,tl) ,
1195 S_V_I(t_j,j)) >= 2L)
1196 return(1L);
1197 }
1198 }
1199 else return(0L);
1200 }
1201 return(0L);
1202 }
1203
1204 /***********************************************************************/
1205 /* */
1206 /* Routine: Pcut_col_row */
1207 /* Bildet den Schnitt einer Zeile und einer Spalte zweier Tabelaux. */
1208 /* */
1209 /***********************************************************************/
1210 /* RH 011091 */
Pcut_col_row(col,l_col,row)1211 static INT Pcut_col_row(col,l_col,row)
1212 OP col; INT l_col; OP row;
1213 {
1214 static INT erg;
1215 static INT i;
1216 static INT j;
1217 static OP zi,zj;
1218
1219 erg = 0L;
1220 for(i=0L,zi = S_V_S(col);i<l_col;++i,zi++)
1221 for(j=0L,zj = S_V_S(row);j<S_V_LI(row);++j,zj++)
1222 if(S_I_I(zi) == S_I_I(zj))
1223 {
1224 erg++;
1225 if(erg>1L) return(erg);
1226 }
1227 return(erg);
1228 }
1229
1230 /***********************************************************************/
1231 /* */
1232 /* Routine: ziffern_existieren (Variante mit transponiertem Tabl.) */
1233 /* Prueft, ob in zwei Tableaux zwei Ziffern gemeinsam in einer */
1234 /* Zeile und Spalte vorkommen */
1235 /* */
1236 /***********************************************************************/
1237 /* RH 011091 */
1238
ziffern_existieren(Tt,t_j)1239 static INT ziffern_existieren(Tt,t_j) OP Tt; OP t_j;
1240 {
1241 static INT i;
1242 static INT j;
1243
1244 for(i=0L;i<S_V_LI(Tt);++i)
1245 for(j=0L;j<S_V_LI(t_j);++j)
1246 if(S_V_LI(S_V_I(Tt,i))>0L)
1247 {
1248 if(cut_col_row(
1249 S_V_I(Tt,i),
1250 S_V_I(t_j,j)) >= 2L)
1251 return(1L);
1252 }
1253 else return(0L);
1254 return(0L);
1255 }
1256
1257 /***********************************************************************/
1258 /* */
1259 /* Routine: cut_col_row (Variante mit transponiertem Tabl.) */
1260 /* Bildet den Schnitt einer Zeile und einer Spalte zweier Tabelaux. */
1261 /* */
1262 /***********************************************************************/
1263 /* RH 011091 */
cut_col_row(col,row)1264 static INT cut_col_row(col,row) OP col; OP row;
1265 {
1266 static INT erg;
1267 static INT i;
1268 static INT j;
1269 static OP zi,zj;
1270
1271 erg = 0L;
1272 for(i=0L,zi = S_V_S(col);i<S_V_LI(col);++i,zi++)
1273 for(j=0L,zj = S_V_S(row);j<S_V_LI(row);++j,zj++)
1274 if(S_I_I(zi) == S_I_I(zj))
1275 {
1276 erg++;
1277 if(erg>1L) return(erg);
1278 }
1279 return(erg);
1280 }
1281
1282 /***********************************************************************/
1283 /* */
1284 /* Routine: get_H_perm */
1285 /* Gibt die Horizontalpermutation zwischen zwei Tableaux zurueck. */
1286 /* */
1287 /***********************************************************************/
1288 /* RH 011091 */
get_H_perm(t_i,t_j,perm)1289 static INT get_H_perm(t_i,t_j,perm) OP t_i; OP t_j; OP perm;
1290 {
1291 INT i;
1292 INT j;
1293
1294 INT sp;
1295
1296 for(i=0L;i<S_V_LI(t_i);++i)
1297 {
1298 for(j=0L;j<S_V_LI(S_V_I(t_i,i));++j)
1299 {
1300 sp = get_col(t_i,S_V_II(S_V_I(t_j,i),j));
1301 m_i_i(S_V_II(S_V_I(t_j,i),j),
1302 S_P_I(perm,S_V_II(S_V_I(t_j,i),sp)-1L));
1303 }
1304 }
1305 return OK;
1306 }
1307
1308
1309 /***********************************************************************/
1310 /* */
1311 /* Routine: get_H_perm */
1312 /* Gibt die Vertikalpermutation zwischen zwei Tableaux zurueck. */
1313 /* */
1314 /***********************************************************************/
1315 /* RH 011091 */
get_Q_perm(T,t_i,perm)1316 static INT get_Q_perm(T,t_i,perm) OP T,t_i,perm;
1317 {
1318 INT i;
1319 INT j;
1320 INT row;
1321
1322 for(i=0L;i<S_V_LI(t_i);++i)
1323 {
1324 for(j=0L;j<S_V_LI(S_V_I(t_i,i));++j)
1325 {
1326 row = get_row(t_i,S_V_II(S_V_I(T,i),j));
1327 m_i_i(S_V_II(S_V_I(T,i),j),
1328 S_P_I(perm,S_V_II(S_V_I(T,row),j)-1L));
1329 }
1330 }
1331 return OK;
1332 }
1333
1334 /***********************************************************************/
1335 /* */
1336 /* Routine: get_T */
1337 /* Gibt das Tableau T:= perm*t_i zurueck. */
1338 /* */
1339 /***********************************************************************/
1340 /* RH 011091 */
1341
get_T(T,t_i,perm)1342 static INT get_T(T,t_i,perm) OP T; OP t_i; OP perm;
1343 {
1344 INT i;
1345 INT j;
1346
1347 for(i=0L;i<S_V_LI(t_i);++i)
1348 for(j=0L;j<S_V_LI(S_V_I(t_i,i));++j)
1349 m_i_i(S_P_II(perm,S_V_II(S_V_I(t_i,i),j)-1L),
1350 S_V_I(S_V_I(T,i),j));
1351 return OK;
1352 }
1353
1354 /***********************************************************************/
1355 /* */
1356 /* Routine: transponiere (fuer Vektoren von Vektoren) */
1357 /* Transponiert das Tableaux T in Tt */
1358 /* */
1359 /***********************************************************************/
1360 /* RH 011091 */
1361
transponiere(T,Tt)1362 static INT transponiere(T,Tt) OP T; OP Tt;
1363 {
1364 static INT i;
1365 static INT j;
1366
1367 for(i=0L;i<S_V_LI(T);++i)
1368 for(j=0;j<S_V_LI(S_V_I(T,i));++j)
1369 m_i_i(
1370 S_V_II(S_V_I(T,i),j),
1371 S_V_I(S_V_I(Tt,j),i));
1372 return OK;
1373 }
1374
1375 /***********************************************************************/
1376 /* */
1377 /* Routine: get_col */
1378 /* Gibt die Spaltennr. an, in der zahl in tab vorkommt. */
1379 /* */
1380 /***********************************************************************/
1381 /* RH 011091 */
get_col(tab,zahl)1382 static INT get_col(tab,zahl) OP tab; INT zahl;
1383 {
1384 INT i;
1385 INT j;
1386 INT erg = OK;
1387
1388 for(i=0L;i<S_V_LI(tab);++i)
1389 for(j=0L;j<S_V_LI(S_V_I(tab,i));++j)
1390 if(S_V_II(S_V_I(tab,i),j) == zahl)
1391 return(j);
1392 SYMCHECK(1,"should never be here");
1393 ENDR("sab.c:get_col");
1394 }
1395
1396 /***********************************************************************/
1397 /* */
1398 /* Routine: get_row */
1399 /* Gibt die Zeilennr. an, in der zahl in tab vorkommt. */
1400 /* */
1401 /***********************************************************************/
1402 /* RH 011091 */
get_row(tab,zahl)1403 static INT get_row(tab,zahl) OP tab; INT zahl;
1404 {
1405 INT i,j,erg=OK;
1406
1407 for(i=0L;i<S_V_LI(tab);++i)
1408 {
1409 for(j=0L;j<S_V_LI(S_V_I(tab,i));++j)
1410 {
1411 if(S_V_II(S_V_I(tab,i),j) == zahl) return(i);
1412 }
1413 }
1414 SYMCHECK(1,"should never be here");
1415 ENDR("get_row");
1416 }
1417
1418 /***********************************************************************/
1419 /* */
1420 /* Routine: D_row_calc */
1421 /* Berechnet den Beitrag eines Tableaux T_i^(m) zur Zeile i */
1422 /* */
1423 /***********************************************************************/
1424 /* RH 011091 */
1425
1426 #ifdef PERMTRUE
D_row_calc(T,t,sgn,D_row,n,Tt)1427 static INT D_row_calc(T,t,sgn,D_row,n,Tt) OP T,t,D_row,Tt; INT sgn,n;
1428 {
1429 INT i;
1430 INT erg = OK;
1431 OP z,p1;
1432 CTO(VECTOR,"D_row_calc(2)",t);
1433
1434 z = callocobject();
1435 p1 = callocobject();
1436
1437 erg += b_ks_p(VECTOR,callocobject(),p1);
1438 erg += m_il_v(n,S_P_S(p1));
1439
1440 for(i=0L;i<S_V_LI(t);++i)
1441 if(ziffern_existieren(Tt,S_V_I(t,i)) == 0L)
1442 {
1443 erg += get_Q_perm(T,S_V_I(t,i),p1);
1444 erg += signum_permutation(p1,z);
1445 m_i_i(S_V_II(D_row,i)+sgn*S_I_I(z),S_V_I(D_row,i));
1446 }
1447
1448 erg += freeall(p1);
1449 erg += freeall(z);
1450 ENDR("internal:SAB1");
1451 }
1452
1453
1454 /***********************************************************************/
1455 /* */
1456 /* Routine: D_calc */
1457 /* Berechnet die i-te Zeile der darstellenden Matrix D^part(perm) */
1458 /* Berechnung der i-ten Zeile der darstellenden Matrix. Dabei wird */
1459 /* P in dfs- Weise durchlaufen, und dann werden in D_row_calc() */
1460 /* die Anteile der Tableaux T_i^(m) berechnet und aufaddiert. */
1461 /* */
1462 /***********************************************************************/
1463 /* RH 011091 */
1464
D_calc(t,D,P,perm,p_inv,sgn,zeile,i,k,Tt)1465 static INT D_calc(t,D,P,perm,p_inv,sgn,zeile,i,k,Tt)
1466 OP t; OP D; OP P; OP perm;
1467 OP p_inv; INT sgn; INT zeile; INT i;
1468 INT k; OP Tt;
1469 {
1470 INT l;
1471 OP T = callocobject();
1472 OP HP;
1473 OP p1 = callocobject();
1474
1475
1476 b_ks_p(VECTOR,callocobject(),p1);
1477 m_il_v(S_P_LI(perm),S_P_S(p1));
1478
1479 m_il_v(S_V_LI(S_V_I(t,0L)),T);
1480 for(l=0L;l<S_V_LI(T);++l)
1481 m_il_v(S_V_LI(S_V_I(S_V_I(t,0L),l)),S_V_I(T,l));
1482
1483 if(i == S_V_LI(t))
1484 goto DC_ende;
1485 else
1486 {
1487 if(k==i)
1488 {
1489 get_T(T,S_V_I(t,zeile),p_inv);
1490 transponiere(T,Tt);
1491 D_row_calc(T,t,sgn,D,S_P_LI(p_inv),Tt);
1492 }
1493 if (!emptyp(S_V_I(P,i)))
1494 {
1495 HP = S_V_I(P,i);
1496 for(l=0L;l<S_V_LI(HP);++l)
1497 {
1498 if(S_V_II(S_V_I(HP,l),0L) > k)
1499 {
1500 sgn*= -1;
1501 invers(p_inv,p1);
1502 mult(p1,perm,perm);
1503 invers(S_V_I(S_V_I(HP,l),1L),p1);
1504 mult(p1,perm,perm);
1505 mult(p_inv,perm,perm);
1506 get_T(T,S_V_I(t,zeile),perm);
1507 transponiere(T,Tt);
1508 D_row_calc(T,t,sgn,D,S_P_LI(p_inv),Tt);
1509 D_calc(t,D,P,perm,p_inv,sgn,zeile,
1510 S_V_II(S_V_I(HP,l),0L),
1511 S_V_II(S_V_I(HP,l),0L)+1L,Tt);
1512 invers(p_inv,p1);
1513 mult(p1,perm,perm);
1514 mult(S_V_I(S_V_I(HP,l),1L),perm,perm);
1515 mult(p_inv,perm,perm);
1516 sgn*= -1;
1517 }
1518 }
1519 }
1520 }
1521 DC_ende:
1522 freeall(T);
1523 freeall(p1);
1524 return OK;
1525 }
1526
1527 #endif /* PERMTRUE */
1528 /***********************************************************************/
1529 /***********************************************************************/
1530 /* */
1531 /* Programme zur Berechnung der seminormalen und der orthogonalen */
1532 /* Matrixdarstellungen von S_n */
1533 /* */
1534 /* Written by: Ralf Hager September 1991 */
1535 /* */
1536 /***********************************************************************/
1537 /***********************************************************************/
1538
1539 /***********************************************************************/
1540 /* */
1541 /* Routine: sdg */
1542 /* Gibt zu gegebener Partition part und Permutation perm die */
1543 /* zugehoerige seminormale irreduzible Matrixdarstellung in D zurueck. */
1544 /* */
1545 /***********************************************************************/
1546 /* RH 011091 */
1547
1548 #ifdef DGTRUE
sdg(part,perm,D)1549 INT sdg(part,perm,D) OP part,perm,D;
1550 {
1551 INT i;
1552 INT j;
1553 INT erg = OK;
1554 INT lls_vgl();
1555 OP n,lehmer,HD,dim,t,inh;
1556
1557 CTO(PARTITION,"sdg(1)",part);
1558 CTO(PERMUTATION,"sdg(2)",perm);
1559
1560 n = callocobject();
1561 lehmer = callocobject();
1562 HD = callocobject();
1563 dim = callocobject();
1564 t = callocobject();
1565 inh = callocobject();
1566
1567
1568 erg += dimension(part,dim);
1569 erg += m_lh_nm(dim,dim,D);
1570 erg += m_lh_nm(dim,dim,HD);
1571 erg += m_i_i(S_P_LI(perm),n);
1572 erg += m_l_v(n,inh);
1573
1574 for(i=0;i<S_I_I(n);++i)
1575 erg += m_i_i(1L,S_V_I(inh,i));
1576
1577 erg += rz_perm(perm,lehmer);
1578 erg += alpha_beta_tabs(part,inh,t);
1579 erg += sortieren(0L,S_V_LI(t)-1L,t,lls_vgl);
1580
1581 if(!einsp(perm))
1582 {
1583 i = S_V_II(lehmer,S_V_LI(lehmer)-1L);
1584 erg += _xdg(i,t,D,0L);
1585
1586 for(j=S_V_LI(lehmer)-2L;j>=0L;--j)
1587 {
1588 erg += _xdg(S_V_II(lehmer,j),t,HD,0L);
1589 MULT_APPLY(HD,D);
1590 }
1591 }
1592 else
1593 {
1594 for(i=0L;i<S_I_I(dim);++i)
1595 M_I_I(1L,S_M_IJ(D,i,i));
1596 }
1597 if(S_PA_LI(part) == 1L)
1598 {
1599 erg += m_i_i(1L,S_M_IJ(D,0L,0L));
1600 }
1601 if(S_PA_LI(part) == S_P_LI(perm))
1602 {
1603 erg += signum_permutation(perm,dim);
1604 erg += m_i_i(S_I_I(dim),S_M_IJ(D,0L,0L));
1605 }
1606
1607 FREEALL6(n,inh,lehmer,HD,dim,t);
1608
1609 ENDR("sdg");
1610 }
1611
1612
1613 /***********************************************************************/
1614 /* */
1615 /* Routine: _xdg */
1616 /* Gibt zu gegebener Partition und Elementartransposition nr */
1617 /* zugehoerige irreduzible Matrixdarstellung in P zurueck. */
1618 /* In t stehen die Standardtableaux zum gegebenen Umriss */
1619 /* */
1620 /***********************************************************************/
1621 /* RH 011091 */
1622
_xdg(nr,t,P,var)1623 static INT _xdg(nr,t,P,var)
1624 INT nr; OP t; OP P; INT var;
1625 {
1626 INT i,j,k;
1627 INT main_diag = 0L;
1628 INT get_main_diag();
1629 INT erg = OK;
1630
1631 for(i=0L;i<S_V_LI(t);++i)
1632 for(j=0L;j<S_V_LI(t);++j)
1633 m_i_i(0L,S_M_IJ(P,i,j));
1634
1635 for(i=0L;i<S_V_LI(t);++i)
1636 {
1637 /* Hauptdiagonale wird besetzt */
1638
1639 main_diag = get_main_diag(S_V_I(t,i),nr);
1640 m_i_i(main_diag,S_M_IJ(P,i,i));
1641 }
1642
1643 /* quadratische Bloecke wernden besetzt */
1644
1645 for(i=0L;i<S_V_LI(t);++i)
1646 {
1647 for(k=i+1L;k<S_V_LI(t);++k)
1648 {
1649 get_submatrix(
1650 nr,i,k,S_V_I(t,i),S_V_I(t,k),P,var);
1651 }
1652 }
1653 CTO(MATRIX,"_xdg(e3)",P);
1654 ENDR("sab.c:_xdg");
1655 }
1656 #endif /* DGTRUE */
1657 /***********************************************************************/
1658 /* */
1659 /* Routine: get_main_diag */
1660 /* Erstellt einen Eintrag in der Hauptdiagonalen der darst. Mat. */
1661 /* zum Tableau tab und der Elementartranspositon nr */
1662 /* */
1663 /***********************************************************************/
1664 /* RH 011091 */
1665
get_main_diag(tab,nr)1666 static INT get_main_diag(tab,nr) OP tab; INT nr;
1667 {
1668 INT row_nr;
1669 INT row_nr_plus_1;
1670
1671
1672 row_nr = get_row(tab,nr);
1673 row_nr_plus_1 = get_row(tab,nr+1L);
1674
1675 if(row_nr == row_nr_plus_1) return(1L);
1676
1677 row_nr = get_col(tab,nr);
1678 row_nr_plus_1 = get_col(tab,nr+1L);
1679
1680 if(row_nr == row_nr_plus_1) return(-1L);
1681
1682 return(0L);
1683 }
1684
1685 /***********************************************************************/
1686 /* */
1687 /* Routine: get_submatrix */
1688 /* Erstellt, je nachdem ob var = 0 oder var = 1 gesetzt ist, */
1689 /* den Diagonalblock zu den Tableaux t_i und t_k fuer die */
1690 /* darstellende Matrix zur Elementartransposition (nr nr+1) */
1691 /* */
1692 /***********************************************************************/
1693 /* RH 011091 */
1694
1695 #ifdef MATRIXTRUE
get_submatrix(nr,i,k,t_i,t_k,P,var)1696 static INT get_submatrix(nr,i,k,t_i,t_k,P,var)
1697 INT nr; INT i; INT k;
1698 OP t_i; OP t_k; OP P; INT var;
1699 {
1700 INT erg = OK;
1701
1702 OP dist;
1703 OP e_ii;
1704 OP e_ik;
1705 OP e_ki;
1706 OP e_kk;
1707
1708 INT ti_eq_transpo_tk();
1709 INT get_axial_distance();
1710
1711 SYMCHECK( (var != 0)&&(var != 1), "sab:get_submatrix:wrong value of var");
1712
1713 if(ti_eq_transpo_tk(nr,t_i,t_k) == 1L)
1714 {
1715 dist = callocobject();
1716 e_ii = callocobject();
1717 e_ik = callocobject();
1718 e_ki = callocobject();
1719 e_kk = callocobject();
1720 get_axial_distance(t_i,nr,dist);
1721 if(var == 0L)
1722 {
1723 m_i_i(-1L,e_ii);
1724 div(e_ii,dist,e_ii);
1725 copy(e_ii,S_M_IJ(P,i,i));
1726
1727 m_i_i(1L,e_kk);
1728 div(e_ii,dist,e_ik);
1729 add_apply(e_kk,e_ik); /* AK 010692 statt add */
1730 copy(e_ik,S_M_IJ(P,i,k));
1731
1732 m_i_i(1L,e_ki);
1733 copy(e_ki,S_M_IJ(P,k,i));
1734
1735 m_i_i(1L,e_kk);
1736 div(e_kk,dist,e_kk);
1737 copy(e_kk,S_M_IJ(P,k,k));
1738 }
1739 if(var == 1L)
1740 {
1741 m_i_i(-1L,e_ii);
1742 div(e_ii,dist,e_ii);
1743 copy(e_ii,S_M_IJ(P,i,i));
1744
1745 m_i_i(1L,e_kk);
1746 div(e_ii,dist,e_ik);
1747 add_apply(e_kk,e_ik); /* AK 010692 statt add */
1748 squareroot(e_ik,e_ik);
1749 copy(e_ik,S_M_IJ(P,i,k));
1750
1751 copy(e_ik,S_M_IJ(P,k,i));
1752
1753 m_i_i(1L,e_kk);
1754 div(e_kk,dist,e_kk);
1755 copy(e_kk,S_M_IJ(P,k,k));
1756 }
1757 freeall(dist);
1758 freeall(e_ii);
1759 freeall(e_ik);
1760 freeall(e_ki);
1761 freeall(e_kk);
1762 }
1763
1764
1765 ENDR("sab:get_submatrix");
1766 }
1767 #endif /* MATRIXTRUE */
1768
1769 /***********************************************************************/
1770 /* */
1771 /* Routine: ti_eq_transpo_tk */
1772 /* Prueft ob gilt: t_i = transpo*t_k */
1773 /* */
1774 /***********************************************************************/
1775 /* RH 011091 */
1776
ti_eq_transpo_tk(nr,t_i,t_k)1777 static INT ti_eq_transpo_tk(nr,t_i,t_k)
1778 INT nr; OP t_i; OP t_k;
1779 {
1780 INT i_r;
1781 INT i_s;
1782 INT j_r;
1783 INT j_s;
1784 INT erg = 0L;
1785 INT ti_eq_tk();
1786
1787 i_r = get_row(t_k,nr);
1788 j_r = get_col(t_k,nr);
1789 i_s = get_row(t_k,nr+1L);
1790 j_s = get_col(t_k,nr+1L);
1791
1792 swap(S_V_I(S_V_I(t_k,i_r),j_r),
1793 S_V_I(S_V_I(t_k,i_s),j_s));
1794 if(ti_eq_tk(t_i,t_k) == 1L) erg = 1L;
1795 swap(S_V_I(S_V_I(t_k,i_r),j_r),
1796 S_V_I(S_V_I(t_k,i_s),j_s));
1797 return(erg);
1798 }
1799
1800
1801 /***********************************************************************/
1802 /* */
1803 /* Routine: get_axial_distanz */
1804 /* Berechnet die Axialdistanz von nr und nr+1 in t_i und schreibt */
1805 /* sie in dist. */
1806 /* */
1807 /***********************************************************************/
1808 /* RH 011091 */
1809
get_axial_distance(t_i,nr,dist)1810 static INT get_axial_distance(t_i,nr,dist)
1811 OP t_i; INT nr; OP dist;
1812 {
1813 INT get_row();
1814 INT get_col();
1815
1816 INT i_r;
1817 INT j_r;
1818 INT i_s;
1819 INT j_s;
1820
1821 i_r = get_row(t_i,nr);
1822 i_s = get_row(t_i,nr+1L);
1823
1824 j_r = get_col(t_i,nr);
1825 j_s = get_col(t_i,nr+1L);
1826
1827 m_i_i(((j_r-j_s)+(i_s-i_r)),dist);
1828 return OK;
1829 }
1830
1831 /***********************************************************************/
1832 /* */
1833 /* Routine: ti_eq_tk */
1834 /* Prueft ob gilt: t_i == t_k */
1835 /* */
1836 /***********************************************************************/
1837 /* RH 011091 */
1838
ti_eq_tk(t_i,t_k)1839 static INT ti_eq_tk(t_i,t_k) OP t_i; OP t_k;
1840 {
1841 INT i;
1842 INT j;
1843
1844 for(i=0L;i<S_V_LI(t_i);++i)
1845 for(j=0L;j<S_V_LI(S_V_I(t_i,i));++j)
1846 if(S_V_II(S_V_I(t_i,i),j) !=
1847 S_V_II(S_V_I(t_k,i),j))
1848 return(0L);
1849 return(1L);
1850 }
1851
1852 /***********************************************************************/
1853 /* */
1854 /* Routine: lls_vgl */
1855 /* Vergleich nach last letter sequenz zweier Operanden */
1856 /* */
1857 /***********************************************************************/
1858 /* RH 011091 */
1859
lls_vgl(a,b)1860 static INT lls_vgl(a,b) OP a; OP b;
1861 {
1862 INT i;
1863 INT j;
1864 INT k;
1865 INT merk_a=0;
1866 INT merk_b=0;
1867
1868 INT n = 0L;
1869
1870 for(i=0L;i<S_V_LI(a);++i)n+= S_V_LI(S_V_I(a,i));
1871
1872 for(k=n;k>=1L;--k)
1873 {
1874 for(i=0L;i<S_V_LI(a);++i)
1875 {
1876 for(j=0L;j<S_V_LI(S_V_I(a,i));++j)
1877 {
1878 if(S_V_II(S_V_I(a,i),j) == k)
1879 {
1880 merk_a = i;
1881 break;
1882 }
1883 }
1884 }
1885 for(i=0;i<S_V_LI(a);++i)
1886 {
1887 for(j=0L;j<S_V_LI(S_V_I(a,i));++j)
1888 {
1889 if(S_V_II(S_V_I(b,i),j) == k)
1890 {
1891 merk_b = i;
1892 break;
1893 }
1894 }
1895 }
1896 if(merk_a < merk_b) return(-1L);
1897 if(merk_a > merk_b) return(1L);
1898 }
1899 return(0L);
1900 }
1901
1902 /***********************************************************************/
1903 /* */
1904 /* Routine: odg */
1905 /* Gibt zu gegebener Partition part und Permutation perm die */
1906 /* zugehoerige orthogonale irreduzible Matrixdarstellung in D zurueck. */
1907 /* */
1908 /***********************************************************************/
1909 /* RH 011091 */
1910
cyclo_odg(a,b,c)1911 INT cyclo_odg(a,b,c) OP a,b,c;
1912 /* AK 110202 */
1913 {
1914 INT erg = OK;
1915 OP d;
1916 INT i,j;
1917 d = CALLOCOBJECT();
1918 erg += odg(a,b,d);
1919 m_lh_m(S_M_L(d),S_M_H(d),c);
1920 for(i=0;i<S_M_HI(c);i++)
1921 for(j=0;j<S_M_LI(c);j++)
1922 {
1923 if (S_O_K(S_M_IJ(d,i,j)) == SQ_RADICAL)
1924 convert_sqrad_cyclo(S_M_IJ(d,i,j), S_M_IJ(c,i,j));
1925 else
1926 copy(S_M_IJ(d,i,j), S_M_IJ(c,i,j));
1927 }
1928 FREEALL(d);
1929 ENDR("cyclo_odg");
1930 }
1931
cyclo_an_tafel(a,c)1932 INT cyclo_an_tafel(a,c) OP a,c;
1933 /* AK 130202 */
1934 {
1935 INT erg = OK;
1936 OP d;
1937 INT i,j;
1938 d = CALLOCOBJECT();
1939 erg += an_tafel(a,d);
1940 m_lh_m(S_M_L(d),S_M_H(d),c);
1941 for(i=0;i<S_M_HI(c);i++)
1942 for(j=0;j<S_M_LI(c);j++)
1943 {
1944 if (S_O_K(S_M_IJ(d,i,j)) == SQ_RADICAL)
1945 convert_sqrad_cyclo(S_M_IJ(d,i,j), S_M_IJ(c,i,j));
1946 else
1947 copy(S_M_IJ(d,i,j), S_M_IJ(c,i,j));
1948 }
1949 FREEALL(d);
1950 ENDR("cyclo_an_tafel");
1951 }
1952
odg(part,perm,D)1953 INT odg(part,perm,D) OP part,perm,D;
1954 /* RH 1991 */
1955 /* AK 220704 V3.0 */
1956 {
1957 INT erg = OK;
1958 CTO(PARTITION,"odg(1)",part);
1959 CTO(PERMUTATION,"odg(2)",perm);
1960 CE3(part,perm,D,odg);
1961 {
1962 #ifdef DGTRUE
1963 INT i;
1964 INT j;
1965 INT lls_vgl();
1966 OP rz,HD,dim,t,inh;
1967
1968
1969 HD = callocobject();
1970 dim = callocobject();
1971 t = callocobject();
1972 inh = callocobject();
1973 rz = callocobject();
1974
1975 erg += weight(part,dim);
1976 if (neq(dim, S_P_L(perm)) )
1977 {
1978 error("odg: part and perm have different weights");
1979 erg += ERROR;
1980 goto odge;
1981 }
1982
1983
1984 erg += dimension(part,dim);
1985 erg += m_lh_nm(dim,dim,D);
1986 erg += m_lh_nm(dim,dim,HD);
1987 erg += m_il_v(S_P_LI(perm),inh);
1988
1989 for(i=0;i<S_V_LI(inh);++i)
1990 M_I_I(1L,S_V_I(inh,i));
1991
1992 erg += rz_perm(perm,rz);
1993 CTTO(INTEGERVECTOR,VECTOR,"odg(i)",rz);
1994 erg += alpha_beta_tabs(part,inh,t);
1995 erg += sortieren(0L,S_V_LI(t)-1L,t,lls_vgl);
1996
1997 if(!einsp(perm))
1998 {
1999 i = S_V_II(rz,S_V_LI(rz)-1L);
2000 erg += _xdg(i,t,D,1L);
2001
2002 for(j=S_V_LI(rz)-2L;j>=0L;--j)
2003 {
2004 erg += _xdg(S_V_II(rz,j),t,HD,1L);
2005 MULT_APPLY(HD,D);
2006 }
2007 }
2008 else
2009 {
2010 for(i=0L;i<S_I_I(dim);++i)
2011 erg += m_i_i(1L,S_M_IJ(D,i,i));
2012 }
2013 if(S_PA_LI(part) == 1L)
2014 {
2015 erg += m_i_i(1,S_M_IJ(D,0,0));
2016 }
2017 if(S_PA_LI(part) == S_P_LI(perm))
2018 {
2019 erg += signum_permutation(perm,dim);
2020 CLEVER_COPY(dim,S_M_IJ(D,0,0));
2021 }
2022
2023 odge:
2024 FREEALL5(inh,rz,HD,dim,t);
2025 #endif
2026 }
2027 ENDR("odg");
2028 }
2029 #ifdef DGTRUE
2030 /***********************************************************************/
2031 /***********************************************************************/
2032 /* */
2033 /* Programm zur Berechnung einer polynomialen, irreduziblen */
2034 /* Matrixdarstellung von GL_m(C). */
2035 /* Dabei wird zu einer Partition part ein Satz von orthonormalen */
2036 /* Basisvektoren bestimmt, mit denen dann ein Basiswechsel */
2037 /* durchgefuehrt wird. Es werden in diesem Spezialfall nur die */
2038 /* Projektoren P_11 berechnet, statt einer Matrixinversion kann */
2039 /* man die transponierte Matrix betrachten. */
2040 /* Zum orthonormalisieren wird das Gram-Schmidt-Verfahren verwendet. */
2041 /* Das Ergebnis steht in Matrix D. */
2042 /* */
2043 /* Written by: Ralf Hager September 1991 */
2044 /* */
2045 /***********************************************************************/
2046 /***********************************************************************/
2047
2048
2049 /***********************************************************************/
2050 /* */
2051 /* Routine: glpdg */
2052 /* Diese Routine ruft zunaechst glm_sab auf, wo die Basis B berechnet */
2053 /* wird. Im Unterprogramm glm_B_W wird dann der Basiswechsel */
2054 /* durchgefuehrt. */
2055 /* */
2056 /***********************************************************************/
2057 /* RH 011091 */
2058
glpdg(m,part,D)2059 INT glpdg(m,part,D) OP m; OP part; OP D;
2060 {
2061 INT erg = OK; /* AK 010692 */
2062 OP B = callocobject();
2063 OP n = callocobject();
2064
2065 erg += weight(part,n); /* AK 010692 */
2066
2067 erg += glm_sab(m,n,B,part); /* Symmetrieangepasste Basis errechnen */
2068 erg += glm_B_W(m,n,B,D); /* Basiswechsel durchfuehren */
2069
2070 erg += freeall(n);
2071 erg += freeall(B);
2072 return erg;
2073 }
2074 #endif /* DGTRUE */
2075
2076 /***********************************************************************/
2077 /* */
2078 /* Routine: glm_sab */
2079 /* Die Routine verlaeuft analog zur Routine sab, nur dass hier */
2080 /* ausgenutzt werden kann, dass die S_n nicht vorher erzeugt werden */
2081 /* muss. Ausserdem wird nur P_11 zu part bestimmt, da man sich nur */
2082 /* fuer einen Block der homogenen Komponente interessiert. */
2083 /* */
2084 /***********************************************************************/
2085 /* RH 011091 */
2086
2087 #ifdef SABTRUE
glm_sab(m,n,B,part)2088 INT glm_sab(m,n,B,part) OP m; OP n; OP B; OP part;
2089 {
2090 INT TP();
2091
2092 OP f_D = callocobject();
2093 OP irr_dim = callocobject();
2094 OP faktor = callocobject();
2095 OP dim;
2096 OP perm = callocobject();
2097 OP tupel = callocobject();
2098 OP D = callocobject();
2099 OP X = callocobject();
2100 OP PROJ = callocobject();
2101
2102 NEW_INTEGER(dim,0);
2103 hoch(m,n,f_D);
2104 dimension_symmetrization(m,part,dim);
2105 m_lh_nm(dim,f_D,B);
2106 m_l_nv(n,tupel);
2107 m_il_p(S_I_I(f_D),X);
2108
2109 if(S_PA_LI(part) <= S_I_I(m))
2110 {
2111 dimension_partition(part,irr_dim);
2112 m_lh_nm(f_D,f_D,PROJ);
2113 first_permutation(n,perm);
2114 do
2115 {
2116 invers(perm,perm);
2117 odg(part,perm,D);
2118 invers(perm,perm);
2119 m_l_nv(n,tupel);
2120 TP(tupel,0L,n,m,perm,X);
2121 copy(S_M_IJ(D,0L,0L),faktor);
2122 if(!nullp(faktor))
2123 add_to_PROJ(PROJ,faktor,X);
2124 }while(next(perm,perm));
2125 m_i_i(0L,dim);
2126 glm_get_BV(PROJ,B,dim);
2127 reell_gram_schmidt(B);
2128 }
2129 freeall(PROJ);
2130 freeall(dim);
2131 freeall(f_D);
2132 freeall(irr_dim);
2133 freeall(X);
2134 freeall(faktor);
2135 freeall(tupel);
2136 freeall(perm);
2137 freeall(D);
2138 return OK;
2139 }
2140 #endif /* SABTRUE */
2141
2142 /***********************************************************************/
2143 /* */
2144 /* Routine: reell_gram_schmidt */
2145 /* Dieses Unterprogramm arbeitet analog zur Routine gram_schmidt. */
2146 /* Es wird hier jedoch keine Konjugation durchgefuehrt, da die */
2147 /* betrachteten Vektoren nur relle Eintraege besitzen. */
2148 /* */
2149 /***********************************************************************/
2150 /* RH 011091 */
2151
2152
2153 #ifdef MATRIXTRUE
reell_gram_schmidt(A)2154 INT reell_gram_schmidt(A) OP A;
2155 {
2156 INT i;
2157 INT j;
2158 INT k;
2159
2160 OP R = callocobject();
2161 OP h = callocobject();
2162 OP n = callocobject();
2163 OP p = callocobject();
2164 OP s = callocobject();
2165
2166 m_lh_nm(S_M_H(A),S_M_L(A),R);
2167 m_i_i(S_M_HI(A),n);
2168 m_i_i(S_M_LI(A),p);
2169
2170 for(k=0L;k<S_I_I(p);++k)
2171 {
2172 for(j=0L;j<k;++j)
2173 {
2174 m_i_i(0L,S_M_IJ(R,j,k));
2175 for(i=0L;i<S_I_I(n);++i)
2176 {
2177 mult(S_M_IJ(A,i,j),S_M_IJ(A,i,k),h);
2178 add_apply(h,S_M_IJ(R,j,k));
2179 }
2180 for(i=0L;i<S_I_I(n);++i)
2181 {
2182 mult(S_M_IJ(R,j,k),S_M_IJ(A,i,j),h);
2183 addinvers_apply(h);
2184 add_apply(h,S_M_IJ(A,i,k));
2185 }
2186 }
2187 m_i_i(0L,s);
2188 for(i=0L;i<S_I_I(n);++i)
2189 {
2190 mult(S_M_IJ(A,i,k),S_M_IJ(A,i,k),h);
2191 add_apply(h,s);
2192 }
2193 squareroot(s,S_M_IJ(R,k,k));
2194 invers(S_M_IJ(R,k,k),h);
2195 for(i=0L;i<S_I_I(n);++i)
2196 mult_apply(h,S_M_IJ(A,i,k));
2197 }
2198
2199 freeall(R);
2200 freeall(h);
2201 freeall(n);
2202 freeall(p);
2203 freeall(s);
2204 return OK;
2205 }
2206 #endif /* MATRIXTRUE */
2207
2208 /***********************************************************************/
2209 /* */
2210 /* Routine: get_poly_self */
2211 /* Diese Routine besetzt, ausgehend von einem Index den self-teil */
2212 /* eines polynomialen Eintrages der Matrix "n-faches Tensorprodukt */
2213 /* der generischen GLm(C) matrix. */
2214 /* Der Zeilenindex row und der Spaltenindex col werden dabei m-naer */
2215 /* ausgedrueckt. Dies ergibt die Indices des Monoms. Diese Indices */
2216 /* werden in die Matrix polyself uebertragen. */
2217 /* */
2218 /***********************************************************************/
2219 /* RH 011091 */
2220
2221
2222 #ifdef MATRIXTRUE
get_poly_self(row,col,m,n,polyself)2223 static INT get_poly_self(row,col,m,n,polyself)
2224 INT row; INT col; OP m; OP n; OP polyself;
2225 /* AK 010692 */
2226 {
2227 INT i;
2228 OP dim,zeile,spalte,oprow,opcol,x,y;
2229 INT erg = OK; /* AK 030197 */
2230 CTO(INTEGER,"get_poly_self(4)",n);
2231
2232
2233 dim = callocobject();
2234 zeile = callocobject();
2235 spalte = callocobject();
2236 oprow = callocobject();
2237 opcol = callocobject();
2238 x = callocobject();
2239 y = callocobject();
2240
2241 erg += m_l_nv(n,zeile);
2242 erg += m_l_nv(n,spalte);
2243 erg += hoch(m,n,dim);
2244 erg += m_i_i(row,oprow);
2245 erg += m_i_i(col,opcol);
2246 for(i=0L;i<S_I_I(n);++i)
2247 {
2248 erg += m_i_i(S_I_I(dim)/S_I_I(m),dim);
2249 erg += m_i_i(S_I_I(opcol)/S_I_I(dim),S_V_I(spalte,i));
2250 erg += m_i_i(S_I_I(oprow)/S_I_I(dim),S_V_I(zeile,i));
2251 erg += m_i_i(S_I_I(dim)*S_V_II(zeile,i),x);
2252 erg += m_i_i(S_I_I(dim)*S_V_II(spalte,i),y);
2253 erg += m_i_i(S_I_I(oprow)-S_I_I(x),oprow);
2254 erg += m_i_i(S_I_I(opcol)-S_I_I(y),opcol);
2255 }
2256
2257 erg += freeall(x);
2258 x = cons_null;
2259 for(i=0;i<S_I_I(n);++i)
2260 if(S_V_II(zeile,i) >= S_I_I(x)) /* copy(S_V_I(zeile,i),x); */
2261 x = S_V_I(zeile,i);
2262 freeall(y);
2263 y = cons_null;
2264 for(i=0;i<S_I_I(n);++i)
2265 if(S_V_II(spalte,i) >= S_I_I(y)) /* copy(S_V_I(spalte,i),y); */
2266 y = S_V_I(spalte,i);
2267 erg += m_ilih_nm(S_I_I(y)+1L,S_I_I(x)+1L,polyself);
2268 for(i=0L;i<S_I_I(n);++i)
2269 inc_integer(S_M_IJ(polyself,S_V_II(zeile,i),S_V_II(spalte,i)));
2270
2271 erg += freeall(oprow);
2272 erg += freeall(opcol);
2273 erg += freeall(dim);
2274 erg += freeall(zeile);
2275 erg += freeall(spalte);
2276 ENDR("get_poly_self");
2277 }
2278 #endif /* MATRIXTRUE */
2279 /***********************************************************************/
2280 /* */
2281 /* Routine: glm_B_W */
2282 /* Der Basiswechsel zur erstellung von <<part>> wird durchgefuehrt. */
2283 /* Er hat den Aufwand O(m^2n f^2), ist also extrem aufwendig. */
2284 /* Dabei ist f die Dimension derSymmetrisierung <<part>>. */
2285 /* */
2286 /***********************************************************************/
2287 /* RH 011091 */
2288
2289 #ifdef DGTRUE
glm_B_W(m,n,B,D)2290 INT glm_B_W(m,n,B,D)
2291 OP m; OP n; OP B; OP D;
2292 {
2293 INT i,j,k,l;
2294
2295 OP polyself = callocobject();
2296 OP erg = callocobject();
2297 OP x = callocobject();
2298 OP mneu = callocobject();
2299
2300 m_lh_m(S_M_L(B),S_M_L(B),D);
2301
2302 for(i=0L;i<S_M_LI(B);++i)
2303 {
2304 for(j=0L;j<S_M_LI(B);++j)
2305 {
2306 m_i_i(0L,mneu);
2307 for(k=0L;k<S_M_HI(B);++k)
2308 {
2309 for(l=0L;l<S_M_HI(B);++l)
2310 {
2311 if(not nullp(S_M_IJ(B,k,i)) &&
2312 not nullp(S_M_IJ(B,l,j)))
2313 {
2314
2315 mult(S_M_IJ(B,k,i),S_M_IJ(B,l,j),erg);
2316 get_poly_self(k,l,m,n,polyself);
2317 m_skn_po(polyself,erg,NULL,x);
2318 add_apply(x,mneu);
2319 }
2320 }
2321 }
2322 copy(mneu,S_M_IJ(D,i,j));
2323 }
2324 }
2325
2326 freeall(polyself);
2327 freeall(erg);
2328 freeall(x);
2329 freeall(mneu);
2330 return OK;
2331 }
2332 #endif /* DGTRUE */
2333
2334 /***********************************************************************/
2335 /* */
2336 /* Routine: TP */
2337 /* Rekursive Routine zur Berechnung der Permutationsdarstellung */
2338 /* von S_n auf dem Tensorraum C^(m^n). Es werden dabei alle */
2339 /* n-stelligen Tupel mit Eintraegen aus m durchlaufen. Die */
2340 /* Permutationsmatrix wird als Permutation bperm codiert. */
2341 /* */
2342 /***********************************************************************/
2343 /* RH 011091 */
2344
2345 #ifdef DGTRUE
TP(tupel,st,n,m,sperm,bperm)2346 static INT TP(tupel,st,n,m,sperm,bperm)
2347 OP tupel; INT st; OP n; OP m;
2348 OP sperm; OP bperm;
2349 {
2350 INT i;
2351 INT zeile = 0L;
2352 INT spalte = 0L;
2353 INT TP();
2354 INT get_nr_of_tupel();
2355 static INT count = 0L;
2356 OP htupel;
2357
2358 if(st == S_I_I(n))
2359 {
2360 htupel = callocobject();
2361 count++;
2362 m_l_nv(n,htupel);
2363 bilde_htupel(sperm,tupel,htupel);
2364 spalte = get_nr_of_tupel(tupel,m);
2365 zeile = get_nr_of_tupel(htupel,m);
2366 m_i_i(zeile+1L,S_P_I(bperm,spalte));
2367 freeall(htupel);
2368 }
2369 else
2370 {
2371 for(i=0L;i<S_I_I(m);++i)
2372 {
2373 m_i_i((INT)(S_V_II(tupel,st)+i),S_V_I(tupel,st));
2374 TP(tupel,(INT)(st+1L),n,m,sperm,bperm);
2375 m_i_i((INT)(S_V_II(tupel,st)-i),S_V_I(tupel,st));
2376 }
2377 }
2378 return OK;
2379 }
2380 #endif /* DGTRUE */
2381
2382 /***********************************************************************/
2383 /* */
2384 /* Routine: bilde_htupel */
2385 /* Hilfsroutine zur Tupelbildung. */
2386 /* */
2387 /***********************************************************************/
2388 /* RH 011091 */
2389
bilde_htupel(perm,tupel,htupel)2390 INT bilde_htupel(perm,tupel,htupel)
2391 OP perm; OP tupel; OP htupel;
2392 {
2393 INT i;
2394
2395 invers(perm,perm);
2396 for(i=0L;i<S_P_LI(perm);++i)
2397 {
2398 m_i_i(S_V_II(tupel,S_P_II(perm,i)-1L),S_V_I(htupel,i));
2399 }
2400 invers(perm,perm);
2401 return OK;
2402 }
2403
2404
2405
2406 /***********************************************************************/
2407 /* */
2408 /* Routine: get_nr_of_tupel */
2409 /* Hilfsroutine zur Tupelbildung. */
2410 /* Ein Tupel wird als m-naere Zahl aufgefasst. Diese wird berechnet */
2411 /* und zurueckgegeben. */
2412 /* */
2413 /***********************************************************************/
2414 /* RH 011091 */
2415
get_nr_of_tupel(tupel,m)2416 static INT get_nr_of_tupel(tupel,m) OP tupel; OP m;
2417 {
2418 INT i;
2419 INT expo;
2420 INT nr = 0L;
2421 INT mhochexpo();
2422
2423 for(i=0L;i<S_V_LI(tupel);++i)
2424 {
2425 if(S_V_II(tupel,i) > 0L)
2426 {
2427 expo = S_V_LI(tupel)-i-1L;
2428 expo = mhochexpo(S_I_I(m),expo);
2429 expo*= S_V_II(tupel,i);
2430 nr+= expo;
2431 }
2432 }
2433
2434 return(nr);
2435 }
mhochexpo(x,y)2436 static INT mhochexpo(x,y) INT x; INT y;
2437 {
2438 INT i;
2439 INT erg = 1L;
2440
2441 if(y == 0L)return(1L);
2442 erg = x;
2443 for(i=2L;i<=y;++i)
2444 erg*=x;
2445 return(erg);
2446 }
2447
2448 /***********************************************************************/
2449 /* */
2450 /* Routine: dimension_symmetrization */
2451 /* Die Dimension der Symmetrisierung von part in GL_m(C) wird bestimmt.*/
2452 /* Grundlange ist die in JK79 S. 188 angegebene Formel. */
2453 /* */
2454 /***********************************************************************/
2455 /* RH 011091 */
2456 /* = anzahl der tableaux vom umriss part und max eintrag m */
2457
2458 #ifdef PARTTRUE
dimension_symmetrization(m,part,dim)2459 INT dimension_symmetrization(m,part,dim) OP m,part,dim;
2460 /* RH 011091 */ /* AK 010692 */
2461 {
2462 INT i;
2463 INT j;
2464 INT erg = OK;
2465
2466 OP nfak;
2467 OP f_part;
2468 OP hpart;
2469 OP _i;
2470 OP _j;
2471 OP h;
2472
2473 CTO(INTEGER,"dimension_symmetrization(1)",m);
2474 CTO(PARTITION,"dimension_symmetrization(2)",part);
2475
2476 nfak = callocobject();
2477 f_part = callocobject();
2478 hpart = callocobject();
2479 _i = callocobject();
2480 _j = callocobject();
2481 h = callocobject();
2482
2483 m_i_i(1L,dim);
2484 for(i=0L;i<S_PA_LI(part);++i)
2485 erg += add_apply(S_PA_I(part,i),nfak);
2486 erg += m_l_v(S_PA_L(part),hpart);
2487 for(i=0;i<S_PA_LI(part);++i)
2488 erg += copy(S_PA_I(part,S_PA_LI(part)-i-1L),S_V_I(hpart,i));
2489 erg += fakul(nfak,nfak);
2490 erg += dimension_partition(part,f_part);
2491
2492 for(i=1L;i<=S_V_LI(hpart);++i)
2493 for(j=1L;j<=S_V_II(hpart,i-1L);++j)
2494 {
2495 m_i_i(-i,_i);
2496 m_i_i( j,_j);
2497 erg += add(_i,_j,h);
2498 erg += add_apply(m,h);
2499 erg += mult_apply(h,dim);
2500 }
2501
2502 erg += div(dim,nfak,dim);
2503 erg += mult(f_part,dim,dim);
2504
2505 erg += freeall(nfak);
2506 erg += freeall(f_part);
2507 erg += freeall(hpart);
2508 erg += freeall(h);
2509 erg += freeall(_j);
2510 erg += freeall(_i);
2511 ENDR("dimension_symmetrization");
2512 }
2513 #endif /* PARTTRUE */
2514
2515 /***********************************************************************/
2516 /* */
2517 /* Routine: glm_get_BV */
2518 /* Arbeitet analog zur Routine get_BV */
2519 /* */
2520 /***********************************************************************/
2521 /* RH 011091 */
2522
glm_get_BV(PROJ,B,count)2523 INT glm_get_BV(PROJ,B,count)
2524 OP PROJ; OP B; OP count;
2525 {
2526 INT i;
2527 INT j;
2528 INT k;
2529 INT l;
2530 INT s;
2531 INT f_D;
2532
2533 OP Proj = callocobject();
2534 OP val = callocobject();
2535 OP faktor = callocobject();
2536 OP h = callocobject();
2537
2538 copy(PROJ,Proj);
2539 f_D = S_M_HI(B);
2540 k = 0L;
2541 for(i=0L;i<f_D;++i)
2542 {
2543 for(j=k;j<f_D;++j)
2544 {
2545 if(!nullp(S_M_IJ(Proj,j,i)))
2546 {
2547 for(l=0L;l<f_D;++l)
2548 {
2549 copy(S_M_IJ(PROJ,l,i),S_M_IJ(B,l,S_I_I(count)));
2550 }
2551 inc(count);
2552 if(j != k)
2553 for(l=0L;l<f_D;++l)
2554 {
2555 swap(S_M_IJ(Proj,k,l),S_M_IJ(Proj,j,l));
2556 }
2557 invers(S_M_IJ(Proj,k,i),h);
2558 for(l=k+1L;l<f_D;++l)
2559 if(!nullp(S_M_IJ(Proj,l,i)))
2560 {
2561 mult(S_M_IJ(Proj,l,i),h,faktor);
2562 for(s=i;s<f_D;++s)
2563 {
2564 /*
2565 mult(faktor,S_M_IJ(Proj,k,s),val);
2566 sub(S_M_IJ(Proj,l,s),val,val);
2567 copy(val,S_M_IJ(Proj,l,s));
2568 */
2569 mult(faktor,S_M_IJ(Proj,k,s),val);
2570 addinvers_apply(val);
2571 add_apply(val,S_M_IJ(Proj,l,s));
2572 }
2573 }
2574 }
2575 }
2576 k++;
2577 }
2578
2579 freeall(Proj);
2580 freeall(val);
2581 freeall(faktor);
2582 freeall(h);
2583 return OK;
2584 }
2585
2586 /***********************************************************************/
2587 /***********************************************************************/
2588 /* */
2589 /* Homomorphietest fuer GLm(C) Darstellungen. */
2590 /* Es werden 2 Matrizen A,B mit Zufallszahlen besetzt, und dann in */
2591 /* <<part>> eingesetzt. Dann wird A*B gebildet, eingesetzt und auf */
2592 /* Gleichheit getestet. */
2593 /* */
2594 /* written by: Ralf Hager (August 1991) */
2595 /***********************************************************************/
2596 /***********************************************************************/
2597 /* RH 011091 */
2598
2599 /***********************************************************************/
2600 /* */
2601 /* Routine: glm_homtest */
2602 /* Hauptroutine fuer den Homomorphietest. */
2603 /* */
2604 /***********************************************************************/
2605 /* RH 011091 */
2606
2607 #ifdef DGTRUE
glm_homtest(m,d)2608 INT glm_homtest(m,d) OP m,d;
2609 {
2610 OP a;
2611 OP b;
2612 OP anz;
2613
2614
2615 a = callocobject();
2616 b = callocobject();
2617 anz = callocobject();
2618
2619 m_lh_nm(m,m,a);
2620 m_lh_nm(m,m,b);
2621 bestimme_zufallsmatrizen(m,a,b);
2622 if(_homtest(a,b,d) == (INT)1)
2623 {
2624 printf("Homtest OK\n");
2625 }
2626 else
2627 {
2628 printf("Fehler in Homtest\n");
2629 }
2630
2631 freeall(a);
2632 freeall(b);
2633 freeall(anz);
2634 return OK;
2635 }
2636 #endif /* DGTRUE */
2637
2638 /***********************************************************************/
2639 /* */
2640 /* Routine: bestimme_zufallsmatrizen */
2641 /* Es werden zwei Matrizen A,B mit Zufallszahlen zwischen -5 und 5 */
2642 /* besetzt. */
2643 /* */
2644 /***********************************************************************/
2645 /* RH 011091 */
2646
2647 #ifdef MATRIXTRUE
bestimme_zufallsmatrizen(m,a,b)2648 INT bestimme_zufallsmatrizen(m,a,b) OP m,a,b;
2649 {
2650 INT i;
2651 INT j;
2652 OP zahl = callocobject();
2653 OP ober = callocobject();
2654 OP unter = callocobject();
2655
2656 m_i_i((INT)-5,unter);
2657 m_i_i((INT)5,ober);
2658 for(i=0L;i<S_I_I(m);++i)
2659 {
2660 for(j=0L;j<S_I_I(m);++j)
2661 {
2662 random_integer(S_M_IJ(a,i,j),unter,ober);
2663 random_integer(S_M_IJ(b,i,j),unter,ober);
2664 }
2665 }
2666
2667 freeall(zahl); freeall(unter); freeall(ober);
2668 return OK;
2669 }
2670 #endif /* MATRIXTRUE */
2671 /***********************************************************************/
2672 /* */
2673 /* Routine: _homtest */
2674 /* Die Homomorphieeigenschaft wird nachgerechnet. */
2675 /* */
2676 /***********************************************************************/
2677 /* RH 011091 */
2678
2679 #ifdef DGTRUE
_homtest(a,b,m)2680 INT _homtest(a,b,m) OP a,b,m;
2681 {
2682 OP ab = callocobject();
2683 OP d_a = callocobject();
2684 OP d_b = callocobject();
2685 OP d_ab = callocobject();
2686
2687 m_lh_nm(S_M_L(m),S_M_L(m),d_a);
2688 m_lh_nm(S_M_L(m),S_M_L(m),d_b);
2689 m_lh_nm(S_M_L(m),S_M_L(m),d_ab);
2690
2691 mult(a,b,ab);
2692 bestimme_D(m,a,d_a);
2693 bestimme_D(m,b,d_b);
2694 bestimme_D(m,ab,d_ab);
2695 mult(d_a,d_b,d_a);
2696 if(!EQ(d_a,d_ab))
2697 {
2698 error("Hilfe, keine Darstellung !!!!!!!!!!!!\n");
2699
2700 freeall(ab);
2701 freeall(d_a);
2702 freeall(d_b);
2703 freeall(d_ab);
2704 return((INT)0);
2705 }
2706 else
2707 {
2708 freeall(ab);
2709 freeall(d_a);
2710 freeall(d_b);
2711 freeall(d_ab);
2712 return((INT)1);
2713 }
2714 }
2715 #endif /* DGTRUE */
2716
2717 /***********************************************************************/
2718 /* */
2719 /* Routine: bestimme_D */
2720 /* Die Matrix A wird in M eingesetzt und das Ergebnis auf D_A */
2721 /* geschrieben. */
2722 /* */
2723 /***********************************************************************/
2724 /* RH 011091 */
2725
2726 #ifdef DGTRUE
bestimme_D(m,a,d_a)2727 INT bestimme_D(m,a,d_a) OP m,a,d_a;
2728 {
2729 INT i;
2730 INT j;
2731
2732 for(i=0L;i<S_M_LI(m);++i)
2733 {
2734 for(j=0L;j<S_M_LI(m);++j)
2735 {
2736 werte_Polynom_aus(a,S_M_IJ(m,i,j),S_M_IJ(d_a,i,j));
2737 }
2738 }
2739 return OK;
2740 }
2741 #endif /* DGTRUE */
2742
2743 /***********************************************************************/
2744 /* */
2745 /* Routine: werte_Polynom_aus */
2746 /* Das multivariate Polynom poly wird mit A ausgewertet. */
2747 /* Das Ergebnis steht in erg. */
2748 /* */
2749 /***********************************************************************/
2750 /* RH 011091 */
2751
2752 #ifdef MATRIXTRUE
werte_Polynom_aus(a,poly,erg)2753 INT werte_Polynom_aus(a,poly,erg)
2754 OP a; OP poly; OP erg;
2755 {
2756 INT k;
2757 INT l;
2758
2759 OP z = poly;
2760 OP x = callocobject();
2761 OP y = callocobject();
2762
2763 m_i_i(0L,x);
2764 m_i_i(0L,y);
2765 m_i_i(0L,erg);
2766
2767 while(z != NULL)
2768 {
2769 if(!nullp(S_PO_K(z)) && !emptyp(S_PO_K(z)))
2770 {
2771 copy(S_PO_K(z),x);
2772 for(k=0L;k<S_M_HI(S_PO_S(z));++k)
2773 {
2774 for(l=0L;l<S_M_LI(S_PO_S(z));++l)
2775 if(S_M_IJI(S_PO_S(z),k,l) > 0L)
2776 {
2777 hoch(S_M_IJ(a,k,l),S_M_IJ(S_PO_S(z),k,l),y);
2778 mult_apply(y,x); /* AK 120692 statt mult */
2779 }
2780
2781 }
2782 }
2783 z = S_PO_N(z);
2784 add_apply(x,erg); /* AK 010692 statt add */
2785 }
2786 freeall(x); /* AK 010692 */
2787 freeall(y); /* AK 010692 */
2788 return OK;
2789 }
2790 #endif /* MATRIXTRUE */
2791
2792 /***********************************************************************/
2793 /***********************************************************************/
2794 /* */
2795 /* Programm zur Bestimmung von GLm(C) Matrixdarstellungen ohne */
2796 /* Orthonormalisierungsschritt. */
2797 /* Es kann durch var zwischen der orthogonalen (var = 0L) und der */
2798 /* Boernerschen (var = 1L) Matrixdarstellung gewaehlt werden. */
2799 /* Es wird eine m^n-dimensionale Matrix zerlegt. Sie muss im Fall */
2800 /* von (S_n,GLm(C)) jedoch nicht explizit bekannt sein, sondern */
2801 /* wird elementweise berechnet. */
2802 /* Ansonsten ist das Verfahren identisch zu dem in sab beschriebenen. */
2803 /* */
2804 /* written by: Ralf Hager (August 1991) */
2805 /***********************************************************************/
2806 /***********************************************************************/
2807
2808 /***********************************************************************/
2809 /* */
2810 /* Routine: input_glmn */
2811 /* Ein Erzeugendensystem fuer Sn wird bestimmt zusammen mit den */
2812 /* darstellenden Matrizen der irreduziblen Matrixdarstellungen */
2813 /* part mit part_1' <= m. */
2814 /* */
2815 /***********************************************************************/
2816 /* RH 011091 */
2817
2818 #ifdef DGTRUE
input_glmn(m,n,S,SMat,var)2819 INT input_glmn(m,n,S,SMat,var)
2820 OP m; OP n; OP S; OP SMat; INT var;
2821 {
2822 INT i;
2823
2824 OP part = callocobject();
2825 OP anz_irr = callocobject();
2826 OP perm1 = callocobject();
2827 OP perm2 = callocobject();
2828 OP m_hoch_n= callocobject();
2829 OP tupel = callocobject();
2830
2831 m_i_i(0L,anz_irr);
2832
2833 hoch(m,n,m_hoch_n);
2834 first_partition(n,part);
2835 do
2836 {
2837 if(S_PA_LI(part)<= S_I_I(m)) inc(anz_irr);
2838 } while(next(part,part));
2839
2840 if(S_I_I(n) > 2)
2841 {
2842 m_il_v(2L,S);
2843 m_il_p(S_I_I(m_hoch_n),S_V_I(S,0L));
2844 m_il_p(S_I_I(m_hoch_n),S_V_I(S,1L));
2845 m_il_p(S_I_I(n),perm1);
2846 m_il_p(S_I_I(n),perm2);
2847
2848 for(i=2L;i<=S_I_I(n);++i)
2849 {
2850 m_i_i(i,S_P_I(perm2,i-2L));
2851 }
2852 m_i_i(1L,S_P_I(perm2,S_I_I(n)-1L));
2853
2854 for(i=1L;i<=S_I_I(n);++i)
2855 {
2856 m_i_i(i,S_P_I(perm1,i-1L));
2857 }
2858
2859 m_i_i(2L,S_P_I(perm1,0L));
2860 m_i_i(1L,S_P_I(perm1,1L));
2861 m_l_v(anz_irr,SMat);
2862 for(i=0L;i<S_I_I(anz_irr);++i) m_il_v(2L,S_V_I(SMat,i));
2863
2864 m_l_nv(n,tupel);
2865 TP(tupel,0L,n,m,perm1,S_V_I(S,0L));
2866 m_l_nv(n,tupel);
2867 TP(tupel,0L,n,m,perm2,S_V_I(S,1L));
2868
2869 i = 0L;
2870 first_partition(n,part);
2871 do
2872 {
2873 if(S_PA_LI(part)<=S_I_I(m))
2874 {
2875 if(var == 0L)
2876 {
2877 odg(part,perm1,S_V_I(S_V_I(SMat,i),0L));
2878 odg(part,perm2,S_V_I(S_V_I(SMat,i),1L));
2879 }
2880 if(var == 1L)
2881 {
2882 bdg(part,perm1,S_V_I(S_V_I(SMat,i),0L));
2883 bdg(part,perm2,S_V_I(S_V_I(SMat,i),1L));
2884 }
2885 i++;
2886 }
2887 } while(next(part,part));
2888 }
2889 else
2890 {
2891 m_il_v(1L,S);
2892 m_il_p(S_I_I(m_hoch_n),S_V_I(S,0L));
2893 m_il_p(S_I_I(n),perm1);
2894 for(i=2L;i<=S_I_I(n);++i)
2895 {
2896 m_i_i(i,S_P_I(perm1,i-2L));
2897 }
2898 m_i_i(1L,S_P_I(perm1,S_I_I(n)-1L));
2899 m_l_v(anz_irr,SMat);
2900 for(i=0L;i<S_I_I(anz_irr);++i) m_il_v(1L,S_V_I(SMat,i));
2901 m_l_nv(n,tupel);
2902 TP(tupel,0L,n,m,perm1,S_V_I(S,0L));
2903 i = 0L;
2904 first_partition(n,part);
2905 do
2906 {
2907 if(S_PA_LI(part)<=S_I_I(m))
2908 {
2909 if(var == 0L)
2910 {
2911 odg(part,perm1,S_V_I(S_V_I(SMat,i),0L));
2912 }
2913 if(var == 1L)
2914 {
2915 bdg(part,perm1,S_V_I(S_V_I(SMat,i),0L));
2916 }
2917 i++;
2918 }
2919 }while(next(part,part));
2920 }
2921
2922 freeall(perm1);
2923 freeall(perm2);
2924 freeall(part);
2925 freeall(anz_irr);
2926 freeall(m_hoch_n);
2927 freeall(tupel); /* AK 010692 */
2928 return OK;
2929 }
2930 #endif /* DGTRUE */
2931
2932 /***********************************************************************/
2933 /* */
2934 /* Routine: glmndg */
2935 /* Es wird zunaechst die Eingabe erzeugt, */
2936 /* dann wird die Gruppe konstruiert und schliesslich die Basis */
2937 /* berechnet. */
2938 /* Dann wird der Basiswechsel durchgefuehrt. */
2939 /* */
2940 /***********************************************************************/
2941 /* RH 011091 */
2942
2943 #ifdef DGTRUE
2944 #ifdef SABTRUE
glmndg(m,n,M,var)2945 INT glmndg(m,n,M,var) OP m,n,M; INT var;
2946 {
2947 OP Di;
2948 OP D;
2949 OP B;
2950 OP sym_dim;
2951 OP sn_dim;
2952 OP S;
2953 OP SMat;
2954
2955 INT erg = OK; /* AK 230893 */
2956 CTO(INTEGER,"glmndg(1)",m);
2957 CTO(INTEGER,"glmndg(2)",n);
2958
2959 CALLOCOBJECT7(Di,D,B,sym_dim,sn_dim,S,SMat);
2960
2961
2962 erg += input_glmn(m,n,S,SMat,var);
2963 erg += group_gen(S,SMat,D,Di);
2964 erg += _sab(Di,D,B,sym_dim,sn_dim);
2965 erg += glmn_B_W(m,n,B,sym_dim,sn_dim,M);
2966
2967 FREEALL7(Di,D,B,sym_dim,sn_dim,S,SMat);
2968 ENDR("glmndg");
2969 }
2970 #endif /* SABTRUE */
2971 #endif /* DBTRUE */
2972
2973 /***********************************************************************/
2974 /* */
2975 /* Routine: glmn_B_W */
2976 /* Der Basiswechsel wird durchgefuehrt, die Matrix M wird ein */
2977 /* Vektor aus Matrizen, in denen die irreduziblen Matrixdarstellungen */
2978 /* stehen. */
2979 /* */
2980 /***********************************************************************/
2981 /* RH 011091 */
2982
2983 #ifdef DGTRUE
glmn_B_W(m,n,B,sym_dim,sn_dim,M)2984 static INT glmn_B_W(m,n,B,sym_dim,sn_dim,M)
2985 OP m; OP n; OP B; OP sym_dim;
2986 OP sn_dim; OP M;
2987 {
2988 INT i,j,k,l,s;
2989
2990 OP B_1;
2991 OP polyself;
2992 OP res;
2993 OP x;
2994 OP mneu;
2995 OP count;
2996 OP r;
2997
2998 INT erg = OK; /* AK 230893 */
2999
3000 CTO(VECTOR,"glmn_B_W(4)",sym_dim);
3001
3002 B_1 = callocobject();
3003 polyself= callocobject();
3004 res = callocobject();
3005 x = callocobject();
3006 mneu = callocobject();
3007 count = callocobject();
3008 r = callocobject();
3009
3010 erg += invers(B,B_1);
3011 m_i_i(0L,count);
3012 for(s=0L;s<S_V_LI(sym_dim);++s)
3013 if(S_V_II(sym_dim,s) != 0L)
3014 erg += inc(count);
3015 erg += m_l_v(count,M);
3016 m_i_i(0L,count);
3017 for(i=0L;i<S_V_LI(sym_dim);++i)
3018 if(S_V_II(sym_dim,i) != 0L)
3019 {
3020 erg += m_lh_m(S_V_I(sym_dim,i),S_V_I(sym_dim,i),S_V_I(M,S_I_I(count)));
3021 erg += inc(count);
3022 }
3023 m_i_i(0L,count);
3024 m_i_i(0L,r);
3025
3026 m_i_i(0L,count);
3027 for(s=0L;s<S_V_LI(sym_dim);++s)
3028 if(S_V_II(sym_dim,s) != 0L) m_i_i(S_I_I(count)+1L,count);
3029
3030 m_i_i(0L,count);
3031 for(s=0L;s<S_V_LI(sym_dim);++s)
3032 {
3033 if(S_V_II(sym_dim,s) != 0L)
3034 {
3035 for(i=S_I_I(count);i<S_I_I(count)+S_V_II(sym_dim,s);++i)
3036 {
3037 for(j=S_I_I(count);j<S_I_I(count)+S_V_II(sym_dim,s);++j)
3038 {
3039 m_i_i(0L,mneu);
3040 for(k=0;k<S_M_LI(B);++k)
3041 {
3042 for(l=0L;l<S_M_LI(B);++l)
3043 {
3044 if(!nullp(S_M_IJ(B_1,i,k)) &&
3045 !nullp(S_M_IJ(B,l,j)))
3046 {
3047
3048 erg += mult(S_M_IJ(B_1,i,k),S_M_IJ(B,l,j),res);
3049 erg += get_poly_self(k,l,m,n,polyself);
3050 erg += m_skn_po(polyself,res,NULL,x);
3051 erg += add_apply(x,mneu); /* AK 010692 statt add */
3052 }
3053 }
3054 }
3055 erg += copy(mneu,
3056 S_M_IJ(S_V_I(M,S_I_I(r)),
3057 i-S_I_I(count),j-S_I_I(count)));
3058 }
3059 }
3060 erg += inc(r);
3061 }
3062 erg += mult(S_V_I(sym_dim,s),S_V_I(sn_dim,s),x);
3063 erg += add_apply(x,count); /* AK 010692 statt add */
3064 }
3065
3066 erg += freeall(B_1);
3067 erg += freeall(polyself);
3068 erg += freeall(res);
3069 erg += freeall(x);
3070 erg += freeall(mneu);
3071 erg += freeall(r);
3072 erg += freeall(count);
3073 ENDR("glmn_B_W");
3074 }
3075 #endif /* DGTRUE */
3076
3077