1
2 /* Routinen zur Berechnung von Zykelindikatorpolynomen
3 * Nikolaus Schueler 90/91
4 */
5 #include "def.h"
6 #include "macro.h"
7
8 /* Hier wird eigenes FALSE und TRUE eingefuehrt , weil es sich
9 * hier im Gegensatz zu def.h auf unsigned short bezieht
10 */
11 #define N_FALSE 0
12 #define N_TRUE 1
13
14
15 #define OFFSET 0
16 /* diverse prozeduren zu den zykelindikatorpolynom- programmen */
17
18
19 #ifdef PERMTRUE
20 #ifdef VECTORTRUE
21 #ifdef MATRIXTRUE
22 #ifdef PARTTRUE
23 #ifdef POLYTRUE
24 #define ZYKTRUE
25 #endif
26 #endif
27 #endif
28 #endif
29 #endif
30
31 #ifdef ZYKTRUE
32
33
34 static INT colltypes();
35 static INT strongen();
36 static INT polyasub();
37
38 /* Die beiden naechsten Routinen werden zur Berechnung des
39 * Zykelindikators der zyklischen Gruppe benoetigt
40 */
41
42 /* liefert den gcd (greatest common divisor) von a und b zurueck */
43
gcd(a,b)44 static INT gcd(a,b) INT a,b;
45 /* NS 060891 V1.3 */
46 {
47 INT c=5L; /* AK 240692 ERROR muss initialisiert werden */
48
49 /* a soll immer den groesseren Wert enthalten */
50 if (b > a) { c=b;b=a;a=c; }
51 /* Euklidischer Algorithmus */
52 while(c != 0L)
53 {
54 c=a%b; /* c ist der Rest bei der Division */
55 a=b; b=c;
56 }
57 return(a);
58 }
59
eulerfunc(n)60 static INT eulerfunc(n) INT n;
61 /* NS 060891 V1.3 */
62 {
63 INT i,h=0L;
64
65 if(n == 1L) return(1L);
66 for(i=1L; i < n; i++)
67 {
68 if(gcd(n,i) == 1L)
69 {
70 h++;
71 }
72 }
73 return(h);
74 }
75 /* Berechnet das Zykelindikatorpolynom der Cn */
76
77 #ifdef BRUCHTRUE
zykelind_Cn(l,pol)78 INT zykelind_Cn(l,pol) OP l; OP pol;
79 /* NS 060891 V1.3 */
80 {
81 INT d,li;
82 INT erg = OK;
83 OP b;
84
85 CTO(INTEGER,"zykelind_Cn",l);
86 if (S_I_I(l) < 1L) /* AK 060792 */
87 {
88 error("zykelind_Cn: input < 1");
89 goto endr_ende;
90 }
91
92
93 init(POLYNOM,pol);
94
95 if (einsp(l)) /* AK 060792 */
96 {
97 erg += m_iindex_monom(0L,pol);
98 goto endr_ende;
99 }
100
101 b=callocobject();
102 li=S_I_I(l);
103 for(d=1L; (d <= li) ; d++)
104 if(li%d == 0L)
105 {
106 /* stopf den Koeffizienten in den Typ BRUCH von symchar */
107 /* stopf das Ergebnis in den Typ POLYNOM von symchar */
108 erg += b_skn_po(callocobject(),callocobject(),NULL,b);
109 erg += m_ioiu_b(eulerfunc(d),li,S_PO_K(b));
110 erg += kuerzen(S_PO_K(b));
111 erg += m_il_nv(li,S_PO_S (b));
112 erg += m_i_i(li/d,S_PO_SI(b,d-1L));
113 erg += add_apply(b,pol);
114 }
115 erg += freeall(b);
116 ENDR("zykelind_Cn");
117 }
118 #endif /* BRUCHTRUE */
119
120 /* Berechnet das Zykelindikatorpolynom der Dn */
121 /* Beruht auf dem Programm fuer den Zykelindex der Cn, da fuer
122 * den Zykelindex der Diedergruppe Dn nur ein Summand dazukommt
123 * , ruft zykelind_Cn auf.
124 */
125
zykelind_Dn(l,pol)126 INT zykelind_Dn(l,pol) OP l; OP pol;
127 /* NS 060891 V1.3 */
128 {
129 #ifdef BRUCHTRUE
130 INT len,erg=OK;
131 OP b,halb,hilf;
132
133 CTO(INTEGER,"zykelind_Dn(1)",l);
134 SYMCHECK( S_I_I(l) < 1L,"zykelind_Dn: input < 1");
135
136
137 len=S_I_I(l);
138
139 init(POLYNOM,pol);
140
141
142 if (einsp(l)) /* AK 060792 */
143 return m_iindex_monom(0L,pol);
144
145 /* Berechne den Zykelindiktor der Cn */
146
147 zykelind_Cn(l,pol);
148
149 b=callocobject();
150 halb=callocobject();
151 hilf=callocobject();
152
153 /* Vorfaktor 1/2 */
154 div(pol,cons_zwei,pol);
155
156 /* Anhaengen der zusaetzlichen Summanden */
157
158 b_skn_po(callocobject(),callocobject(),NULL,b);
159 m_l_nv(l,S_PO_S (b));
160
161 /* Wenn m gerade ist ..*/
162
163 if((long)len%2L == 0L)
164 {
165 erg += m_ioiu_b(1L,4L,S_PO_K(b));
166 erg += kuerzen(S_PO_K(b));
167 m_i_i(len/2L,S_PO_SI(b,1L));
168 add_apply(b,pol); /* addiere die zusaetzlichen Summanden */
169 erg += m_ioiu_b(1L,4L,S_PO_K(b));
170 erg += kuerzen(S_PO_K(b));
171 m_i_i(2L,S_PO_SI(b,0L));
172 m_i_i((len-2L)/2L,S_PO_SI(b,1L));
173 erg += add_apply(b,pol); /* addiere die zusaetzlichen Summanden */
174 }
175
176 /* Wenn m ungerade ist .. */
177
178 if(len%2L == 1L)
179 {
180 m_ioiu_b(1L,2L,S_PO_K(b));
181 kuerzen(S_PO_K(b));
182
183 /* y1 in das Polynom eintragen */
184
185 m_i_i(1L,S_PO_SI(b,0L));
186
187 /* y2 hoch (n-1L)/2 in das Polynom eintragen */
188
189 m_i_i(((long)len-1L)/2L,S_PO_SI(b,1L));
190 add_apply(b,pol); /* addiere die zusaetzlichen Summanden */
191 }
192 freeall(b);
193 freeall(halb);
194 freeall(hilf);
195 ENDR("zykelind_Dn");
196 #endif /*BRUCHTRUE*/
197 }
198 /* Berechnet das Zykelindikatorpolynom der An */
199
200
zykelind_An(l,pol)201 INT zykelind_An(l,pol) OP l; OP pol;
202 /* NS 060891 V1.3 */
203 {
204 #ifdef BRUCHTRUE
205 INT i,j,veklen,veklen2;
206
207 OP a;
208 OP hilf;
209 OP k;
210 OP n;
211 OP v;
212 OP party;
213 OP zahl;
214 OP zwisch;
215
216 if (S_O_K(l) != INTEGER) /* AK 060792 */
217 return error("zykelind_An: input not INTEGER");
218 if (S_I_I(l) < 1L) /* AK 060792 */
219 return error("zykelind_An: input < 1");
220
221 if (einsp(l)) /* AK 040692 */
222 {
223 return m_iindex_monom(0L,pol);
224 }
225
226 init(POLYNOM,pol);
227
228 a = callocobject();
229 hilf=callocobject();
230 k=callocobject();
231 n = callocobject();
232 v = callocobject();
233 party = callocobject();
234 zahl = callocobject();
235 zwisch = callocobject();
236
237 b_skn_po(callocobject(),callocobject(),NULL,a);
238
239 /* lasse alle partitionen berechnen */
240
241 makevectorofpart(l,v);
242 veklen=S_V_LI(v);
243 m_l_nv(l,S_PO_S (a)); /* AK 040692 Macro */
244 for(i=0L; i < veklen ; i++)
245 {
246 /* umwandeln in Exponentenschreibweise */
247 t_VECTOR_EXPONENT(S_V_I(v,i),party);
248 /* und umwandeln in ein Monom */
249 copy(S_PA_S(party),S_PO_S (a)); /* AK 040692 Macro */
250
251 veklen2=S_V_LI(S_PO_S (a)); /* AK 040692 Macro */
252 m_i_i(0L,zwisch); /* Variable entleeren */
253
254 for(j=1L; j < veklen2; j+=2L)
255 {
256 /* addiere a[2], a[4], ... auf */
257 add_apply(S_PO_SI(a,j),zwisch); /* AK 040692 statt add */
258 }
259 /* Nur wenn a[2]+a[4]+... ungerade ist, dann gibts einen Koeff-
260 * izienten
261 */
262 if(even(zwisch))
263 {
264 /* Berechnen der Koeffizienten */
265
266 m_i_i(1L,k);
267 for(j=0L; j < veklen2; j++)
268 {
269 fakul(S_PO_SI(a,j),zwisch);
270 mult(k,zwisch,k);
271 m_i_i(j+1L,zahl);
272 hoch(zahl,S_PO_SI(a,j),zwisch);
273 mult(k,zwisch,k);
274 }
275 m_i_i(2L,zwisch);
276 m_ou_b(zwisch,k,S_PO_K(a));
277 kuerzen(S_PO_K(a));
278 add_apply(a,pol); /* AK 040692 statt add */
279 }
280 }
281 freeall(a);
282 freeall(hilf);
283 freeall(k);
284 freeall(n);
285 freeall(party);
286 freeall(v);
287 freeall(zahl);
288 freeall(zwisch);
289 return OK;
290 #else /* BRUCHTRUE */
291 return error("zykelind_An: BRUCH not available");
292 #endif /* BRUCHTRUE */
293 }
294 /*
295 * Berechnet das Zykelindikatorpolynom der Sn
296 */
297
298
zykelind_Sn(l,pol)299 INT zykelind_Sn(l,pol) OP l; OP pol;
300 /* NS 060891 V1.3 */
301 /* AK 300998 V2.0 */
302 /* l and pol may be equal */
303 {
304 INT erg = OK;
305
306 CTO(INTEGER,"zykelind_Sn(1)",l);
307 SYMCHECK(S_I_I(l)<1,"zykelind_Sn(1): input < 1");
308 {
309 INT i,j,veklen,veklen2;
310 OP a,hilf,k,v,party,zahl,zwisch;
311
312 CALLOCOBJECT4(a,k,hilf,party);
313 CALLOCOBJECT3(v,zahl,zwisch);
314
315 erg += b_skn_po(CALLOCOBJECT(),CALLOCOBJECT(),NULL,a);
316
317 erg += makevectorofpart(l,v);
318 veklen=S_V_LI(v);
319 erg += m_l_nv(l,S_PO_S (a));
320 erg += init(POLYNOM,pol);
321 for(i=0L; i < veklen ; i++)
322 {
323 /* umwandeln in Exponentenschreibweise */
324 erg += t_VECTOR_EXPONENT(S_V_I(v,i),party);
325 /* und umwandeln in ein Monom */
326 CLEVER_COPY(S_PA_S(party),S_PO_S(a));
327
328 /* Berechnen der Koeffizienten */
329
330 m_i_i(1,k);
331 veklen2=S_V_LI(S_PO_S (a));
332 for(j=0L; j < veklen2; j++)
333 {
334 erg += fakul(S_V_I(S_PO_S(a),j),zwisch);
335 MULT_APPLY(zwisch,k);
336 M_I_I(j+1,zahl);
337 erg += hoch(zahl,S_V_I(S_PO_S(a),j),zwisch);
338 MULT_APPLY(zwisch,k);
339 }
340 erg += invers(k,S_PO_K(a));
341 ADD_APPLY(a,pol);
342 }
343 FREEALL4(a,k,hilf,party);
344 FREEALL3(v,zahl,zwisch);
345 }
346 ENDR("zykelind_Sn");
347 }
348
349 /* Hier folgen die Routinen fuer die Berechnung der Zykel-
350 * indizes von beliebigen Permutationsgruppen, die durch
351 * erzeugende Permutationen gegeben sind
352 */
353 /* routinen zur berechnung eines starken Erzeugers nach Hoffmann
354 */
355
356
357 struct treecomp {
358 unsigned short atr; /* 0=FALSE, 1=TRUE */
359 INT gen;
360 INT point;
361 };
362
stabilizer(i,vec,stabi)363 static INT stabilizer(i,vec,stabi) INT i; OP vec; OP stabi;
364 /* NS 060891 V1.3 */
365 /* AK 101291 vec ist VECTOR of PERMUTATION
366 stabi wird VECTOR of PERMUTATION */
367 {
368 unsigned short is_stab;
369 INT j,k; /* Schleifenzaehler */
370 INT veclen;
371
372 m_il_v(0L,stabi); /* stabi is freed first */
373
374 veclen=S_V_LI(vec);
375 for(j=0L; j < veclen; j++)
376 {
377 is_stab=N_TRUE;
378 for(k=0 ; k < i-1L; k++)
379 {
380 if(S_P_II(S_V_I(vec,j),k) != k+1L)
381 {
382 is_stab=N_FALSE;
383 break;
384 }
385 }
386 if(is_stab)
387 {
388 inc(stabi);
389 copy(S_V_I(vec,j),S_V_I(stabi,S_V_LI(stabi)-1L));
390 }
391 }
392 return OK;
393 }
394
updatemat(degree,i,tree,stabi,repma)395 static INT updatemat(degree,i,tree,stabi,repma)
396 INT degree,i; struct treecomp *tree; OP stabi,repma;
397 /* NS 060891 V1.3 */
398 {
399 INT l,next;
400 OP id=callocobject();
401 OP ob_degree=callocobject();
402 OP operm=callocobject();
403
404 m_i_i(degree,ob_degree);
405 first_permutation(ob_degree,id);
406
407 /* Untersuche die Bahn von i, fuer jedes k aus der Bahn von
408 * von i mache einen Eintrag i,k in der Repraesentations
409 * matrix und zwar das Produkt aller Erzeuger, die i
410 * sukzessive nach k bewegt haben, dazu benutzte alle
411 * Eintrage in tree.gen auf dem Weg zurueck von k nach i (=wort)
412 * und multipliziere die entsprechenden Erzeuger
413 */
414 for(l=i; l < degree; l++)
415 {
416 if(tree[l].atr == N_TRUE)
417 {
418 /* suche rueckwaerts */
419 next=l;
420 copy(id,operm);
421 while(tree[next].point)
422 {
423 mult(operm,S_V_I(stabi,tree[next].gen),operm);
424 next=(tree[next].point)-1L;
425 }
426 copy(operm,S_M_IJ(repma,i-1L,l));
427 }
428 }
429 freeall(id);
430 freeall(ob_degree);
431 freeall(operm);
432 return OK;
433 }
434
435
436 /* sift() sieht nach, ob perm in der Repraesentationsmatrix
437 * repma enthalten ist. sift() wird von strongen benoetigt.
438 */
439
sift(degree,insrow,perm,repma)440 static INT sift(degree,insrow,perm,repma) INT degree, *insrow; OP perm,repma;
441 /* NS 060891 V1.3 */
442 {
443 register unsigned short ismember=N_TRUE;
444 INT i=0L,j;
445 OP invperm=callocobject();
446
447 while((i < degree) && ismember)
448 {
449 i++;
450 j=S_P_II(perm,i-1L);
451 if(not EMPTYP (S_M_IJ(repma,i-1L,j-1L)))
452 {
453 invers(S_M_IJ(repma,i-1L,j-1L),invperm);
454 mult(invperm,perm,perm);
455 }
456 else
457 {
458 ismember=N_FALSE;
459 copy(perm,S_M_IJ(repma,i-1L,j-1L));
460 /* In diese Zeile wurde die Permutation eingefuegt */
461 (*insrow)=i-1L;
462 break;
463 }
464 }
465 freeall(invperm);
466 return(ismember);
467 }
468
469 /* porbit(degree,i,stabi,tree) berechnet die Bahn eines Punktes
470 * i und gibt
471 * folgenden Baum (tree) zurueck:
472 *
473 *
474 * tree[r].atr: TRUE, wenn r in der Bahn von i liegt, sonst FALSE
475 * tree[r].gen: Nummer k des Generators g[k], der s-> r abb.
476 * tree[r].point: Bahnpunkt s, der von dem generator g[k] nach
477 * r abbgebildet wird, tree[i].point = 0L, um
478 * anzuzeigen, dass i die Wurzel des Baumes ist.
479 *
480 * porbit() wird von strongen() benoetigt.
481 */
482
483
porbit(degree,i,stabi,tree)484 static INT porbit(degree,i,stabi,tree) INT degree,i;
485 OP stabi; /* Stabilisator von 1L,...,i-1 */
486 struct treecomp tree[];
487 /* NS 060891 V1.3 */
488 {
489 INT pos=0L;
490 INT stablen;
491 INT g,j, /* Schleifenzaehler */ r,s; /* Bahnpunkte */
492 INT *points;
493
494 stablen=S_V_LI(stabi);
495 points=(INT *) SYM_malloc(sizeof(INT)*degree+OFFSET);
496 points[pos] = i; /* initialisiere points mit dem Punkt,
497 * dessen Bahn gesucht wird
498 */
499 for(j=0L; j < degree; j++)
500 {
501 tree[j].atr=N_FALSE;
502 tree[i-1].atr=N_TRUE;
503 tree[i-1].point=0L; /* i ist Wurzel des Baumes */
504 }
505
506 while(pos >= 0L)
507 {
508 s=points[pos];
509 points[pos--]=0L;
510 for(g=0L; g < stablen; g++)
511 {
512 /* Bestimme das Bild r von s unter dem g-ten Stabilisator */
513 r=S_P_II(S_V_I(stabi,g),s-1L);
514
515 if(not tree[r-1].atr)
516 {
517 points[++pos]=r;
518 tree[r-1].atr=N_TRUE;
519 tree[r-1].gen=g;
520 tree[r-1].point=s;
521 }
522 }
523 }
524 SYM_free(points);
525 return OK;
526 }
527
528 /* Eigentliche Prozedur zur Berechnung eines starken Erzeugers
529 * beziehungsweise einer Repraesentationsmatrix repma
530 */
531
532 #ifdef MATRIXTRUE
strong_generators(a,b)533 INT strong_generators(a,b) OP a,b;
534 /* AK 290192 */
535 /* a VECTOR of generators
536 b becomes MATRIX of stronggenerators */
537 {
538 INT degree, numgen;
539 INT erg = OK;
540 degree=S_P_LI (S_V_I(a,0L));
541 numgen=S_V_LI(a);
542 erg += m_ilih_m(degree+1L,degree+1L,b);
543 erg += strongen(degree,numgen,a,b);
544 ENDR("strong_generators");
545 }
546
547
548
strongen(degree,numgen,genvec,repma)549 static INT strongen(degree,numgen,genvec,repma)
550 INT degree; /* numgen= Anzahl der Erzeuger wird von
551 strongen an sift() * weitergereicht. */
552 INT numgen;
553 OP genvec; /* Vektor der die Erzeuger einer Gruppe enthaelt */
554 OP repma; /* Repraesentationsmatrix fuer die Gruppe */
555 /* NS 060891 V1.3 */
556 {
557 INT i,j,k,l; /* Schleifenzaehler */
558 INT row;
559 INT stablen;
560 INT erg = OK; /* AK 290192 */
561 struct treecomp *tree;
562
563 OP id=callocobject();
564 OP queue=callocobject();
565 OP ob_degree=callocobject();
566 OP perm_eins=callocobject();
567 /* stabi ist Vektor von Permutationen.
568 */
569 OP stabi=callocobject();
570 OP strgset=callocobject();
571
572 tree=(struct treecomp*) SYM_malloc(degree*sizeof(struct treecomp)+OFFSET);
573 m_i_i(degree,ob_degree);
574 erg +=first_permutation(ob_degree,id);
575
576
577 /* Stabilisator ist am Anfang der ganze Erzeuger */
578 erg +=m_il_v(numgen,strgset);
579 erg +=m_il_v(0L,stabi);
580 for(k=0L; k < numgen; k++)
581 {
582 erg +=copy(S_V_I(genvec,k),S_V_I(strgset,k));
583 }
584
585 /* Diagonale der Repraesentationsmatrix mit id besetzen */
586 for(k=0L; k < degree; k++)
587 erg +=copy(id,S_M_IJ(repma,k,k));
588
589 for(i=1L; i <= degree; i++)
590 {
591
592 erg +=stabilizer(i,strgset,stabi);
593 erg +=porbit(degree,i,stabi,tree);
594 erg +=updatemat(degree,i,tree,stabi,repma);
595 }
596 m_il_v(0L,queue);
597
598 for(i=1L; i <= degree; i++)
599 {
600 erg +=stabilizer(i,strgset,stabi);
601 stablen=S_V_LI(stabi);
602
603
604 for(l=0L; l < stablen; l++)
605 {
606 /* for(k=i-1L; k <= degree; k++) */
607 for(k=i-1L; k < degree; k++)
608 /* statt <= degree < degree */
609 if(not EMPTYP(S_M_IJ(repma,i-1L,k)))
610 {
611 erg +=mult(S_V_I(stabi,l),S_M_IJ(repma,i-1L,k),perm_eins);
612 erg +=inc(queue);
613 erg +=copy(perm_eins,S_V_I(queue,S_V_LI(queue)-1L));
614 }
615 }
616
617 while(S_V_LI(queue))
618 {
619 erg +=copy(S_V_I(queue,S_V_LI(queue)-1L),perm_eins);
620 erg +=dec(queue);
621 if(not sift(degree,&row,perm_eins,repma)) /* ismember == 0 */
622 {
623 erg +=inc(strgset);
624 erg +=copy(perm_eins,S_V_I(strgset,S_V_LI(strgset)-1L));
625
626 for(j=1L; j <= row; j++)
627
628 {
629 erg +=stabilizer(j,strgset,stabi);
630 erg +=porbit(degree,j,stabi,tree);
631 erg +=updatemat(degree,j,tree,stabi,repma);
632
633 /*for(l=j-1L; l <= degree; l++)*/
634 for(l=j-1L; l < degree; l++)
635 /* < degree statt <= degree */
636 {
637 if(not EMPTYP(S_M_IJ(repma,j-1L,l)))
638 {
639 erg +=mult(perm_eins,S_M_IJ(repma,j-1L,l),perm_eins);
640 erg +=inc(queue);
641 erg +=copy(perm_eins,S_V_I(queue,S_V_LI(queue)-1L));
642 }
643 }
644 }
645 } /* end if */
646 } /* end while */
647 } /* AK end for 110292 */
648 erg +=freeall(id);
649 erg +=freeall(queue);
650 erg +=freeall(ob_degree);
651 erg +=freeall(perm_eins);
652 erg +=freeall(stabi);
653 erg +=freeall(strgset);
654 SYM_free(tree);
655 return erg;
656 } /* end all :-) */
657 #endif /* MATRIXTRUE */
658
659
recu(degree,start,numnontriv,ztvec,numztvec,perm,repma)660 static INT recu(degree,start,numnontriv,ztvec,numztvec,perm,repma)
661 INT degree; INT start; INT numnontriv;
662 OP ztvec; OP numztvec; OP perm; OP repma;
663 /* NS 060891 V1.3 */
664 {
665 INT i,j;
666
667 OP saveperm=callocobject();
668
669
670 if(start == numnontriv-1L)
671 {
672 for(i=start; i< degree; i++)
673 {
674 if(not EMPTYP(S_M_IJ(repma,start,i)))
675 {
676 mult(perm,S_M_IJ(repma,start,i),saveperm);
677 colltypes(saveperm,ztvec,numztvec);
678 }
679 }
680
681
682 }
683 else
684 for(j=start; j < degree; j++)
685 {
686 if(not EMPTYP(S_M_IJ(repma,start,j)))
687 {
688 mult(perm,S_M_IJ(repma,start,j),saveperm);
689 recu(degree,start+1L,numnontriv,ztvec,numztvec,saveperm,repma);
690 }
691 }
692
693 freeall(saveperm);
694 return OK;
695 }
696
697
698
699
700
callrecu(grad,ztvec,numztvec,repma)701 static INT callrecu(grad,ztvec,numztvec,repma) INT grad; OP ztvec,numztvec,repma;
702 /* NS 060891 V1.3 */
703 {
704 unsigned short trivrow;
705 INT i,j;
706 INT numnontriv=1L; /* AK 021291 */
707 OP id=callocobject();
708 OP ob_grad=callocobject();
709 OP perm=callocobject();
710
711 /* Weil die unteren Zeilen der Matrix meistens bis auf
712 * die Identitaet in der Diagonalen leer sind, wird
713 * hier erstmal festgestellt, ab wo die Matrix leer
714 * ist
715 */
716
717 for(i=grad-1L; i > 0L; i--)
718 {
719 trivrow=N_TRUE;
720
721 for(j=i+1L; j < grad; j++)
722 {
723 if(not EMPTYP(S_M_IJ(repma,i,j)))
724 {
725 trivrow=N_FALSE;
726 break;
727 }
728 }
729 if(trivrow == N_FALSE)
730 {
731 numnontriv=i+1L;
732 break;
733 }
734 }
735 m_i_i(grad,ob_grad);
736 first_permutation(ob_grad,id);
737 copy(id,perm);
738 recu(grad,0L,numnontriv,ztvec,numztvec,perm,repma);
739 freeall(id);
740 freeall(ob_grad);
741 freeall(perm);
742 return OK;
743 } /* end all */
744
745
746
colltypes(perm,ztvec,numztvec)747 static INT colltypes(perm,ztvec,numztvec) OP perm, ztvec, numztvec;
748 /* NS 060891 V1.3 */
749 {
750 INT i;
751 INT ztveclen;
752 OP ztperm=callocobject();
753 OP expztperm=callocobject();
754
755 ztveclen=S_V_LI(ztvec);
756 zykeltyp(perm,ztperm);
757 t_VECTOR_EXPONENT(ztperm,expztperm);
758
759 for(i=0L; i < ztveclen; i++)
760 {
761 if(comp(expztperm,S_V_I(ztvec,i)) == 0L)
762 {
763 inc(S_V_I(numztvec,i));
764 goto ende;
765 }
766 }
767
768 inc(ztvec);
769 copy(expztperm,S_V_I(ztvec,S_V_LI(ztvec)-1L));
770 inc(numztvec);
771 m_i_i(1L,S_V_I(numztvec,S_V_LI(numztvec)-1L));
772 ende:
773 freeall(ztperm);
774 freeall(expztperm);
775 return OK;
776 }
777 /* berechnet das Zykelindikatorpolynom einer beliebigen Permutations-
778 * gruppe, benutzt dazu die uebergebenen Vektoren, die die
779 * Zykeltypen (expztvec), bzw deren Anzahlen (numztvec)
780 * enthalten.
781 */
782
783 #ifdef BRUCHTRUE
zykelind_arb_co(expztvec,numztvec,pol)784 static INT zykelind_arb_co(expztvec,numztvec,pol)
785 OP expztvec;
786 OP numztvec;
787 OP pol; /* enhaelt nach Ablauf der Routine das Zykelindikator-
788 * polynom (noch nicht Polyasubstituiert)
789 */
790 /* NS 060891 V1.3 */
791 {
792 INT i,order,veklen;
793 INT erg = OK;
794
795 OP a,hilf,k,party,zahl,zwisch,zykeltypvec;
796 OP ak_order;
797
798 a = CALLOCOBJECT();
799 hilf=CALLOCOBJECT();
800 k=CALLOCOBJECT();
801 party = CALLOCOBJECT();
802 zahl = CALLOCOBJECT();
803 zwisch = CALLOCOBJECT();
804 zykeltypvec = CALLOCOBJECT();
805 ak_order = CALLOCOBJECT();
806
807 SYM_sum(numztvec,ak_order); /* AK 060295 */
808
809
810
811 b_skn_po(CALLOCOBJECT(),CALLOCOBJECT(),NULL,a);
812
813 veklen=S_V_LI(expztvec);
814 m_il_nv(5,S_PO_S(a));
815 init(POLYNOM,pol);
816 for(i=0L; i < veklen ; i++)
817 {
818 COPY(S_PA_S(S_V_I(expztvec,i)),S_PO_S (a));
819 /* m_i_i(order,k); */
820 m_ou_b(S_V_I(numztvec,i),ak_order,S_PO_K(a));
821 /* m_ou_b(S_V_I(numztvec,i),k,S_PO_K(a)); */
822 kuerzen(S_PO_K(a));
823 add_apply(a,pol);
824 }
825 FREEALL3(a,hilf,k);
826 FREEALL5(ak_order,party,zahl,zwisch,zykeltypvec);
827 ENDR("zykelind_arb_co:internal routine");
828 }
829 #endif /* BRUCHTRUE */
830
831 /*
832 * Die Funktion zykelind_arb fasst die Funktionen strongen,
833 * callrecu und zykelind_arb_co zusammen. Eingabeparameter
834 * der Vektor von Permutationen genvec, er enthaelt die Erzeuger
835 * der Gruppe, in pol wird dann das Zykelindikator polynom
836 * geliefert.
837 */
838
zykelind_arb(genvec,pol)839 INT zykelind_arb(genvec,pol) OP genvec; OP pol;
840 /* NS 060891 V1.3 */
841 /* AK 180998: now the generating permutations may have different
842 degree */
843 {
844 #ifdef BRUCHTRUE
845 INT degree,i,j;
846 INT numgen;
847 INT erg = OK;
848 OP mat=callocobject();
849 OP numztvec=callocobject();
850 OP ztvec=callocobject();
851 OP axl = callocobject();
852 OP mygenvec = callocobject();
853
854 erg += m_l_v(cons_null,numztvec);
855 erg += m_l_v(cons_null,ztvec);
856
857 /* degree und numgen bestimmen */
858 degree=S_P_LI (S_V_I(genvec,0L));
859 for (i=1;i<S_V_LI(genvec);i++) /* AK 180998 */
860 if (S_P_LI (S_V_I(genvec,i)) > degree)
861 degree=S_P_LI (S_V_I(genvec,i));
862 numgen=S_V_LI(genvec);
863 erg += m_il_v(numgen,mygenvec); /* AK 180998 */
864 for (i=0;i<numgen;i++) /* AK 180998 */
865 {
866 if (degree==S_P_LI (S_V_I(genvec,i)))
867 erg += copy_permutation(S_V_I(genvec,i), S_V_I(mygenvec,i));
868 else {
869 erg += m_il_p(degree,S_V_I(mygenvec,i));
870 for (j=0;j<S_P_LI(S_V_I(genvec,i));j++)
871 M_I_I(S_P_II(S_V_I(genvec,i),j),S_P_I(S_V_I(mygenvec,i),j));
872 for(;j<S_P_LI(S_V_I(mygenvec,i));j++)
873 M_I_I(j+1,S_P_I(S_V_I(mygenvec,i),j));
874 }
875 }
876 erg += m_i_i(degree+1L,axl);
877 erg += m_lh_m(axl,axl,mat);
878 erg += strongen(degree,numgen,mygenvec,mat);
879 erg += callrecu(degree,ztvec,numztvec,mat);
880 erg += zykelind_arb_co(ztvec,numztvec,pol);
881 FREEALL5(axl,numztvec,ztvec,mygenvec,mat);
882 ENDR("zykelind_arb");
883 #else /* BRUCHTRUE */
884 return error("zkelind_arb: BRUCH not available");
885 #endif /* BRUCHTRUE */
886 }
887 /* routine zur Polyasubstitution */
888
889
polya_n_sub(p,n,e)890 INT polya_n_sub(p,n,e) OP p,n,e;
891 /* AK 060792 */
892 {
893 return polyasub(S_I_I(n),p,e);
894 }
895
polyasub(numcol,pol,pattpol)896 static INT polyasub(numcol,pol,pattpol)
897 INT numcol; /* Anzahl der Farben */
898 OP pol; /* Zykelindikatorpolynom */
899 OP pattpol; /* Ergebnis: Musterpolynom */
900 /* NS 060891 V1.3 */
901 {
902 INT i,j; /* Schleifenzaehler */
903 INT degree;
904 INT erg = OK;
905 OP colvec=callocobject();
906 OP compcolpol=callocobject();
907 OP hpol=callocobject();
908 OP exp=callocobject();
909
910 /* Farbenvector herstellen:
911 * Monom abc...... -> [11...1][22...2]...[nn...n]
912 * wobei die Laenge von [xx...x] gleich der Anzahl der
913 * Farben ist und n die Maechtigkeit der Menge X,
914 * auf der die Gruppe operiert.
915 */
916 erg += numberofvariables(pol,exp); /* AK 211194 */
917 degree=S_I_I(exp);
918 erg += m_il_v(degree,colvec);
919 erg += b_skn_po(callocobject(),callocobject(),NULL,compcolpol);
920 erg += b_skn_po(callocobject(),callocobject(),NULL,hpol);
921
922 FREESELF(pattpol);
923
924 for(i=0L; i < degree; i++)
925 {
926 erg += init(POLYNOM,compcolpol);
927
928 for(j=0L; j < numcol; j++)
929 {
930 if (not EMPTYP(hpol))
931 erg += freeself(hpol);
932 erg += m_iindex_iexponent_monom(j,i+1L,hpol);
933 erg += add_apply(hpol,compcolpol);
934 }
935 erg += copy(compcolpol,S_V_I(colvec,i));
936 }
937 erg += eval_polynom(pol,colvec,pattpol);
938
939 FREEALL4(colvec,compcolpol,hpol,exp);
940 ENDR("zyk:internal function polyasub");
941 }
942
943 /* Algorithmus von dimino, berechnet eine Liste mit allen
944 * Elementen einer Gruppe aus den erzeugenden Permutationen.
945 */
946
dimino(elm)947 INT dimino(elm) OP elm;
948 /* enthaelt am Anfang die Erzeuger der Gruppe, nach
949 * Ablauf der Routine dann alle Gruppenelemente */
950 /* NS 060891 V1.3 */
951 {
952 INT i,j,k,
953 cosetlen, elt_not_elm, numgen,
954 order=0L,
955 rep_pos, s_count, si_not_elm;
956 INT erg = OK;
957
958
959 OP elt,g,genvec,id;
960 CTO(VECTOR,"dimino(1)",elm);
961
962 elt=callocobject(); /* Hilfsvariable fuer Test auf Enhaltensein
963 * in elm
964 */
965 g=callocobject(); /* eine Permutation */
966 genvec=callocobject(); /* enhaelt die Erzeuger */
967 id=callocobject(); /* die identische Permutation */
968
969 numgen=S_V_LI(elm);
970 erg += m_il_v(numgen,genvec);
971 /* Kopiere die Erzeuger in einen eigenen Vektor genvec */
972 for(i=0L; i < numgen; i++)
973 COPY(S_V_I(elm,i),S_V_I(genvec,i));
974 eins(S_V_I(genvec,0),id);
975
976 /* Liste der Elemente anlegen, laenge ist erstmal = 1 */
977 erg += m_il_v(1L,elm);
978
979 /* Spezialfall G= <S1> */
980
981 COPY(id,S_V_I(elm,order)); /* 1. Element ist id */
982 COPY(S_V_I(genvec,0L),g); /* g:=s1 */
983 // while(comp(g,id)) /* Solange g ungleich id */
984 while(not einsp(g)) /* Solange g ungleich id */
985 {
986 /* Elementevektor muss jedesmal erst um 1 verlaengert werden */
987 INC(elm);
988 ++order; /* AK 060891 */
989 COPY(g,S_V_I(elm,order)); /* elm[order]=g */
990 CLEVER_MULT(S_V_I(elm,order),S_V_I(genvec,0L),g); /* g:=g*s1 */
991 }
992
993 /* Laenge der Nebenklassen feststellen, muss man nur einmal machen,
994 * da alle Nebenklassen gleiche Laenge haben
995 */
996 cosetlen=S_V_LI(elm);
997
998 /* Falls es mehr als einen Erzeuger gibt */
999 for(i=1L; i < numgen; i++)
1000 {
1001 si_not_elm=1L;
1002 for(k=0L; k <= order; k++) /* s(i) in elm ? */
1003 if((si_not_elm=comp(S_V_I(genvec,i),S_V_I(elm,k))) == 0L)
1004 break;
1005
1006 /* Wenn s[i] nicht in elm:
1007 * s[i] und seine Nebenklasse g*s[i]
1008 * zu elm hinzufuegen
1009 */
1010 if(si_not_elm) /* Wenn s(i) nicht in elm */
1011 {
1012 /* s[i] hinzufuegen */
1013 /* Elementevektor muss jedesmal erst um 1 verlaengert werden */
1014 inc(elm);
1015 ++order;
1016 COPY(S_V_I(genvec,i),S_V_I(elm,order));
1017 /* Nebenklasse zu elm hinzufuegen */
1018 for(j=1L; j < cosetlen; j++)
1019 {
1020 /* ++order,elm[order]:=elm[j]*s[i] */
1021 /* Elementevektor muss jedesmal erst
1022 um 1 verlaengert werden */
1023 inc(elm);
1024 ++order;
1025 MULT(S_V_I(elm,j),S_V_I(genvec,i),S_V_I(elm,order));
1026 } /* end for */
1027
1028 rep_pos=cosetlen;
1029
1030 do {
1031 for(s_count=0L; s_count <= i; s_count++)
1032 {
1033 /* elt=elm[rep_pos]*s[s_count] */
1034 MULT(S_V_I(elm,rep_pos),
1035 S_V_I(genvec,s_count),elt);
1036
1037 elt_not_elm=1L;
1038 for(k=0L; k <= order; k++)
1039 /* elt in elm ? */
1040 if((elt_not_elm=comp(elt,S_V_I(elm,k))) == 0L)
1041 break;
1042
1043 /* Wenn elt nicht in elm:
1044 * elt und seine Nebenklasse g*elt
1045 * zu elm hinzufuegen
1046 */
1047 if(elt_not_elm)
1048 {
1049 /* elt hinzufuegen */
1050 /* Elementevektor muss jedesmal erst
1051 * um 1 verlaengert werden
1052 */
1053 INC(elm);
1054 ++order;
1055 COPY(elt,S_V_I(elm,order));
1056 /* Nebenklasse zu elm hinzufuegen */
1057 for(j=1L; j < cosetlen; j++)
1058 {
1059 /* ++order,elm[order]:=elm[j]*s[i] */
1060 /* Elementevektor muss jedesmal erst
1061 * um 1 verlaengert werden
1062 */
1063 INC(elm);
1064 ++order;
1065 MULT(S_V_I(elm,j),elt,S_V_I(elm,order));
1066 } /* end for */
1067 } /* end if */
1068 } /* end for */
1069 rep_pos+=cosetlen;
1070 } while(rep_pos <= order);
1071 } /* end if */
1072 cosetlen=order+1L;
1073 } /* end for */
1074
1075 FREEALL4(elt,g,genvec,id);
1076 CTO(VECTOR,"dimino(1e)",elm);
1077 ENDR("dimino");
1078 }
1079
1080
grf_arb(gr,n,res)1081 INT grf_arb(gr,n,res) OP gr,n,res;
1082 /* AK 220998 V2.0 */
1083 /* AK 091204 V3.0 */
1084 {
1085 INT erg = OK;
1086 CTO(INTEGER,"grf_arb(2)",n);
1087 CTO(VECTOR,"grf_arb(1)",gr);
1088 CE3(gr,n,res,grf_arb);
1089 {
1090 OP zw;
1091 zw = CALLOCOBJECT();
1092 erg += zykelind_arb(gr,zw);
1093 erg += polya_n_sub(zw,n,res);
1094 FREEALL(zw);
1095 }
1096 ENDR("grf_arb");
1097 }
1098
1099
grf_Sn(gr,n,res)1100 INT grf_Sn(gr,n,res) OP gr,n,res;
1101 /* AK 220998 V2.0 */
1102 {
1103 INT erg = OK;
1104 CTO(INTEGER,"grf_Sn",n);
1105 CTO(INTEGER,"grf_Sn",gr);
1106 CE3(gr,n,res,grf_Sn);
1107 {
1108 OP zw;
1109 zw = callocobject();
1110 erg += zykelind_Sn(gr,zw);
1111 erg += polya_n_sub(zw,n,res);
1112 erg += freeall(zw);
1113 }
1114 ENDR("grf_Sn");
1115 }
1116
grf_An(gr,n,res)1117 INT grf_An(gr,n,res) OP gr,n,res;
1118 /* AK 220998 V2.0 */
1119 {
1120 OP zw;
1121 INT erg = OK;
1122 CTO(INTEGER,"grf_An",n);
1123 CTO(INTEGER,"grf_An",gr);
1124 CE3(gr,n,res,grf_An);
1125 zw = callocobject();
1126 erg += zykelind_An(gr,zw);
1127 erg += polya_n_sub(zw,n,res);
1128 erg += freeall(zw);
1129 ENDR("grf_An");
1130 }
1131
grf_Cn(gr,n,res)1132 INT grf_Cn(gr,n,res) OP gr,n,res;
1133 /* AK 220998 V2.0 */
1134 {
1135 OP zw;
1136 INT erg = OK;
1137 CTO(INTEGER,"grf_Cn",n);
1138 CTO(INTEGER,"grf_Cn",gr);
1139 CE3(gr,n,res,grf_Cn);
1140 zw = callocobject();
1141 erg += zykelind_Cn(gr,zw);
1142 erg += polya_n_sub(zw,n,res);
1143 erg += freeall(zw);
1144 ENDR("grf_Cn");
1145 }
1146
grf_Dn(gr,n,res)1147 INT grf_Dn(gr,n,res) OP gr,n,res;
1148 /* AK 220998 V2.0 */
1149 {
1150 OP zw;
1151 INT erg = OK;
1152 CTO(INTEGER,"grf_Dn",n);
1153 CTO(INTEGER,"grf_Dn",gr);
1154 CE3(gr,n,res,grf_Dn);
1155 zw = callocobject();
1156 erg += zykelind_Dn(gr,zw);
1157 erg += polya_n_sub(zw,n,res);
1158 erg += freeall(zw);
1159 ENDR("grf_Dn");
1160 }
1161
1162
1163
1164
no_orbits_arb(a,b,c)1165 INT no_orbits_arb(a,b,c) OP a,b,c;
1166 /* AK 071098 V2.0 */
1167 {
1168 OP d,e;
1169 OP z;
1170 INT erg = OK;
1171 CE3(a,b,c,no_orbits_arb);
1172 d = callocobject();
1173 e = callocobject();
1174 erg += zykelind_arb(a,d);
1175 z = d;
1176 erg += m_i_i(0,c);
1177 while (z!=NULL)
1178 {
1179 erg += SYM_sum(S_PO_S(z),e);
1180 erg += hoch(b,e,e);
1181 erg += mult_apply(S_PO_K(z),e);
1182 erg += add_apply(e,c);
1183 z = S_PO_N(z);
1184 }
1185
1186 erg += freeall(d);
1187 erg += freeall(e);
1188 ENDR("no_orbits_arb");
1189 }
1190
1191 #endif /* ZYKTRUE */
1192
1193