Line data Source code
1 : /* Copyright (C) 2000-2003 The PARI group.
2 :
3 : This file is part of the PARI/GP package.
4 :
5 : PARI/GP is free software; you can redistribute it and/or modify it under the
6 : terms of the GNU General Public License as published by the Free Software
7 : Foundation; either version 2 of the License, or (at your option) any later
8 : version. It is distributed in the hope that it will be useful, but WITHOUT
9 : ANY WARRANTY WHATSOEVER.
10 :
11 : Check the License for details. You should have received a copy of it, along
12 : with the package; see the file 'COPYING'. If not, write to the Free Software
13 : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
14 :
15 : #include "pari.h"
16 : #include "paripriv.h"
17 :
18 : /*************************************************************************/
19 : /** **/
20 : /** Routines for handling VEC/COL **/
21 : /** **/
22 : /*************************************************************************/
23 : int
24 1841 : vec_isconst(GEN v)
25 : {
26 1841 : long i, l = lg(v);
27 : GEN w;
28 1841 : if (l==1) return 1;
29 1841 : w = gel(v,1);
30 6342 : for(i=2; i<l; i++)
31 5768 : if (!gequal(gel(v,i), w)) return 0;
32 574 : return 1;
33 : }
34 :
35 : int
36 17424 : vecsmall_isconst(GEN v)
37 : {
38 17424 : long i, l = lg(v);
39 : ulong w;
40 17424 : if (l==1) return 1;
41 17424 : w = uel(v,1);
42 30092 : for(i=2; i<l; i++)
43 24267 : if (uel(v,i) != w) return 0;
44 5825 : return 1;
45 : }
46 :
47 : /* Check if all the elements of v are different.
48 : * Use a quadratic algorithm. Could be done in n*log(n) by sorting. */
49 : int
50 0 : vec_is1to1(GEN v)
51 : {
52 0 : long i, j, l = lg(v);
53 0 : for (i=1; i<l; i++)
54 : {
55 0 : GEN w = gel(v,i);
56 0 : for(j=i+1; j<l; j++)
57 0 : if (gequal(gel(v,j), w)) return 0;
58 : }
59 0 : return 1;
60 : }
61 :
62 : GEN
63 98084 : vec_insert(GEN v, long n, GEN x)
64 : {
65 98084 : long i, l=lg(v);
66 98084 : GEN V = cgetg(l+1,t_VEC);
67 711340 : for(i=1; i<n; i++) gel(V,i) = gel(v,i);
68 98084 : gel(V,n) = x;
69 471646 : for(i=n+1; i<=l; i++) gel(V,i) = gel(v,i-1);
70 98084 : return V;
71 : }
72 : /*************************************************************************/
73 : /** **/
74 : /** Routines for handling VECSMALL **/
75 : /** **/
76 : /*************************************************************************/
77 : /* Sort v[0]...v[n-1] and put result in w[0]...w[n-1].
78 : * We accept v==w. w must be allocated. */
79 : static void
80 181375578 : vecsmall_sortspec(GEN v, long n, GEN w)
81 : {
82 181375578 : pari_sp ltop=avma;
83 181375578 : long nx=n>>1, ny=n-nx;
84 : long m, ix, iy;
85 : GEN x, y;
86 181375578 : if (n<=2)
87 : {
88 107442454 : if (n==1)
89 23124836 : w[0]=v[0];
90 84317618 : else if (n==2)
91 : {
92 88673266 : long v0=v[0], v1=v[1];
93 88673266 : if (v0<=v1) { w[0]=v0; w[1]=v1; }
94 2208537 : else { w[0]=v1; w[1]=v0; }
95 : }
96 107442454 : return;
97 : }
98 73933124 : x=new_chunk(nx); y=new_chunk(ny);
99 80510669 : vecsmall_sortspec(v,nx,x);
100 84249563 : vecsmall_sortspec(v+nx,ny,y);
101 392605563 : for (m=0, ix=0, iy=0; ix<nx && iy<ny; )
102 300981751 : if (x[ix]<=y[iy])
103 253567754 : w[m++]=x[ix++];
104 : else
105 47413997 : w[m++]=y[iy++];
106 94094489 : for(;ix<nx;) w[m++]=x[ix++];
107 338909840 : for(;iy<ny;) w[m++]=y[iy++];
108 91623812 : set_avma(ltop);
109 : }
110 :
111 : static long
112 34753450 : vecsmall_sort_max(GEN v)
113 : {
114 34753450 : long i, l = lg(v), max = -1;
115 120000256 : for (i = 1; i < l; i++)
116 118290771 : if (v[i] > max) { max = v[i]; if (max >= l) return -1; }
117 21068514 : else if (v[i] < 0) return -1;
118 1709485 : return max;
119 : }
120 : /* assume 0 <= v[i] <= M. In place. */
121 : void
122 1843058 : vecsmall_counting_sort(GEN v, long M)
123 : {
124 : pari_sp av;
125 : long i, j, k, l;
126 : GEN T;
127 1843058 : if (M == 0) return;
128 1843058 : av = avma; T = new_chunk(M + 1); l = lg(v);
129 8978409 : for (i = 0; i <= M; i++) T[i] = 0;
130 7193256 : for (i = 1; i < l; i++) T[v[i]]++; /* T[j] is # keys = j */
131 8978421 : for (j = 0, k = 1; j <= M; j++)
132 12485570 : for (i = 1; i <= T[j]; i++) v[k++] = j;
133 1843057 : set_avma(av);
134 : }
135 : /* not GC-clean, suitable for gerepileupto */
136 : GEN
137 16115 : vecsmall_counting_uniq(GEN v, long M)
138 : {
139 16115 : long i, k, l = lg(v);
140 : GEN T, U;
141 16115 : if (l == 1) return cgetg(1, t_VECSMALL);
142 16115 : if (M == 0) return mkvecsmall(0);
143 16115 : if (l == 2) return leafcopy(v);
144 15982 : U = new_chunk(M + 2);
145 15982 : T = U+1; /* allows to rewrite result over T also if T[0] = 1 */
146 124248 : for (i = 0; i <= M; i++) T[i] = 0;
147 201760 : for (i = 1; i < l; i++) T[v[i]] = 1;
148 124248 : for (i = 0, k = 1; i <= M; i++)
149 108266 : if (T[i]) U[k++] = i;
150 15982 : U[0] = evaltyp(t_VECSMALL) | _evallg(k); return U;
151 : }
152 : GEN
153 10631 : vecsmall_counting_indexsort(GEN v, long M)
154 : {
155 : pari_sp av;
156 10631 : long i, l = lg(v);
157 : GEN T, p;
158 10631 : if (M == 0 || l <= 2) return identity_zv(l - 1);
159 10617 : p = cgetg(l, t_VECSMALL); av = avma; T = new_chunk(M + 1);
160 48859 : for (i = 0; i <= M; i++) T[i] = 0;
161 15317063 : for (i = 1; i < l; i++) T[v[i]]++; /* T[j] is # keys = j */
162 38242 : for (i = 1; i <= M; i++) T[i] += T[i-1]; /* T[j] is # keys <= j */
163 15317063 : for (i = l-1; i > 0; i--) { p[T[v[i]]] = i; T[v[i]]--; }
164 10617 : return gc_const(av, p);
165 : }
166 :
167 : /* in place sort */
168 : void
169 39161312 : vecsmall_sort(GEN v)
170 : {
171 39161312 : long n = lg(v) - 1, max;
172 39161312 : if (n <= 1) return;
173 32438962 : if ((max = vecsmall_sort_max(v)) >= 0)
174 1843058 : vecsmall_counting_sort(v, max);
175 : else
176 30493448 : vecsmall_sortspec(v+1, n, v+1);
177 : }
178 :
179 : /* cf gen_sortspec */
180 : static GEN
181 7914453 : vecsmall_indexsortspec(GEN v, long n)
182 : {
183 : long nx, ny, m, ix, iy;
184 : GEN x, y, w;
185 7914453 : switch(n)
186 : {
187 148941 : case 1: return mkvecsmall(1);
188 3597815 : case 2: return (v[1] <= v[2])? mkvecsmall2(1,2): mkvecsmall2(2,1);
189 1398281 : case 3:
190 1398281 : if (v[1] <= v[2]) {
191 668801 : if (v[2] <= v[3]) return mkvecsmall3(1,2,3);
192 201891 : return (v[1] <= v[3])? mkvecsmall3(1,3,2)
193 611052 : : mkvecsmall3(3,1,2);
194 : } else {
195 729480 : if (v[1] <= v[3]) return mkvecsmall3(2,1,3);
196 273883 : return (v[2] <= v[3])? mkvecsmall3(2,3,1)
197 828578 : : mkvecsmall3(3,2,1);
198 : }
199 : }
200 2769416 : nx = n>>1; ny = n-nx;
201 2769416 : w = cgetg(n+1,t_VECSMALL);
202 2769415 : x = vecsmall_indexsortspec(v,nx);
203 2769417 : y = vecsmall_indexsortspec(v+nx,ny);
204 33591344 : for (m=1, ix=1, iy=1; ix<=nx && iy<=ny; )
205 30821927 : if (v[x[ix]] <= v[y[iy]+nx])
206 14936175 : w[m++] = x[ix++];
207 : else
208 15885752 : w[m++] = y[iy++]+nx;
209 5284240 : for(;ix<=nx;) w[m++] = x[ix++];
210 5298647 : for(;iy<=ny;) w[m++] = y[iy++]+nx;
211 2769417 : set_avma((pari_sp)w); return w;
212 : }
213 :
214 : /*indirect sort.*/
215 : GEN
216 2386314 : vecsmall_indexsort(GEN v)
217 : {
218 2386314 : long n = lg(v) - 1, max;
219 2386314 : if (n==0) return cgetg(1, t_VECSMALL);
220 2386251 : if ((max = vecsmall_sort_max(v)) >= 0)
221 10631 : return vecsmall_counting_indexsort(v, max);
222 : else
223 2375621 : return vecsmall_indexsortspec(v,n);
224 : }
225 :
226 : /* assume V sorted */
227 : GEN
228 31340 : vecsmall_uniq_sorted(GEN v)
229 : {
230 : long i, j, l;
231 31340 : GEN w = cgetg_copy(v, &l);
232 31340 : if (l == 1) return w;
233 31286 : w[1] = v[1];
234 34259 : for(i = j = 2; i < l; i++)
235 2973 : if (v[i] != w[j-1]) w[j++] = v[i];
236 31286 : stackdummy((pari_sp)(w + l), (pari_sp)(w + j));
237 31286 : setlg(w, j); return w;
238 : }
239 :
240 : GEN
241 16446 : vecsmall_uniq(GEN v)
242 : {
243 16446 : pari_sp av = avma;
244 : long max;
245 16446 : if ((max = vecsmall_sort_max(v)) >= 0)
246 16115 : v = vecsmall_counting_uniq(v, max);
247 : else
248 331 : { v = zv_copy(v); vecsmall_sort(v); v = vecsmall_uniq_sorted(v); }
249 16446 : return gerepileuptoleaf(av, v);
250 : }
251 :
252 : /* assume x sorted */
253 : long
254 0 : vecsmall_duplicate_sorted(GEN x)
255 : {
256 0 : long i,k,l=lg(x);
257 0 : if (l==1) return 0;
258 0 : for (k=x[1],i=2; i<l; k=x[i++])
259 0 : if (x[i] == k) return i;
260 0 : return 0;
261 : }
262 :
263 : long
264 19437 : vecsmall_duplicate(GEN x)
265 : {
266 19437 : pari_sp av=avma;
267 19437 : GEN p=vecsmall_indexsort(x);
268 19437 : long k,i,r=0,l=lg(x);
269 19437 : if (l==1) return 0;
270 26608 : for (k=x[p[1]],i=2; i<l; k=x[p[i++]])
271 7171 : if (x[p[i]] == k) { r=p[i]; break; }
272 19437 : set_avma(av);
273 19437 : return r;
274 : }
275 :
276 : static int
277 54173 : vecsmall_is1to1spec(GEN v, long n, GEN w)
278 : {
279 54173 : pari_sp ltop=avma;
280 54173 : long nx=n>>1, ny=n-nx;
281 : long m, ix, iy;
282 : GEN x, y;
283 54173 : if (n<=2)
284 : {
285 32592 : if (n==1)
286 11158 : w[0]=v[0];
287 21434 : else if (n==2)
288 : {
289 21434 : long v0=v[0], v1=v[1];
290 21434 : if (v0==v1) return 0;
291 21406 : else if (v0<v1) { w[0]=v0; w[1]=v1; }
292 4605 : else { w[0]=v1; w[1]=v0; }
293 : }
294 32564 : return 1;
295 : }
296 21581 : x = new_chunk(nx);
297 21581 : if (!vecsmall_is1to1spec(v,nx,x)) return 0;
298 21483 : y = new_chunk(ny);
299 21483 : if (!vecsmall_is1to1spec(v+nx,ny,y)) return 0;
300 84073 : for (m=0, ix=0, iy=0; ix<nx && iy<ny; )
301 62688 : if (x[ix]==y[iy]) return 0;
302 62639 : else if (x[ix]<y[iy])
303 38003 : w[m++]=x[ix++];
304 : else
305 24636 : w[m++]=y[iy++];
306 23499 : for(;ix<nx;) w[m++]=x[ix++];
307 52959 : for(;iy<ny;) w[m++]=y[iy++];
308 21385 : set_avma(ltop);
309 21385 : return 1;
310 : }
311 :
312 : int
313 11207 : vecsmall_is1to1(GEN V)
314 : {
315 11207 : pari_sp av = avma;
316 : long l;
317 11207 : GEN W = cgetg_copy(V, &l);
318 11207 : if (l <= 2) return 1;
319 11109 : return gc_bool(av, vecsmall_is1to1spec(V+1,l,W+1));
320 : }
321 :
322 : /*************************************************************************/
323 : /** **/
324 : /** Routines for handling vectors of VECSMALL **/
325 : /** **/
326 : /*************************************************************************/
327 :
328 : GEN
329 7 : vecvecsmall_sort(GEN x)
330 7 : { return gen_sort(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
331 : GEN
332 361041 : vecvecsmall_sort_shallow(GEN x)
333 361041 : { return gen_sort_shallow(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
334 :
335 : void
336 119 : vecvecsmall_sort_inplace(GEN x, GEN *perm)
337 119 : { gen_sort_inplace(x, (void*)&vecsmall_lexcmp, cmp_nodata, perm); }
338 :
339 : GEN
340 462 : vecvecsmall_sort_uniq(GEN x)
341 462 : { return gen_sort_uniq(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
342 :
343 : GEN
344 441 : vecvecsmall_indexsort(GEN x)
345 441 : { return gen_indexsort(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
346 :
347 : long
348 22934359 : vecvecsmall_search(GEN x, GEN y)
349 22934359 : { return gen_search(x,y,(void*)vecsmall_prefixcmp, cmp_nodata); }
350 :
351 : /* assume x non empty */
352 : long
353 0 : vecvecsmall_max(GEN x)
354 : {
355 0 : long i, l = lg(x), m = vecsmall_max(gel(x,1));
356 0 : for (i = 2; i < l; i++)
357 : {
358 0 : long t = vecsmall_max(gel(x,i));
359 0 : if (t > m) m = t;
360 : }
361 0 : return m;
362 : }
363 :
364 : /*************************************************************************/
365 : /** **/
366 : /** Routines for handling permutations **/
367 : /** **/
368 : /*************************************************************************/
369 :
370 : /* Permutations may be given by
371 : * perm (VECSMALL): a bijection from 1...n to 1...n i-->perm[i]
372 : * cyc (VEC of VECSMALL): a product of disjoint cycles. */
373 :
374 : /* Multiply (compose) two permutations, putting the result in the second one. */
375 : static void
376 21 : perm_mul_inplace2(GEN s, GEN t)
377 : {
378 21 : long i, l = lg(s);
379 525 : for (i = 1; i < l; i++) t[i] = s[t[i]];
380 21 : }
381 :
382 : GEN
383 0 : vecperm_extendschreier(GEN C, GEN v, long n)
384 : {
385 0 : pari_sp av = avma;
386 0 : long mj, lv = lg(v), m = 1, mtested = 1;
387 0 : GEN bit = const_vecsmall(n, 0);
388 0 : GEN cy = cgetg(n+1, t_VECSMALL);
389 0 : GEN sh = const_vec(n, gen_0);
390 0 : for(mj=1; mj<=n; mj++)
391 : {
392 0 : if (isintzero(gel(C,mj))) continue;
393 0 : gel(sh,mj) = gcopy(gel(C,mj));
394 0 : if (bit[mj]) continue;
395 0 : cy[m++] = mj;
396 0 : bit[mj] = 1;
397 : for(;;)
398 0 : {
399 0 : long o, mold = m;
400 0 : for (o = 1; o < lv; o++)
401 : {
402 0 : GEN vo = gel(v,o);
403 : long p;
404 0 : for (p = mtested; p < mold; p++) /* m increases! */
405 : {
406 0 : long j = vo[ cy[p] ];
407 0 : if (!bit[j])
408 : {
409 0 : gel(sh,j) = perm_mul(vo, gel(sh, cy[p]));
410 0 : cy[m++] = j;
411 : }
412 0 : bit[j] = 1;
413 : }
414 : }
415 0 : mtested = mold;
416 0 : if (m == mold) break;
417 : }
418 : }
419 0 : return gerepileupto(av, sh);
420 : }
421 :
422 : /* Orbits of the subgroup generated by v on {1,..,n} */
423 : static GEN
424 1420658 : vecperm_orbits_i(GEN v, long n)
425 : {
426 1420658 : long mj = 1, lv = lg(v), k, l;
427 1420658 : GEN cycle = cgetg(n+1, t_VEC), bit = const_vecsmall(n, 0);
428 8715801 : for (k = 1, l = 1; k <= n;)
429 : {
430 7295063 : pari_sp ltop = avma;
431 7295063 : long m = 1;
432 7295063 : GEN cy = cgetg(n+1, t_VECSMALL);
433 8643749 : for ( ; bit[mj]; mj++) /*empty*/;
434 7294989 : k++; cy[m++] = mj;
435 7294989 : bit[mj++] = 1;
436 : for(;;)
437 2350641 : {
438 9645630 : long o, mold = m;
439 19305533 : for (o = 1; o < lv; o++)
440 : {
441 9659903 : GEN vo = gel(v,o);
442 : long p;
443 30311861 : for (p = 1; p < m; p++) /* m increases! */
444 : {
445 20651958 : long j = vo[ cy[p] ];
446 20651958 : if (!bit[j]) cy[m++] = j;
447 20651958 : bit[j] = 1;
448 : }
449 : }
450 9645630 : if (m == mold) break;
451 2350641 : k += m - mold;
452 : }
453 7294989 : setlg(cy, m);
454 7295058 : gel(cycle,l++) = gerepileuptoleaf(ltop, cy);
455 : }
456 1420738 : setlg(cycle, l); return cycle;
457 : }
458 : /* memory clean version */
459 : GEN
460 4746 : vecperm_orbits(GEN v, long n)
461 : {
462 4746 : pari_sp av = avma;
463 4746 : return gerepilecopy(av, vecperm_orbits_i(v, n));
464 : }
465 :
466 : static int
467 2667 : isperm(GEN v)
468 : {
469 2667 : pari_sp av = avma;
470 2667 : long i, n = lg(v)-1;
471 : GEN w;
472 2667 : if (typ(v) != t_VECSMALL) return 0;
473 2667 : w = zero_zv(n);
474 26411 : for (i=1; i<=n; i++)
475 : {
476 23779 : long d = v[i];
477 23779 : if (d < 1 || d > n || w[d]) return gc_bool(av,0);
478 23744 : w[d] = 1;
479 : }
480 2632 : return gc_bool(av,1);
481 : }
482 :
483 : /* Compute the cyclic decomposition of a permutation */
484 : GEN
485 13314 : perm_cycles(GEN v)
486 : {
487 13314 : pari_sp av = avma;
488 13314 : return gerepilecopy(av, vecperm_orbits_i(mkvec(v), lg(v)-1));
489 : }
490 :
491 : GEN
492 259 : permcycles(GEN v)
493 : {
494 259 : if (!isperm(v)) pari_err_TYPE("permcycles",v);
495 252 : return perm_cycles(v);
496 : }
497 :
498 : /* Output the order of p */
499 : ulong
500 433449 : perm_orderu(GEN v)
501 : {
502 433449 : pari_sp av = avma;
503 433449 : GEN c = vecperm_orbits_i(mkvec(v), lg(v)-1);
504 : long i, d;
505 3088445 : for(i=1, d=1; i<lg(c); i++) d = ulcm(d, lg(gel(c,i))-1);
506 433478 : return gc_ulong(av,d);
507 : }
508 :
509 : static GEN
510 2002 : _domul(void *data, GEN x, GEN y)
511 : {
512 2002 : GEN (*mul)(GEN,GEN)=(GEN (*)(GEN,GEN)) data;
513 2002 : return mul(x,y);
514 : }
515 :
516 : /* Output the order of p */
517 : GEN
518 427 : perm_order(GEN v)
519 : {
520 427 : pari_sp av = avma;
521 427 : GEN c = vecperm_orbits_i(mkvec(v), lg(v)-1);
522 427 : long i, l = lg(c);
523 427 : GEN V = cgetg(l, t_VEC);
524 2856 : for (i = 1; i < l; i++)
525 2429 : gel(V,i) = utoi(lg(gel(c,i))-1);
526 427 : return gerepileuptoint(av, gen_product(V, (void *)lcmii, _domul));
527 : }
528 :
529 : GEN
530 434 : permorder(GEN v)
531 : {
532 434 : if (!isperm(v)) pari_err_TYPE("permorder",v);
533 427 : return perm_order(v);
534 : }
535 :
536 : /* sign of a permutation */
537 : long
538 968724 : perm_sign(GEN v)
539 : {
540 968724 : pari_sp av = avma;
541 968724 : GEN c = vecperm_orbits_i(mkvec(v), lg(v)-1);
542 968766 : long i, l = lg(c), s = 1;
543 5550387 : for (i = 1; i < l; i++)
544 4581627 : if (odd(lg(gel(c, i)))) s = -s;
545 968760 : return gc_long(av,s);
546 : }
547 :
548 : long
549 273 : permsign(GEN v)
550 : {
551 273 : if (!isperm(v)) pari_err_TYPE("permsign",v);
552 259 : return perm_sign(v);
553 : }
554 :
555 : GEN
556 5915 : Z_to_perm(long n, GEN x)
557 : {
558 : pari_sp av;
559 : ulong i, r;
560 5915 : GEN v = cgetg(n+1, t_VECSMALL);
561 5915 : if (n==0) return v;
562 5908 : uel(v,n) = 1; av = avma;
563 5908 : if (signe(x) <= 0) x = modii(x, mpfact(n));
564 27146 : for (r=n-1; r>=1; r--)
565 : {
566 : ulong a;
567 21238 : x = absdiviu_rem(x, n+1-r, &a);
568 71687 : for (i=r+1; i<=(ulong)n; i++)
569 50449 : if (uel(v,i) > a) uel(v,i)++;
570 21238 : uel(v,r) = a+1;
571 : }
572 5908 : set_avma(av); return v;
573 : }
574 : GEN
575 5915 : numtoperm(long n, GEN x)
576 : {
577 5915 : if (n < 0) pari_err_DOMAIN("numtoperm", "n", "<", gen_0, stoi(n));
578 5915 : if (typ(x) != t_INT) pari_err_TYPE("numtoperm",x);
579 5915 : return Z_to_perm(n, x);
580 : }
581 :
582 : /* destroys v */
583 : static GEN
584 1701 : perm_to_Z_inplace(GEN v)
585 : {
586 1701 : long l = lg(v), i, r;
587 1701 : GEN x = gen_0;
588 1701 : if (!isperm(v)) return NULL;
589 10143 : for (i = 1; i < l; i++)
590 : {
591 8449 : long vi = v[i];
592 8449 : if (vi <= 0) return NULL;
593 8449 : x = i==1 ? utoi(vi-1): addiu(muliu(x,l-i), vi-1);
594 25396 : for (r = i+1; r < l; r++)
595 16947 : if (v[r] > vi) v[r]--;
596 : }
597 1694 : return x;
598 : }
599 : GEN
600 1680 : perm_to_Z(GEN v)
601 : {
602 1680 : pari_sp av = avma;
603 1680 : GEN x = perm_to_Z_inplace(leafcopy(v));
604 1680 : if (!x) pari_err_TYPE("permtonum",v);
605 1680 : return gerepileuptoint(av, x);
606 : }
607 : GEN
608 1708 : permtonum(GEN p)
609 : {
610 1708 : pari_sp av = avma;
611 : GEN v, x;
612 1708 : switch(typ(p))
613 : {
614 1680 : case t_VECSMALL: return perm_to_Z(p);
615 21 : case t_VEC: case t_COL:
616 21 : if (RgV_is_ZV(p)) { v = ZV_to_zv(p); break; }
617 7 : default: pari_err_TYPE("permtonum",p);
618 : return NULL;/*LCOV_EXCL_LINE*/
619 : }
620 21 : x = perm_to_Z_inplace(v);
621 21 : if (!x) pari_err_TYPE("permtonum",p);
622 14 : return gerepileuptoint(av, x);
623 : }
624 :
625 : GEN
626 7301 : cyc_pow(GEN cyc, long exp)
627 : {
628 : long i, j, k, l, r;
629 : GEN c;
630 22169 : for (r = j = 1; j < lg(cyc); j++)
631 : {
632 14868 : long n = lg(gel(cyc,j)) - 1;
633 14868 : r += cgcd(n, exp);
634 : }
635 7301 : c = cgetg(r, t_VEC);
636 22169 : for (r = j = 1; j < lg(cyc); j++)
637 : {
638 14868 : GEN v = gel(cyc,j);
639 14868 : long n = lg(v) - 1, e = umodsu(exp,n), g = (long)ugcd(n, e), m = n / g;
640 31654 : for (i = 0; i < g; i++)
641 : {
642 16786 : GEN p = cgetg(m+1, t_VECSMALL);
643 16786 : gel(c,r++) = p;
644 54719 : for (k = 1, l = i; k <= m; k++)
645 : {
646 37933 : p[k] = v[l+1];
647 37933 : l += e; if (l >= n) l -= n;
648 : }
649 : }
650 : }
651 7301 : return c;
652 : }
653 :
654 : /* Compute the power of a permutation given by product of cycles
655 : * Ouput a perm, not a cyc */
656 : GEN
657 0 : cyc_pow_perm(GEN cyc, long exp)
658 : {
659 : long e, j, k, l, n;
660 : GEN p;
661 0 : for (n = 0, j = 1; j < lg(cyc); j++) n += lg(gel(cyc,j))-1;
662 0 : p = cgetg(n + 1, t_VECSMALL);
663 0 : for (j = 1; j < lg(cyc); j++)
664 : {
665 0 : GEN v = gel(cyc,j);
666 0 : n = lg(v) - 1; e = umodsu(exp, n);
667 0 : for (k = 1, l = e; k <= n; k++)
668 : {
669 0 : p[v[k]] = v[l+1];
670 0 : if (++l == n) l = 0;
671 : }
672 : }
673 0 : return p;
674 : }
675 :
676 : GEN
677 49 : perm_pow(GEN perm, GEN exp)
678 : {
679 49 : long i, r = lg(perm)-1;
680 49 : GEN p = zero_zv(r);
681 49 : pari_sp av = avma;
682 49 : GEN v = cgetg(r+1, t_VECSMALL);
683 168 : for (i=1; i<=r; i++)
684 : {
685 : long e, n, k, l;
686 119 : if (p[i]) continue;
687 49 : v[1] = i;
688 119 : for (n=1, k=perm[i]; k!=i; k=perm[k], n++) v[n+1] = k;
689 49 : e = umodiu(exp, n);
690 168 : for (k = 1, l = e; k <= n; k++)
691 : {
692 119 : p[v[k]] = v[l+1];
693 119 : if (++l == n) l = 0;
694 : }
695 : }
696 49 : set_avma(av); return p;
697 : }
698 :
699 : GEN
700 18585 : perm_powu(GEN perm, ulong exp)
701 : {
702 18585 : ulong i, r = lg(perm)-1;
703 18585 : GEN p = zero_zv(r);
704 18585 : pari_sp av = avma;
705 18585 : GEN v = cgetg(r+1, t_VECSMALL);
706 245994 : for (i=1; i<=r; i++)
707 : {
708 : ulong e, n, k, l;
709 227409 : if (p[i]) continue;
710 84525 : v[1] = i;
711 227409 : for (n=1, k=perm[i]; k!=i; k=perm[k], n++) v[n+1] = k;
712 84525 : e = exp % n;
713 311934 : for (k = 1, l = e; k <= n; k++)
714 : {
715 227409 : p[v[k]] = v[l+1];
716 227409 : if (++l == n) l = 0;
717 : }
718 : }
719 18585 : set_avma(av); return p;
720 : }
721 :
722 : GEN
723 21 : perm_to_GAP(GEN p)
724 : {
725 21 : pari_sp ltop=avma;
726 : GEN gap;
727 : GEN x;
728 : long i;
729 21 : long nb, c=0;
730 : char *s;
731 : long sz;
732 21 : long lp=lg(p)-1;
733 21 : if (typ(p) != t_VECSMALL) pari_err_TYPE("perm_to_GAP",p);
734 21 : x = perm_cycles(p);
735 21 : sz = (long) ((bfffo(lp)+1) * LOG10_2 + 1);
736 : /*Dry run*/
737 133 : for (i = 1, nb = 1; i < lg(x); ++i)
738 : {
739 112 : GEN z = gel(x,i);
740 112 : long lz = lg(z)-1;
741 112 : nb += 1+lz*(sz+2);
742 : }
743 21 : nb++;
744 : /*Real run*/
745 21 : gap = cgetg(nchar2nlong(nb) + 1, t_STR);
746 21 : s = GSTR(gap);
747 133 : for (i = 1; i < lg(x); ++i)
748 : {
749 : long j;
750 112 : GEN z = gel(x,i);
751 112 : if (lg(z) > 2)
752 : {
753 112 : s[c++] = '(';
754 364 : for (j = 1; j < lg(z); ++j)
755 : {
756 252 : if (j > 1)
757 : {
758 140 : s[c++] = ','; s[c++] = ' ';
759 : }
760 252 : sprintf(s+c,"%ld",z[j]);
761 567 : while(s[c++]) /* empty */;
762 252 : c--;
763 : }
764 112 : s[c++] = ')';
765 : }
766 : }
767 21 : if (!c) { s[c++]='('; s[c++]=')'; }
768 21 : s[c] = '\0';
769 21 : return gerepileupto(ltop,gap);
770 : }
771 :
772 : int
773 572495 : perm_commute(GEN s, GEN t)
774 : {
775 572495 : long i, l = lg(t);
776 40373487 : for (i = 1; i < l; i++)
777 39820382 : if (t[ s[i] ] != s[ t[i] ]) return 0;
778 553105 : return 1;
779 : }
780 :
781 : /*************************************************************************/
782 : /** **/
783 : /** Routines for handling groups **/
784 : /** **/
785 : /*************************************************************************/
786 : /* A Group is a t_VEC [gen,orders]
787 : * gen (vecvecsmall): list of generators given by permutations
788 : * orders (vecsmall): relatives orders of generators. */
789 936900 : INLINE GEN grp_get_gen(GEN G) { return gel(G,1); }
790 1586121 : INLINE GEN grp_get_ord(GEN G) { return gel(G,2); }
791 :
792 : /* A Quotient Group is a t_VEC [gen,coset]
793 : * gen (vecvecsmall): coset generators
794 : * coset (vecsmall): gen[coset[p[1]]] generate the p-coset.
795 : */
796 140950 : INLINE GEN quo_get_gen(GEN C) { return gel(C,1); }
797 29708 : INLINE GEN quo_get_coset(GEN C) { return gel(C,2); }
798 :
799 : static GEN
800 52346 : trivialsubgroups(void)
801 52346 : { GEN L = cgetg(2, t_VEC); gel(L,1) = trivialgroup(); return L; }
802 :
803 : /* Compute the order of p modulo the group given by a set */
804 : long
805 220108 : perm_relorder(GEN p, GEN set)
806 : {
807 220108 : pari_sp ltop = avma;
808 220108 : long n = 1, q = p[1];
809 653863 : while (!F2v_coeff(set,q)) { q = p[q]; n++; }
810 220108 : return gc_long(ltop,n);
811 : }
812 :
813 : GEN
814 13048 : perm_generate(GEN S, GEN H, long o)
815 : {
816 13048 : long i, n = lg(H)-1;
817 13048 : GEN L = cgetg(n*o + 1, t_VEC);
818 45829 : for(i=1; i<=n; i++) gel(L,i) = vecsmall_copy(gel(H,i));
819 50617 : for( ; i <= n*o; i++) gel(L,i) = perm_mul(gel(L,i-n), S);
820 13048 : return L;
821 : }
822 :
823 : /*Return the order (cardinality) of a group */
824 : long
825 705584 : group_order(GEN G)
826 : {
827 705584 : return zv_prod(grp_get_ord(G));
828 : }
829 :
830 : /* G being a subgroup of S_n, output n */
831 : long
832 26628 : group_domain(GEN G)
833 : {
834 26628 : GEN gen = grp_get_gen(G);
835 26628 : if (lg(gen) < 2) pari_err_DOMAIN("group_domain", "#G", "=", gen_1,G);
836 26628 : return lg(gel(gen,1)) - 1;
837 : }
838 :
839 : /*Left coset of g mod G: gG*/
840 : GEN
841 301621 : group_leftcoset(GEN G, GEN g)
842 : {
843 301621 : GEN gen = grp_get_gen(G), ord = grp_get_ord(G);
844 301619 : GEN res = cgetg(group_order(G)+1, t_VEC);
845 : long i, j, k;
846 301619 : gel(res,1) = vecsmall_copy(g);
847 301625 : k = 1;
848 554644 : for (i = 1; i < lg(gen); i++)
849 : {
850 253022 : long c = k * (ord[i] - 1);
851 695857 : for (j = 1; j <= c; j++) gel(res,++k) = perm_mul(gel(res,j), gel(gen,i));
852 : }
853 301622 : return res;
854 : }
855 : /*Right coset of g mod G: Gg*/
856 : GEN
857 179842 : group_rightcoset(GEN G, GEN g)
858 : {
859 179842 : GEN gen = grp_get_gen(G), ord = grp_get_ord(G);
860 179842 : GEN res = cgetg(group_order(G)+1, t_VEC);
861 : long i, j, k;
862 179841 : gel(res,1) = vecsmall_copy(g);
863 179845 : k = 1;
864 310750 : for (i = 1; i < lg(gen); i++)
865 : {
866 130907 : long c = k * (ord[i] - 1);
867 411899 : for (j = 1; j <= c; j++) gel(res,++k) = perm_mul(gel(gen,i), gel(res,j));
868 : }
869 179843 : return res;
870 : }
871 : /*Elements of a group from the generators, cf group_leftcoset*/
872 : GEN
873 140867 : group_elts(GEN G, long n)
874 : {
875 140867 : GEN gen = grp_get_gen(G), ord = grp_get_ord(G);
876 140867 : GEN res = cgetg(group_order(G)+1, t_VEC);
877 : long i, j, k;
878 140867 : gel(res,1) = identity_perm(n);
879 140867 : k = 1;
880 285528 : for (i = 1; i < lg(gen); i++)
881 : {
882 144661 : long c = k * (ord[i] - 1);
883 : /* j = 1, use res[1] = identity */
884 144661 : gel(res,++k) = vecsmall_copy(gel(gen,i));
885 384425 : for (j = 2; j <= c; j++) gel(res,++k) = perm_mul(gel(res,j), gel(gen,i));
886 : }
887 140867 : return res;
888 : }
889 :
890 : GEN
891 14448 : groupelts_conj_set(GEN elts, GEN p)
892 : {
893 14448 : long i, j, l = lg(elts), n = lg(p)-1;
894 14448 : GEN res = zero_F2v(n);
895 241465 : for(j = 1; j < n; j++)
896 241465 : if (p[j]==1) break;
897 101136 : for(i = 1; i < l; i++)
898 86688 : F2v_set(res, p[mael(elts,i,j)]);
899 14448 : return res;
900 : }
901 :
902 : GEN
903 28014 : groupelts_set(GEN elts, long n)
904 : {
905 28014 : GEN res = zero_F2v(n);
906 28014 : long i, l = lg(elts);
907 137410 : for(i=1; i<l; i++)
908 109396 : F2v_set(res,mael(elts,i,1));
909 28014 : return res;
910 : }
911 :
912 : /*Elements of a group from the generators, returned as a set (bitmap)*/
913 : GEN
914 90671 : group_set(GEN G, long n)
915 : {
916 90671 : GEN res = zero_F2v(n);
917 90671 : pari_sp av = avma;
918 90671 : GEN elts = group_elts(G, n);
919 90671 : long i, l = lg(elts);
920 284214 : for(i=1; i<l; i++)
921 193543 : F2v_set(res,mael(elts,i,1));
922 90671 : set_avma(av);
923 90671 : return res;
924 : }
925 :
926 : static int
927 17353 : sgcmp(GEN a, GEN b) { return vecsmall_lexcmp(gel(a,1),gel(b,1)); }
928 :
929 : GEN
930 497 : subgroups_tableset(GEN S, long n)
931 : {
932 497 : long i, l = lg(S);
933 497 : GEN v = cgetg(l, t_VEC);
934 5411 : for(i=1; i<l; i++)
935 4914 : gel(v,i) = mkvec2(group_set(gel(S,i), n), mkvecsmall(i));
936 497 : gen_sort_inplace(v,(void*)sgcmp,cmp_nodata, NULL);
937 497 : return v;
938 : }
939 :
940 : long
941 2002 : tableset_find_index(GEN tbl, GEN set)
942 : {
943 2002 : long i = tablesearch(tbl,mkvec2(set,mkvecsmall(0)),sgcmp);
944 2002 : if (!i) return 0;
945 2002 : return mael3(tbl,i,2,1);
946 : }
947 :
948 : GEN
949 52367 : trivialgroup(void) { retmkvec2(cgetg(1,t_VEC), cgetg(1,t_VECSMALL)); }
950 :
951 : /*Cyclic group generated by g of order s*/
952 : GEN
953 26670 : cyclicgroup(GEN g, long s)
954 26670 : { retmkvec2(mkvec( vecsmall_copy(g) ), mkvecsmall(s)); }
955 :
956 : /*Return the group generated by g1,g2 of relative orders s1,s2*/
957 : GEN
958 1085 : dicyclicgroup(GEN g1, GEN g2, long s1, long s2)
959 1085 : { retmkvec2( mkvec2(vecsmall_copy(g1), vecsmall_copy(g2)),
960 : mkvecsmall2(s1, s2) ); }
961 :
962 : /* return the quotient map G --> G/H */
963 : /*The ouput is [gen,hash]*/
964 : /* gen (vecvecsmall): coset generators
965 : * coset (vecsmall): vecsmall of coset number) */
966 : GEN
967 11522 : groupelts_quotient(GEN elt, GEN H)
968 : {
969 11522 : pari_sp ltop = avma;
970 : GEN p2, p3;
971 11522 : long i, j, a = 1;
972 11522 : long n = lg(gel(elt,1))-1, o = group_order(H);
973 : GEN el;
974 11522 : long le = lg(elt)-1;
975 11522 : GEN used = zero_F2v(le+1);
976 11522 : long l = le/o;
977 11522 : p2 = cgetg(l+1, t_VEC);
978 11522 : p3 = zero_zv(n);
979 11522 : el = zero_zv(n);
980 148722 : for (i = 1; i<=le; i++)
981 137200 : el[mael(elt,i,1)]=i;
982 68075 : for (i = 1; i <= l; ++i)
983 : {
984 : GEN V;
985 150003 : while(F2v_coeff(used,a)) a++;
986 56560 : V = group_leftcoset(H,gel(elt,a));
987 56560 : gel(p2,i) = gel(V,1);
988 193655 : for(j=1;j<lg(V);j++)
989 : {
990 137102 : long b = el[mael(V,j,1)];
991 137102 : if (b==0) pari_err_IMPL("group_quotient for a non-WSS group");
992 137095 : F2v_set(used,b);
993 : }
994 193641 : for (j = 1; j <= o; j++)
995 137088 : p3[mael(V, j, 1)] = i;
996 : }
997 11515 : return gerepilecopy(ltop,mkvec2(p2,p3));
998 : }
999 :
1000 : GEN
1001 10157 : group_quotient(GEN G, GEN H)
1002 : {
1003 10157 : return groupelts_quotient(group_elts(G, group_domain(G)), H);
1004 : }
1005 :
1006 : /*Compute the image of a permutation by a quotient map.*/
1007 : GEN
1008 29708 : quotient_perm(GEN C, GEN p)
1009 : {
1010 29708 : GEN gen = quo_get_gen(C);
1011 29708 : GEN coset = quo_get_coset(C);
1012 29708 : long j, l = lg(gen);
1013 29708 : GEN p3 = cgetg(l, t_VECSMALL);
1014 282163 : for (j = 1; j < l; ++j)
1015 : {
1016 252455 : p3[j] = coset[p[mael(gen,j,1)]];
1017 252455 : if (p3[j]==0) pari_err_IMPL("quotient_perm for a non-WSS group");
1018 : }
1019 29708 : return p3;
1020 : }
1021 :
1022 : /* H is a subgroup of G, C is the quotient map G --> G/H
1023 : *
1024 : * Lift a subgroup S of G/H to a subgroup of G containing H */
1025 : GEN
1026 50547 : quotient_subgroup_lift(GEN C, GEN H, GEN S)
1027 : {
1028 50547 : GEN genH = grp_get_gen(H);
1029 50547 : GEN genS = grp_get_gen(S);
1030 50547 : GEN genC = quo_get_gen(C);
1031 50547 : long l1 = lg(genH)-1;
1032 50547 : long l2 = lg(genS)-1, j;
1033 50547 : GEN p1 = cgetg(3, t_VEC), L = cgetg(l1+l2+1, t_VEC);
1034 101262 : for (j = 1; j <= l1; ++j) gel(L,j) = gel(genH,j);
1035 117719 : for (j = 1; j <= l2; ++j) gel(L,l1+j) = gel(genC, mael(genS,j,1));
1036 50547 : gel(p1,1) = L;
1037 50547 : gel(p1,2) = vecsmall_concat(grp_get_ord(H), grp_get_ord(S));
1038 50547 : return p1;
1039 : }
1040 :
1041 : /* Let G a group and C a quotient map G --> G/H
1042 : * Assume H is normal, return the group G/H */
1043 : GEN
1044 10150 : quotient_group(GEN C, GEN G)
1045 : {
1046 10150 : pari_sp ltop = avma;
1047 : GEN Qgen, Qord, Qelt, Qset, Q;
1048 10150 : GEN Cgen = quo_get_gen(C);
1049 10150 : GEN Ggen = grp_get_gen(G);
1050 10150 : long i,j, n = lg(Cgen)-1, l = lg(Ggen);
1051 10150 : Qord = cgetg(l, t_VECSMALL);
1052 10150 : Qgen = cgetg(l, t_VEC);
1053 10150 : Qelt = mkvec(identity_perm(n));
1054 10150 : Qset = groupelts_set(Qelt, n);
1055 31052 : for (i = 1, j = 1; i < l; ++i)
1056 : {
1057 20902 : GEN g = quotient_perm(C, gel(Ggen,i));
1058 20902 : long o = perm_relorder(g, Qset);
1059 20902 : gel(Qgen,j) = g;
1060 20902 : Qord[j] = o;
1061 20902 : if (o != 1)
1062 : {
1063 13048 : Qelt = perm_generate(g, Qelt, o);
1064 13048 : Qset = groupelts_set(Qelt, n);
1065 13048 : j++;
1066 : }
1067 : }
1068 10150 : setlg(Qgen,j);
1069 10150 : setlg(Qord,j); Q = mkvec2(Qgen, Qord);
1070 10150 : return gerepilecopy(ltop,Q);
1071 : }
1072 :
1073 : GEN
1074 1365 : quotient_groupelts(GEN C)
1075 : {
1076 1365 : GEN G = quo_get_gen(C);
1077 1365 : long i, l = lg(G);
1078 1365 : GEN Q = cgetg(l, t_VEC);
1079 10171 : for (i = 1; i < l; ++i)
1080 8806 : gel(Q,i) = quotient_perm(C, gel(G,i));
1081 1365 : return Q;
1082 : }
1083 :
1084 : /* Return 1 if g normalizes N, 0 otherwise */
1085 : long
1086 179843 : group_perm_normalize(GEN N, GEN g)
1087 : {
1088 179843 : pari_sp ltop = avma;
1089 179843 : long r = gequal(vecvecsmall_sort_shallow(group_leftcoset(N, g)),
1090 : vecvecsmall_sort_shallow(group_rightcoset(N, g)));
1091 179843 : return gc_long(ltop, r);
1092 : }
1093 :
1094 : /* L is a list of subgroups, C is a coset and r a relative order.*/
1095 : static GEN
1096 65219 : liftlistsubgroups(GEN L, GEN C, long r)
1097 : {
1098 65219 : pari_sp ltop = avma;
1099 65219 : long c = lg(C)-1, l = lg(L)-1, n = lg(gel(C,1))-1, i, k;
1100 : GEN R;
1101 65219 : if (!l) return cgetg(1,t_VEC);
1102 58527 : R = cgetg(l*c+1, t_VEC);
1103 142968 : for (i = 1, k = 1; i <= l; ++i)
1104 : {
1105 84441 : GEN S = gel(L,i), Selt = group_set(S,n);
1106 84441 : GEN gen = grp_get_gen(S);
1107 84441 : GEN ord = grp_get_ord(S);
1108 : long j;
1109 279202 : for (j = 1; j <= c; ++j)
1110 : {
1111 194761 : GEN p = gel(C,j);
1112 194761 : if (perm_relorder(p, Selt) == r && group_perm_normalize(S, p))
1113 108731 : gel(R,k++) = mkvec2(vec_append(gen, p),
1114 : vecsmall_append(ord, r));
1115 : }
1116 : }
1117 58527 : setlg(R, k);
1118 58527 : return gerepilecopy(ltop, R);
1119 : }
1120 :
1121 : /* H is a normal subgroup, C is the quotient map G -->G/H,
1122 : * S is a subgroup of G/H, and G is embedded in Sym(l)
1123 : * Return all the subgroups K of G such that
1124 : * S= K mod H and K inter H={1} */
1125 : static GEN
1126 49182 : liftsubgroup(GEN C, GEN H, GEN S)
1127 : {
1128 49182 : pari_sp ltop = avma;
1129 49182 : GEN V = trivialsubgroups();
1130 49182 : GEN Sgen = grp_get_gen(S);
1131 49182 : GEN Sord = grp_get_ord(S);
1132 49182 : GEN Cgen = quo_get_gen(C);
1133 49182 : long n = lg(Sgen), i;
1134 114401 : for (i = 1; i < n; ++i)
1135 : { /*loop over generators of S*/
1136 65219 : GEN W = group_leftcoset(H, gel(Cgen, mael(Sgen, i, 1)));
1137 65219 : V = liftlistsubgroups(V, W, Sord[i]);
1138 : }
1139 49182 : return gerepilecopy(ltop,V);
1140 : }
1141 :
1142 : /* 1:A4, 2:S4, 3:F36, 0: other */
1143 : long
1144 9982 : group_isA4S4(GEN G)
1145 : {
1146 9982 : GEN elt = grp_get_gen(G);
1147 9982 : GEN ord = grp_get_ord(G);
1148 9982 : long n = lg(ord);
1149 9982 : if (n != 4 && n != 5) return 0;
1150 2219 : if (n==4 && ord[1]==3 && ord[2]==3 && ord[3]==4)
1151 : {
1152 : long i;
1153 7 : GEN p = gel(elt,1), q = gel(elt,2), r = gel(elt,3);
1154 259 : for(i=1; i<=36; i++)
1155 252 : if (p[r[i]]!=r[q[i]]) return 0;
1156 7 : return 3;
1157 : }
1158 2212 : if (ord[1]!=2 || ord[2]!=2 || ord[3]!=3) return 0;
1159 42 : if (perm_commute(gel(elt,1),gel(elt,3))) return 0;
1160 42 : if (n==4) return 1;
1161 21 : if (ord[4]!=2) return 0;
1162 21 : if (perm_commute(gel(elt,3),gel(elt,4))) return 0;
1163 21 : return 2;
1164 : }
1165 : /* compute all the subgroups of a group G */
1166 : GEN
1167 13146 : group_subgroups(GEN G)
1168 : {
1169 13146 : pari_sp ltop = avma;
1170 : GEN p1, H, C, Q, M, sg1, sg2, sg3;
1171 13146 : GEN gen = grp_get_gen(G);
1172 13146 : GEN ord = grp_get_ord(G);
1173 13146 : long lM, i, j, n = lg(gen);
1174 : long t;
1175 13146 : if (n == 1) return trivialsubgroups();
1176 9982 : t = group_isA4S4(G);
1177 9982 : if (t == 3)
1178 : {
1179 7 : GEN H = mkvec2(mkvec3(gel(gen,1), gel(gen,2), perm_sqr(gel(gen,3))),
1180 : mkvecsmall3(3, 3, 2));
1181 7 : GEN S = group_subgroups(H);
1182 7 : GEN V = cgetg(11,t_VEC);
1183 7 : gel(V,1) = cyclicgroup(gel(gen,3),4);
1184 63 : for (i=2; i<10; i++)
1185 56 : gel(V,i) = cyclicgroup(perm_mul(gmael3(V,i-1,1,1),gel(gen,i%3==1 ? 2:1)),4);
1186 7 : gel(V,10) = G;
1187 7 : return gerepilecopy(ltop,shallowconcat(S,V));
1188 : }
1189 9975 : else if (t)
1190 : {
1191 42 : GEN s = gel(gen,1); /*s = (1,2)(3,4) */
1192 42 : GEN t = gel(gen,2); /*t = (1,3)(2,4) */
1193 42 : GEN st = perm_mul(s, t); /*st = (1,4)(2,3) */
1194 42 : H = dicyclicgroup(s, t, 2, 2);
1195 : /* sg3 is the list of subgroups intersecting only partially with H*/
1196 42 : sg3 = cgetg((n==4)?4: 10, t_VEC);
1197 42 : gel(sg3,1) = cyclicgroup(s, 2);
1198 42 : gel(sg3,2) = cyclicgroup(t, 2);
1199 42 : gel(sg3,3) = cyclicgroup(st, 2);
1200 42 : if (n==5)
1201 : {
1202 21 : GEN u = gel(gen,3);
1203 21 : GEN v = gel(gen,4), w, u2;
1204 21 : if (zv_equal(perm_conj(u,s), t)) /*u=(2,3,4)*/
1205 21 : u2 = perm_sqr(u);
1206 : else
1207 : {
1208 0 : u2 = u;
1209 0 : u = perm_sqr(u);
1210 : }
1211 21 : if (perm_orderu(v)==2)
1212 : {
1213 21 : if (!perm_commute(s,v)) /*v=(1,2)*/
1214 : {
1215 0 : v = perm_conj(u,v);
1216 0 : if (!perm_commute(s,v)) v = perm_conj(u,v);
1217 : }
1218 21 : w = perm_mul(v,t); /*w=(1,4,2,3)*/
1219 : }
1220 : else
1221 : {
1222 0 : w = v;
1223 0 : if (!zv_equal(perm_sqr(w), s)) /*w=(1,4,2,3)*/
1224 : {
1225 0 : w = perm_conj(u,w);
1226 0 : if (!zv_equal(perm_sqr(w), s)) w = perm_conj(u,w);
1227 : }
1228 0 : v = perm_mul(w,t); /*v=(1,2)*/
1229 : }
1230 21 : gel(sg3,4) = dicyclicgroup(s,v,2,2);
1231 21 : gel(sg3,5) = dicyclicgroup(t,perm_conj(u,v),2,2);
1232 21 : gel(sg3,6) = dicyclicgroup(st,perm_conj(u2,v),2,2);
1233 21 : gel(sg3,7) = dicyclicgroup(s,w,2,2);
1234 21 : gel(sg3,8) = dicyclicgroup(t,perm_conj(u,w),2,2);
1235 21 : gel(sg3,9) = dicyclicgroup(st,perm_conj(u2,w),2,2);
1236 : }
1237 : }
1238 : else
1239 : {
1240 9933 : ulong osig = mael(factoru(ord[1]), 1, 1);
1241 9933 : GEN sig = perm_powu(gel(gen,1), ord[1]/osig);
1242 9933 : H = cyclicgroup(sig,osig);
1243 9933 : sg3 = NULL;
1244 : }
1245 9975 : C = group_quotient(G,H);
1246 9968 : Q = quotient_group(C,G);
1247 9968 : M = group_subgroups(Q); lM = lg(M);
1248 : /* sg1 is the list of subgroups containing H*/
1249 9961 : sg1 = cgetg(lM, t_VEC);
1250 59143 : for (i = 1; i < lM; ++i) gel(sg1,i) = quotient_subgroup_lift(C,H,gel(M,i));
1251 : /*sg2 is a list of lists of subgroups not intersecting with H*/
1252 9961 : sg2 = cgetg(lM, t_VEC);
1253 : /* Loop over all subgroups of G/H */
1254 59143 : for (j = 1; j < lM; ++j) gel(sg2,j) = liftsubgroup(C, H, gel(M,j));
1255 9961 : p1 = gconcat(sg1, shallowconcat1(sg2));
1256 9961 : if (sg3)
1257 : {
1258 42 : p1 = gconcat(p1, sg3);
1259 42 : if (n==5) /*ensure that the D4 subgroups of S4 are in supersolvable format*/
1260 84 : for(j = 3; j <= 5; j++)
1261 : {
1262 63 : GEN c = gmael(p1,j,1);
1263 63 : if (!perm_commute(gel(c,1),gel(c,3)))
1264 : {
1265 42 : if (perm_commute(gel(c,2),gel(c,3))) { swap(gel(c,1), gel(c,2)); }
1266 : else
1267 21 : perm_mul_inplace2(gel(c,2), gel(c,1));
1268 : }
1269 : }
1270 : }
1271 9961 : return gerepileupto(ltop,p1);
1272 : }
1273 :
1274 : /*return 1 if G is abelian, else 0*/
1275 : long
1276 9240 : group_isabelian(GEN G)
1277 : {
1278 9240 : GEN g = grp_get_gen(G);
1279 9240 : long i, j, n = lg(g);
1280 13160 : for(i=2; i<n; i++)
1281 13034 : for(j=1; j<i; j++)
1282 9114 : if (!perm_commute(gel(g,i), gel(g,j))) return 0;
1283 4298 : return 1;
1284 : }
1285 :
1286 : /*If G is abelian, return its HNF matrix*/
1287 : GEN
1288 385 : group_abelianHNF(GEN G, GEN S)
1289 : {
1290 385 : GEN M, g = grp_get_gen(G), o = grp_get_ord(G);
1291 385 : long i, j, k, n = lg(g);
1292 385 : if (!group_isabelian(G)) return NULL;
1293 315 : if (n==1) return cgetg(1,t_MAT);
1294 301 : if (!S) S = group_elts(G, group_domain(G));
1295 301 : M = cgetg(n,t_MAT);
1296 980 : for(i=1; i<n; i++)
1297 : {
1298 679 : GEN P, C = cgetg(n,t_COL);
1299 679 : pari_sp av = avma;
1300 679 : gel(M,i) = C;
1301 679 : P = perm_inv(perm_powu(gel(g,i), o[i]));
1302 959 : for(j=1; j<lg(S); j++)
1303 959 : if (zv_equal(P, gel(S,j))) break;
1304 679 : set_avma(av);
1305 679 : if (j==lg(S)) pari_err_BUG("galoisisabelian [inconsistent group]");
1306 679 : j--;
1307 1218 : for(k=1; k<i; k++)
1308 : {
1309 539 : long q = j / o[k];
1310 539 : gel(C,k) = stoi(j - q*o[k]);
1311 539 : j = q;
1312 : }
1313 679 : gel(C,k) = stoi(o[i]);
1314 1218 : for (k++; k<n; k++) gel(C,k) = gen_0;
1315 : }
1316 301 : return M;
1317 : }
1318 :
1319 : /*If G is abelian, return its abstract SNF matrix*/
1320 : GEN
1321 336 : group_abelianSNF(GEN G, GEN L)
1322 : {
1323 336 : pari_sp ltop = avma;
1324 336 : GEN H = group_abelianHNF(G,L);
1325 336 : if (!H) return NULL;
1326 266 : return gerepileupto(ltop, smithclean( ZM_snf(H) ));
1327 : }
1328 :
1329 : GEN
1330 469 : abelian_group(GEN v)
1331 : {
1332 469 : long card = zv_prod(v), i, d = 1, l = lg(v);
1333 469 : GEN G = cgetg(3,t_VEC), gen = cgetg(l,t_VEC);
1334 469 : gel(G,1) = gen;
1335 469 : gel(G,2) = vecsmall_copy(v);
1336 980 : for(i=1; i<l; i++)
1337 : {
1338 511 : GEN p = cgetg(card+1, t_VECSMALL);
1339 511 : long o = v[i], u = d*(o-1), j, k, l;
1340 511 : gel(gen, i) = p;
1341 : /* The following loop is over-optimized. Remember that I wrote it for
1342 : * testpermutation. Something has survived... BA */
1343 1162 : for(j=1;j<=card;)
1344 : {
1345 2394 : for(k=1;k<o;k++)
1346 4613 : for(l=1;l<=d; l++,j++) p[j] = j+d;
1347 2121 : for (l=1; l<=d; l++,j++) p[j] = j-u;
1348 : }
1349 511 : d += u;
1350 : }
1351 469 : return G;
1352 : }
1353 :
1354 : static long
1355 13783 : groupelts_subgroup_isnormal(GEN G, GEN H)
1356 : {
1357 13783 : long i, n = lg(G);
1358 61922 : for(i = 1; i < n; i++)
1359 60494 : if (!group_perm_normalize(H, gel(G,i))) return 0;
1360 1428 : return 1;
1361 : }
1362 :
1363 : /*return 1 if H is a normal subgroup of G*/
1364 : long
1365 336 : group_subgroup_isnormal(GEN G, GEN H)
1366 : {
1367 336 : if (lg(grp_get_gen(H)) > 1 && group_domain(G) != group_domain(H))
1368 0 : pari_err_DOMAIN("group_subgroup_isnormal","domain(H)","!=",
1369 : strtoGENstr("domain(G)"), H);
1370 336 : return groupelts_subgroup_isnormal(grp_get_gen(G), H);
1371 : }
1372 :
1373 : static GEN
1374 4816 : group_subgroup_kernel_set(GEN G, GEN H)
1375 : {
1376 : pari_sp av;
1377 4816 : GEN g = grp_get_gen(G);
1378 4816 : long i, n = lg(g);
1379 : GEN S, elts;
1380 4816 : long d = group_domain(G);
1381 4816 : if (lg(grp_get_gen(H)) > 1 && group_domain(G) != group_domain(H))
1382 0 : pari_err_DOMAIN("group_subgroup_isnormal","domain(H)","!=",
1383 : strtoGENstr("domain(G)"), H);
1384 4816 : elts = group_elts(H,d);
1385 4816 : S = groupelts_set(elts, d);
1386 4816 : av = avma;
1387 19264 : for(i=1; i<n; i++)
1388 : {
1389 14448 : F2v_and_inplace(S, groupelts_conj_set(elts,gel(g,i)));
1390 14448 : set_avma(av);
1391 : }
1392 4816 : return S;
1393 : }
1394 :
1395 : int
1396 4816 : group_subgroup_is_faithful(GEN G, GEN H)
1397 : {
1398 4816 : pari_sp av = avma;
1399 4816 : GEN K = group_subgroup_kernel_set(G,H);
1400 4816 : F2v_clear(K,1);
1401 4816 : return gc_long(av, F2v_equal0(K));
1402 : }
1403 :
1404 : long
1405 0 : groupelts_exponent(GEN elts)
1406 : {
1407 0 : long i, n = lg(elts)-1, expo = 1;
1408 0 : for(i=1; i<=n; i++) expo = ulcm(expo, perm_orderu(gel(elts,i)));
1409 0 : return expo;
1410 : }
1411 :
1412 : GEN
1413 700 : groupelts_center(GEN S)
1414 : {
1415 700 : pari_sp ltop = avma;
1416 700 : long i, j, n = lg(S)-1, l = n;
1417 700 : GEN V, elts = zero_F2v(n+1);
1418 25732 : for(i=1; i<=n; i++)
1419 : {
1420 25032 : if (F2v_coeff(elts,i)) { l--; continue; }
1421 573384 : for(j=1; j<=n; j++)
1422 563192 : if (!perm_commute(gel(S,i),gel(S,j)))
1423 : {
1424 14322 : F2v_set(elts,i);
1425 14322 : F2v_set(elts,j); l--; break;
1426 : }
1427 : }
1428 700 : V = cgetg(l+1,t_VEC);
1429 25732 : for (i=1, j=1; i<=n ;i++)
1430 25032 : if (!F2v_coeff(elts,i)) gel(V,j++) = vecsmall_copy(gel(S,i));
1431 700 : return gerepileupto(ltop,V);
1432 : }
1433 :
1434 : GEN
1435 4270 : groupelts_conjclasses(GEN elts, long *pnbcl)
1436 : {
1437 4270 : long i, j, cl = 0, n = lg(elts)-1;
1438 4270 : GEN c = const_vecsmall(n,0);
1439 4270 : pari_sp av = avma;
1440 52850 : for (i=1; i<=n; i++)
1441 : {
1442 48580 : GEN g = gel(elts,i);
1443 48580 : if (c[i]) continue;
1444 34965 : c[i] = ++cl;
1445 486871 : for(j=1; j<=n; j++)
1446 451906 : if (j != i)
1447 : {
1448 416941 : GEN h = perm_conj(gel(elts,j), g);
1449 416941 : long i2 = gen_search(elts,h,(void*)&vecsmall_lexcmp,&cmp_nodata);
1450 416941 : c[i2] = cl; set_avma(av);
1451 : }
1452 : }
1453 4270 : if (pnbcl) *pnbcl = cl;
1454 4270 : return c;
1455 : }
1456 :
1457 : GEN
1458 4270 : conjclasses_repr(GEN conj, long nb)
1459 : {
1460 4270 : long i, l = lg(conj);
1461 4270 : GEN e = const_vecsmall(nb, 0);
1462 52850 : for(i=1; i<l; i++)
1463 : {
1464 48580 : long ci = conj[i];
1465 48580 : if (!e[ci]) e[ci] = i;
1466 : }
1467 4270 : return e;
1468 : }
1469 :
1470 : /* elts of G sorted wrt vecsmall_lexcmp order: g in G is determined by g[1]
1471 : * so sort by increasing g[1] */
1472 : static GEN
1473 3885 : galois_elts_sorted(GEN gal)
1474 : {
1475 : long i, l;
1476 3885 : GEN elts = gal_get_group(gal), v = cgetg_copy(elts, &l);
1477 43141 : for (i = 1; i < l; i++) { GEN g = gel(elts,i); gel(v, g[1]) = g; }
1478 3885 : return v;
1479 : }
1480 : GEN
1481 4291 : group_to_cc(GEN G)
1482 : {
1483 4291 : GEN elts = checkgroupelts(G), z = cgetg(5,t_VEC);
1484 4270 : long n, flag = 1;
1485 4270 : if (typ(gel(G,1)) == t_POL)
1486 3885 : elts = galois_elts_sorted(G); /* galoisinit */
1487 : else
1488 : {
1489 385 : long i, l = lg(elts);
1490 385 : elts = gen_sort_shallow(elts,(void*)vecsmall_lexcmp,cmp_nodata);
1491 5824 : for (i = 1; i < l; i++)
1492 5586 : if (gel(elts,i)[1] != i) { flag = 0; break; }
1493 : }
1494 4270 : gel(z,1) = elts;
1495 4270 : gel(z,2) = groupelts_conjclasses(elts,&n);
1496 4270 : gel(z,3) = conjclasses_repr(gel(z,2),n);
1497 4270 : gel(z,4) = utoi(flag); return z;
1498 : }
1499 :
1500 : /* S a list of generators */
1501 : GEN
1502 0 : groupelts_abelian_group(GEN S)
1503 : {
1504 0 : pari_sp ltop = avma;
1505 : GEN Qgen, Qord, Qelt;
1506 0 : long i, j, n = lg(gel(S,1))-1, l = lg(S);
1507 0 : Qord = cgetg(l, t_VECSMALL);
1508 0 : Qgen = cgetg(l, t_VEC);
1509 0 : Qelt = mkvec(identity_perm(n));
1510 0 : for (i = 1, j = 1; i < l; ++i)
1511 : {
1512 0 : GEN g = gel(S,i);
1513 0 : long o = perm_relorder(g, groupelts_set(Qelt, n));
1514 0 : gel(Qgen,j) = g;
1515 0 : Qord[j] = o;
1516 0 : if (o != 1) { Qelt = perm_generate(g, Qelt, o); j++; }
1517 : }
1518 0 : setlg(Qgen,j);
1519 0 : setlg(Qord,j);
1520 0 : return gerepilecopy(ltop, mkvec2(Qgen, Qord));
1521 : }
1522 :
1523 : GEN
1524 14 : group_export_GAP(GEN G)
1525 : {
1526 14 : pari_sp av = avma;
1527 14 : GEN s, comma, g = grp_get_gen(G);
1528 14 : long i, k, l = lg(g);
1529 14 : if (l == 1) return strtoGENstr("Group(())");
1530 7 : s = cgetg(2*l, t_VEC);
1531 7 : comma = strtoGENstr(", ");
1532 7 : gel(s,1) = strtoGENstr("Group(");
1533 28 : for (i=1, k=2; i < l; ++i)
1534 : {
1535 21 : if (i > 1) gel(s,k++) = comma;
1536 21 : gel(s,k++) = perm_to_GAP(gel(g,i));
1537 : }
1538 7 : gel(s,k++) = strtoGENstr(")");
1539 7 : return gerepilecopy(av, shallowconcat1(s));
1540 : }
1541 :
1542 : GEN
1543 14 : group_export_MAGMA(GEN G)
1544 : {
1545 14 : pari_sp av = avma;
1546 14 : GEN s, comma, g = grp_get_gen(G);
1547 14 : long i, k, l = lg(g);
1548 14 : if (l == 1) return strtoGENstr("PermutationGroup<1|>");
1549 7 : s = cgetg(2*l, t_VEC);
1550 7 : comma = strtoGENstr(", ");
1551 7 : gel(s,1) = gsprintf("PermutationGroup<%ld|",group_domain(G));
1552 28 : for (i=1, k=2; i < l; ++i)
1553 : {
1554 21 : if (i > 1) gel(s,k++) = comma;
1555 21 : gel(s,k++) = GENtoGENstr( vecsmall_to_vec(gel(g,i)) );
1556 : }
1557 7 : gel(s,k++) = strtoGENstr(">");
1558 7 : return gerepilecopy(av, shallowconcat1(s));
1559 : }
1560 :
1561 : GEN
1562 28 : group_export(GEN G, long format)
1563 : {
1564 28 : switch(format)
1565 : {
1566 14 : case 0: return group_export_GAP(G);
1567 14 : case 1: return group_export_MAGMA(G);
1568 : }
1569 0 : pari_err_FLAG("galoisexport");
1570 0 : return NULL; /*-Wall*/
1571 : }
1572 :
1573 : static GEN
1574 3010 : groupelts_cyclic_subgroups(GEN G)
1575 : {
1576 3010 : pari_sp av = avma;
1577 3010 : long i, j, n = lg(G)-1;
1578 : GEN elts, f, gen, ord;
1579 3010 : if (n==1) return cgetg(1,t_VEC);
1580 3010 : elts = zero_F2v(lg(gel(G,1))-1);
1581 3010 : gen = cgetg(n+1, t_VECSMALL);
1582 3010 : ord = cgetg(n+1, t_VECSMALL);
1583 48958 : for (i=1, j=1; i<=n; i++)
1584 : {
1585 45948 : long k = 1, o, c = 0;
1586 45948 : GEN p = gel(G, i);
1587 45948 : if (F2v_coeff(elts, p[1])) continue;
1588 33439 : o = perm_orderu(p);
1589 33439 : gen[j] = i; ord[j] = o; j++;
1590 : do
1591 : {
1592 90748 : if (cgcd(o, ++c)==1) F2v_set(elts, p[k]);
1593 90748 : k = p[k];
1594 90748 : } while (k!=1);
1595 : }
1596 3010 : setlg(gen, j);
1597 3010 : setlg(ord, j);
1598 3010 : f = vecsmall_indexsort(ord);
1599 3010 : return gerepilecopy(av, mkvec2(vecpermute(gen, f), vecpermute(ord, f)));
1600 : }
1601 :
1602 : GEN
1603 3017 : groupelts_to_group(GEN G)
1604 : {
1605 3017 : pari_sp av = avma;
1606 : GEN L, cyc, ord;
1607 3017 : long i, l, n = lg(G)-1;
1608 3017 : if (n==1) return trivialgroup();
1609 3003 : L = groupelts_cyclic_subgroups(G);
1610 3003 : cyc = gel(L,1); ord = gel(L,2);
1611 3003 : l = lg(cyc);
1612 15085 : for (i = l-1; i >= 2; i--)
1613 : {
1614 14455 : GEN p = gel(G,cyc[i]);
1615 14455 : long o = ord[i];
1616 14455 : GEN H = cyclicgroup(p, o);
1617 14455 : if (o == n) return gerepileupto(av, H);
1618 13447 : if (groupelts_subgroup_isnormal(G, H))
1619 : {
1620 1365 : GEN C = groupelts_quotient(G, H);
1621 1365 : GEN Q = quotient_groupelts(C);
1622 1365 : GEN R = groupelts_to_group(Q);
1623 1365 : if (!R) return gc_NULL(av);
1624 1365 : return gerepilecopy(av, quotient_subgroup_lift(C, H, R));
1625 : }
1626 : }
1627 630 : if (n==12 && l==9 && ord[2]==2 && ord[3]==2 && ord[5]==3)
1628 532 : return gerepilecopy(av,
1629 266 : mkvec2(mkvec3(gel(G,cyc[2]), gel(G,cyc[3]), gel(G,cyc[5])), mkvecsmall3(2,2,3)));
1630 364 : if (n==24 && l==18 && ord[11]==3 && ord[15]==4 && ord[16]==4)
1631 : {
1632 350 : GEN t21 = perm_sqr(gel(G,cyc[15]));
1633 350 : GEN t22 = perm_sqr(gel(G,cyc[16]));
1634 350 : GEN s = perm_mul(t22, gel(G,cyc[15]));
1635 700 : return gerepilecopy(av,
1636 350 : mkvec2(mkvec4(t21,t22, gel(G,cyc[11]), s), mkvecsmall4(2,2,3,2)));
1637 : }
1638 14 : if (n==36 && l==24 && ord[11]==3 && ord[15]==4)
1639 : {
1640 7 : GEN t1 = gel(G,cyc[11]), t3 = gel(G,cyc[15]);
1641 7 : return gerepilecopy(av,
1642 : mkvec2(mkvec3(perm_conj(t3, t1), t1, t3), mkvecsmall3(3,3,4)));
1643 : }
1644 7 : return gc_NULL(av);
1645 : }
1646 :
1647 : static GEN
1648 581 : subg_get_gen(GEN subg) { return gel(subg, 1); }
1649 :
1650 : static GEN
1651 4333 : subg_get_set(GEN subg) { return gel(subg, 2); }
1652 :
1653 : static GEN
1654 399 : groupelt_subg_normalize(GEN elt, GEN subg, GEN cyc)
1655 : {
1656 399 : GEN gen = subg_get_gen(subg), set = subg_get_set(subg);
1657 399 : long i, j, u, n = lg(elt)-1, lgen = lg(gen);
1658 399 : GEN b = F2v_copy(cyc), res = zero_F2v(n);
1659 24339 : for(i = 1; i <= n; i++)
1660 : {
1661 : GEN g;
1662 23940 : if (!F2v_coeff(b, i)) continue;
1663 11186 : g = gel(elt,i);
1664 381731 : for(u=1; u<=n; u++)
1665 381731 : if (g[u]==1) break;
1666 12159 : for(j=1; j<lgen; j++)
1667 : {
1668 11473 : GEN h = gel(elt,gen[j]);
1669 11473 : if (!F2v_coeff(set,g[h[u]])) break;
1670 : }
1671 11186 : if (j < lgen) continue;
1672 686 : F2v_set(res,i);
1673 41846 : for(j=1; j <= n; j++)
1674 41160 : if (F2v_coeff(set, j))
1675 2940 : F2v_clear(b,g[gel(elt,j)[1]]);
1676 : }
1677 399 : return res;
1678 : }
1679 :
1680 : static GEN
1681 7 : triv_subg(GEN elt)
1682 : {
1683 7 : GEN v = cgetg(3, t_VEC);
1684 7 : gel(v,1) = cgetg(1,t_VECSMALL);
1685 7 : gel(v,2) = zero_F2v(lg(elt)-1);
1686 7 : F2v_set(gel(v,2),1);
1687 7 : return v;
1688 : }
1689 :
1690 : static GEN
1691 182 : subg_extend(GEN U, long e, long o, GEN elt)
1692 : {
1693 182 : long i, j, n = lg(elt)-1;
1694 182 : GEN g = gel(elt, e);
1695 182 : GEN gen = vecsmall_append(subg_get_gen(U), e);
1696 182 : GEN set = subg_get_set(U);
1697 182 : GEN Vset = zv_copy(set);
1698 11102 : for(i = 1; i <= n; i++)
1699 10920 : if (F2v_coeff(set, i))
1700 : {
1701 630 : long h = gel(elt, i)[1];
1702 1400 : for(j = 1; j < o; j++)
1703 : {
1704 770 : h = g[h];
1705 770 : F2v_set(Vset, h);
1706 : }
1707 : }
1708 182 : return mkvec2(gen, Vset);
1709 : }
1710 :
1711 : static GEN
1712 217 : cyclic_subg(long e, long o, GEN elt)
1713 : {
1714 217 : long j, n = lg(elt)-1, h = 1;
1715 217 : GEN g = gel(elt, e);
1716 217 : GEN gen = mkvecsmall(e);
1717 217 : GEN set = zero_F2v(n);
1718 217 : F2v_set(set,1);
1719 630 : for(j = 1; j < o; j++)
1720 : {
1721 413 : h = g[h];
1722 413 : F2v_set(set, h);
1723 : }
1724 217 : return mkvec2(gen, set);
1725 : }
1726 :
1727 : static GEN
1728 7 : groupelts_to_regular(GEN elt)
1729 : {
1730 7 : long i, j, n = lg(elt)-1;
1731 7 : GEN V = cgetg(n+1,t_VEC);
1732 427 : for (i=1; i<=n; i++)
1733 : {
1734 420 : pari_sp av = avma;
1735 420 : GEN g = gel(elt, i);
1736 420 : GEN W = cgetg(n+1,t_VEC);
1737 25620 : for(j=1; j<=n; j++)
1738 25200 : gel(W,j) = perm_mul(g, gel(elt,j));
1739 420 : gel(V, i) = gerepileuptoleaf(av,vecvecsmall_indexsort(W));
1740 : }
1741 7 : vecvecsmall_sort_inplace(V, NULL);
1742 7 : return V;
1743 : }
1744 :
1745 : static long
1746 217 : groupelts_pow(GEN elt, long j, long n)
1747 : {
1748 217 : GEN g = gel(elt,j);
1749 217 : long i, h = 1;
1750 847 : for (i=1; i<=n; i++)
1751 630 : h = g[h];
1752 217 : return h;
1753 : }
1754 :
1755 : static GEN
1756 7 : groupelts_cyclic_primepow(GEN elt, GEN *pt_pr, GEN *pt_po)
1757 : {
1758 7 : GEN R = groupelts_cyclic_subgroups(elt);
1759 7 : GEN gen = gel(R,1), ord = gel(R,2);
1760 7 : long i, n = lg(elt)-1, l = lg(gen);
1761 7 : GEN set = zero_F2v(n);
1762 7 : GEN pr = zero_Flv(n);
1763 7 : GEN po = zero_Flv(n);
1764 231 : for (i = 1; i < l; i++)
1765 : {
1766 224 : long h = gen[i];
1767 : ulong p;
1768 224 : if (uisprimepower(ord[i], &p))
1769 : {
1770 217 : F2v_set(set, h);
1771 217 : uel(pr,h) = p;
1772 217 : po[h] = groupelts_pow(elt, h, p);
1773 : }
1774 : }
1775 7 : *pt_pr = pr; *pt_po = po;
1776 7 : return set;
1777 : }
1778 :
1779 : static GEN
1780 7 : all_cyclic_subg(GEN pr, GEN po, GEN elt)
1781 : {
1782 7 : long i, n = lg(pr)-1, m = 0, k = 1;
1783 : GEN W;
1784 427 : for (i=1; i <= n; i++)
1785 420 : m += po[i]==1;
1786 7 : W = cgetg(m+1, t_VEC);
1787 427 : for (i=1; i <= n; i++)
1788 420 : if (po[i]==1)
1789 217 : gel(W, k++) = cyclic_subg(i, pr[i], elt);
1790 7 : return W;
1791 : }
1792 :
1793 : static GEN
1794 7 : groupelts_subgroups_raw(GEN elts)
1795 : {
1796 7 : pari_sp av = avma;
1797 7 : GEN elt = groupelts_to_regular(elts);
1798 7 : GEN pr, po, cyc = groupelts_cyclic_primepow(elt, &pr, &po);
1799 7 : long n = lg(elt)-1;
1800 7 : long i, j, nS = 1;
1801 : GEN S, L;
1802 7 : S = cgetg(1+bigomegau(n)+1, t_VEC);
1803 7 : gel(S, nS++) = mkvec(triv_subg(elt));
1804 7 : gel(S, nS++) = L = all_cyclic_subg(pr, po, elt);
1805 7 : if (DEBUGLEVEL) err_printf("subgroups: level %ld: %ld\n",nS,lg(L)-1);
1806 28 : while (lg(L) > 1)
1807 : {
1808 21 : pari_sp av2 = avma;
1809 21 : long nW = 1, lL = lg(L);
1810 21 : long ng = n;
1811 21 : GEN W = cgetg(1+ng, t_VEC);
1812 420 : for (i=1; i<lL; i++)
1813 : {
1814 399 : GEN U = gel(L, i), set = subg_get_set(U);
1815 399 : GEN G = groupelt_subg_normalize(elt, U, cyc);
1816 3570 : for (j=1; j<nW; j++)
1817 : {
1818 3171 : GEN Wj = subg_get_set(gel(W, j));
1819 3171 : if (F2v_subset(set, Wj))
1820 364 : F2v_negimply_inplace(G, Wj);
1821 : }
1822 24339 : for (j=1; j<=n; j++)
1823 23940 : if(F2v_coeff(G,j))
1824 : {
1825 378 : long p = pr[j];
1826 378 : if (F2v_coeff(set, j)) continue;
1827 182 : if (F2v_coeff(set, po[j]))
1828 : {
1829 182 : GEN U2 = subg_extend(U, j, p, elt);
1830 182 : F2v_negimply_inplace(G, subg_get_set(U2));
1831 182 : if (nW > ng) { ng<<=1; W = vec_lengthen(W, ng); }
1832 182 : gel(W, nW++) = U2;
1833 : }
1834 : }
1835 : }
1836 21 : setlg(W, nW);
1837 21 : if (DEBUGLEVEL) err_printf("subgroups: level %ld: %ld\n",nS,nW-1);
1838 21 : L = W;
1839 21 : if (nW > 1) gel(S, nS++) = L = gerepilecopy(av2, W);
1840 : }
1841 7 : setlg(S, nS);
1842 7 : return gerepilecopy(av, shallowconcat1(S));
1843 : }
1844 :
1845 : static GEN
1846 406 : set_groupelts(GEN S, GEN x)
1847 : {
1848 406 : long i, n = F2v_hamming(x), k=1, m = x[1];
1849 406 : GEN v = cgetg(n+1, t_VEC);
1850 24766 : for (i=1; i<=m; i++)
1851 24360 : if (F2v_coeff(x,i))
1852 2037 : gel(v,k++) = gel(S,i);
1853 406 : return v;
1854 : }
1855 :
1856 : static GEN
1857 7 : subg_to_elts(GEN S, GEN x)
1858 413 : { pari_APPLY_type(t_VEC, set_groupelts(S, gmael(x,i,2))); }
1859 :
1860 : GEN
1861 7 : groupelts_solvablesubgroups(GEN G)
1862 : {
1863 7 : pari_sp av = avma;
1864 7 : GEN S = vecvecsmall_sort(checkgroupelts(G));
1865 7 : GEN L = groupelts_subgroups_raw(S);
1866 7 : return gerepilecopy(av, subg_to_elts(S, L));
1867 : }
|
