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 1848 : vec_isconst(GEN v)
25 : {
26 1848 : long i, l = lg(v);
27 : GEN w;
28 1848 : if (l==1) return 1;
29 1848 : w = gel(v,1);
30 6349 : for(i=2; i<l; i++)
31 5768 : if (!gequal(gel(v,i), w)) return 0;
32 581 : return 1;
33 : }
34 :
35 : int
36 17467 : vecsmall_isconst(GEN v)
37 : {
38 17467 : long i, l = lg(v);
39 : ulong w;
40 17467 : if (l==1) return 1;
41 17467 : w = uel(v,1);
42 30138 : for(i=2; i<l; i++)
43 24310 : if (uel(v,i) != w) return 0;
44 5828 : 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 143927026 : vecsmall_sortspec(GEN v, long n, GEN w)
81 : {
82 143927026 : pari_sp ltop=avma;
83 143927026 : long nx=n>>1, ny=n-nx;
84 : long m, ix, iy;
85 : GEN x, y;
86 143927026 : if (n<=2)
87 : {
88 83071270 : if (n==1)
89 17859329 : w[0]=v[0];
90 65211941 : else if (n==2)
91 : {
92 66974148 : long v0=v[0], v1=v[1];
93 66974148 : if (v0<=v1) { w[0]=v0; w[1]=v1; }
94 2893938 : else { w[0]=v1; w[1]=v0; }
95 : }
96 83071270 : return;
97 : }
98 60855756 : x=new_chunk(nx); y=new_chunk(ny);
99 64055850 : vecsmall_sortspec(v,nx,x);
100 64925145 : vecsmall_sortspec(v+nx,ny,y);
101 290068935 : for (m=0, ix=0, iy=0; ix<nx && iy<ny; )
102 223992664 : if (x[ix]<=y[iy])
103 188259579 : w[m++]=x[ix++];
104 : else
105 35733085 : w[m++]=y[iy++];
106 68293370 : for(;ix<nx;) w[m++]=x[ix++];
107 247901312 : for(;iy<ny;) w[m++]=y[iy++];
108 66076271 : set_avma(ltop);
109 : }
110 :
111 : static long
112 25250598 : vecsmall_sort_max(GEN v)
113 : {
114 25250598 : long i, l = lg(v), max = -1;
115 89804049 : for (i = 1; i < l; i++)
116 88007443 : if (v[i] > max) { max = v[i]; if (max >= l) return -1; }
117 19265235 : else if (v[i] < 0) return -1;
118 1796606 : return max;
119 : }
120 : /* assume 0 <= v[i] <= M. In place. */
121 : void
122 1913002 : vecsmall_counting_sort(GEN v, long M)
123 : {
124 : pari_sp av;
125 : long i, j, k, l;
126 : GEN T;
127 1913002 : if (M == 0) return;
128 1913002 : av = avma; T = new_chunk(M + 1); l = lg(v);
129 9321767 : for (i = 0; i <= M; i++) T[i] = 0;
130 7466665 : for (i = 1; i < l; i++) T[v[i]]++; /* T[j] is # keys = j */
131 9321769 : for (j = 0, k = 1; j <= M; j++)
132 12962442 : for (i = 1; i <= T[j]; i++) v[k++] = j;
133 1913002 : 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 12430 : vecsmall_counting_indexsort(GEN v, long M)
154 : {
155 : pari_sp av;
156 12430 : long i, l = lg(v);
157 : GEN T, p;
158 12430 : if (M == 0 || l <= 2) return identity_zv(l - 1);
159 12416 : p = cgetg(l, t_VECSMALL); av = avma; T = new_chunk(M + 1);
160 55478 : for (i = 0; i <= M; i++) T[i] = 0;
161 14431678 : for (i = 1; i < l; i++) T[v[i]]++; /* T[j] is # keys = j */
162 43062 : for (i = 1; i <= M; i++) T[i] += T[i-1]; /* T[j] is # keys <= j */
163 14431678 : for (i = l-1; i > 0; i--) { p[T[v[i]]] = i; T[v[i]]--; }
164 12416 : return gc_const(av, p);
165 : }
166 :
167 : /* in place sort */
168 : void
169 29649801 : vecsmall_sort(GEN v)
170 : {
171 29649801 : long n = lg(v) - 1, max;
172 29649801 : if (n <= 1) return;
173 22942559 : if ((max = vecsmall_sort_max(v)) >= 0)
174 1913003 : vecsmall_counting_sort(v, max);
175 : else
176 21014033 : vecsmall_sortspec(v+1, n, v+1);
177 : }
178 :
179 : /* cf gen_sortspec */
180 : static GEN
181 7954270 : vecsmall_indexsortspec(GEN v, long n)
182 : {
183 : long nx, ny, m, ix, iy;
184 : GEN x, y, w;
185 7954270 : switch(n)
186 : {
187 152072 : case 1: return mkvecsmall(1);
188 3617625 : case 2: return (v[1] <= v[2])? mkvecsmall2(1,2): mkvecsmall2(2,1);
189 1406058 : case 3:
190 1406058 : if (v[1] <= v[2]) {
191 670678 : if (v[2] <= v[3]) return mkvecsmall3(1,2,3);
192 203037 : return (v[1] <= v[3])? mkvecsmall3(1,3,2)
193 613476 : : mkvecsmall3(3,1,2);
194 : } else {
195 735380 : if (v[1] <= v[3]) return mkvecsmall3(2,1,3);
196 276615 : return (v[2] <= v[3])? mkvecsmall3(2,3,1)
197 836411 : : mkvecsmall3(3,2,1);
198 : }
199 : }
200 2778515 : nx = n>>1; ny = n-nx;
201 2778515 : w = cgetg(n+1,t_VECSMALL);
202 2778523 : x = vecsmall_indexsortspec(v,nx);
203 2778522 : y = vecsmall_indexsortspec(v+nx,ny);
204 33633915 : for (m=1, ix=1, iy=1; ix<=nx && iy<=ny; )
205 30855391 : if (v[x[ix]] <= v[y[iy]+nx])
206 14946861 : w[m++] = x[ix++];
207 : else
208 15908530 : w[m++] = y[iy++]+nx;
209 5302541 : for(;ix<=nx;) w[m++] = x[ix++];
210 5308011 : for(;iy<=ny;) w[m++] = y[iy++]+nx;
211 2778524 : set_avma((pari_sp)w); return w;
212 : }
213 :
214 : /*indirect sort.*/
215 : GEN
216 2409727 : vecsmall_indexsort(GEN v)
217 : {
218 2409727 : long n = lg(v) - 1, max;
219 2409727 : if (n==0) return cgetg(1, t_VECSMALL);
220 2409664 : if ((max = vecsmall_sort_max(v)) >= 0)
221 12430 : return vecsmall_counting_indexsort(v, max);
222 : else
223 2397233 : return vecsmall_indexsortspec(v,n);
224 : }
225 :
226 : /* assume V sorted */
227 : GEN
228 31390 : vecsmall_uniq_sorted(GEN v)
229 : {
230 : long i, j, l;
231 31390 : GEN w = cgetg_copy(v, &l);
232 31390 : if (l == 1) return w;
233 31336 : w[1] = v[1];
234 34309 : for(i = j = 2; i < l; i++)
235 2973 : if (v[i] != w[j-1]) w[j++] = v[i];
236 31336 : stackdummy((pari_sp)(w + l), (pari_sp)(w + j));
237 31336 : setlg(w, j); return w;
238 : }
239 :
240 : GEN
241 16474 : vecsmall_uniq(GEN v)
242 : {
243 16474 : pari_sp av = avma;
244 : long max;
245 16474 : if ((max = vecsmall_sort_max(v)) >= 0)
246 16115 : v = vecsmall_counting_uniq(v, max);
247 : else
248 359 : { v = zv_copy(v); vecsmall_sort(v); v = vecsmall_uniq_sorted(v); }
249 16474 : 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 20193 : vecsmall_duplicate(GEN x)
265 : {
266 20193 : pari_sp av = avma;
267 20193 : GEN p = vecsmall_indexsort(x);
268 20193 : long k, i, r = 0, l = lg(x);
269 20193 : if (l == 1) return gc_long(av, 0);
270 28141 : for (k = x[p[1]], i = 2; i < l; k = x[p[i++]])
271 7948 : if (x[p[i]] == k) { r = p[i]; break; }
272 20193 : return gc_long(av, r);
273 : }
274 :
275 : static int
276 54275 : vecsmall_is1to1spec(GEN v, long n, GEN w)
277 : {
278 54275 : pari_sp av = avma;
279 54275 : long nx = n>>1, ny = n - nx, m, ix, iy;
280 : GEN x, y;
281 54275 : if (n <= 2)
282 : {
283 32660 : if (n == 1) w[0] = v[0];
284 21475 : else if (n==2)
285 : {
286 21475 : long v0 = v[0], v1 = v[1];
287 21475 : if (v0 == v1) return 0;
288 21447 : else if (v0 < v1) { w[0] = v0; w[1] = v1; }
289 4644 : else { w[0] = v1; w[1] = v0; }
290 : }
291 32632 : return 1;
292 : }
293 21615 : x = new_chunk(nx); if (!vecsmall_is1to1spec(v, nx, x)) return 0;
294 21517 : y = new_chunk(ny); if (!vecsmall_is1to1spec(v+nx, ny, y)) return 0;
295 84193 : for (m = ix = iy = 0; ix < nx && iy < ny; )
296 62774 : if (x[ix] == y[iy]) return 0;
297 62725 : else if (x[ix] < y[iy])
298 38053 : w[m++] = x[ix++];
299 : else
300 24672 : w[m++] = y[iy++];
301 23524 : while (ix < nx) w[m++] = x[ix++];
302 53025 : while (iy < ny) w[m++] = y[iy++];
303 21419 : return gc_bool(av, 1);
304 : }
305 :
306 : int
307 11241 : vecsmall_is1to1(GEN V)
308 : {
309 11241 : pari_sp av = avma;
310 : long l;
311 11241 : GEN W = cgetg_copy(V, &l);
312 11241 : return gc_bool(av, l <= 2? 1: vecsmall_is1to1spec(V+1, l, W+1));
313 : }
314 :
315 : /*************************************************************************/
316 : /** **/
317 : /** Routines for handling vectors of VECSMALL **/
318 : /** **/
319 : /*************************************************************************/
320 :
321 : GEN
322 14 : vecvecsmall_sort(GEN x)
323 14 : { return gen_sort(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
324 : GEN
325 365848 : vecvecsmall_sort_shallow(GEN x)
326 365848 : { return gen_sort_shallow(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
327 :
328 : void
329 126 : vecvecsmall_sort_inplace(GEN x, GEN *perm)
330 126 : { gen_sort_inplace(x, (void*)&vecsmall_lexcmp, cmp_nodata, perm); }
331 :
332 : GEN
333 462 : vecvecsmall_sort_uniq(GEN x)
334 462 : { return gen_sort_uniq(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
335 :
336 : GEN
337 861 : vecvecsmall_indexsort(GEN x)
338 861 : { return gen_indexsort(x, (void*)&vecsmall_lexcmp, cmp_nodata); }
339 :
340 : long
341 22934359 : vecvecsmall_search(GEN x, GEN y)
342 22934359 : { return gen_search(x,y,(void*)vecsmall_prefixcmp, cmp_nodata); }
343 :
344 : /* assume x non empty */
345 : long
346 0 : vecvecsmall_max(GEN x)
347 : {
348 0 : long i, l = lg(x), m = vecsmall_max(gel(x,1));
349 0 : for (i = 2; i < l; i++)
350 : {
351 0 : long t = vecsmall_max(gel(x,i));
352 0 : if (t > m) m = t;
353 : }
354 0 : return m;
355 : }
356 :
357 : /*************************************************************************/
358 : /** **/
359 : /** Routines for handling permutations **/
360 : /** **/
361 : /*************************************************************************/
362 :
363 : /* Permutations may be given by
364 : * perm (VECSMALL): a bijection from 1...n to 1...n i-->perm[i]
365 : * cyc (VEC of VECSMALL): a product of disjoint cycles. */
366 :
367 : /* Multiply (compose) two permutations, putting the result in the second one. */
368 : static void
369 21 : perm_mul_inplace2(GEN s, GEN t)
370 : {
371 21 : long i, l = lg(s);
372 525 : for (i = 1; i < l; i++) t[i] = s[t[i]];
373 21 : }
374 :
375 : GEN
376 0 : vecperm_extendschreier(GEN C, GEN v, long n)
377 : {
378 0 : pari_sp av = avma;
379 0 : long mj, lv = lg(v), m = 1, mtested = 1;
380 0 : GEN bit = const_vecsmall(n, 0);
381 0 : GEN cy = cgetg(n+1, t_VECSMALL);
382 0 : GEN sh = const_vec(n, gen_0);
383 0 : for(mj=1; mj<=n; mj++)
384 : {
385 0 : if (isintzero(gel(C,mj))) continue;
386 0 : gel(sh,mj) = gcopy(gel(C,mj));
387 0 : if (bit[mj]) continue;
388 0 : cy[m++] = mj;
389 0 : bit[mj] = 1;
390 : for(;;)
391 0 : {
392 0 : long o, mold = m;
393 0 : for (o = 1; o < lv; o++)
394 : {
395 0 : GEN vo = gel(v,o);
396 : long p;
397 0 : for (p = mtested; p < mold; p++) /* m increases! */
398 : {
399 0 : long j = vo[ cy[p] ];
400 0 : if (!bit[j])
401 : {
402 0 : gel(sh,j) = perm_mul(vo, gel(sh, cy[p]));
403 0 : cy[m++] = j;
404 : }
405 0 : bit[j] = 1;
406 : }
407 : }
408 0 : mtested = mold;
409 0 : if (m == mold) break;
410 : }
411 : }
412 0 : return gerepileupto(av, sh);
413 : }
414 :
415 : /* Orbits of the subgroup generated by v on {1,..,n} */
416 : static GEN
417 1410976 : vecperm_orbits_i(GEN v, long n)
418 : {
419 1410976 : long mj = 1, lv = lg(v), k, l;
420 1410976 : GEN cycle = cgetg(n+1, t_VEC), bit = const_vecsmall(n, 0);
421 8780123 : for (k = 1, l = 1; k <= n;)
422 : {
423 7369097 : pari_sp ltop = avma;
424 7369097 : long m = 1;
425 7369097 : GEN cy = cgetg(n+1, t_VECSMALL);
426 8762453 : for ( ; bit[mj]; mj++) /*empty*/;
427 7369019 : k++; cy[m++] = mj;
428 7369019 : bit[mj++] = 1;
429 : for(;;)
430 2405294 : {
431 9774313 : long o, mold = m;
432 19562914 : for (o = 1; o < lv; o++)
433 : {
434 9788601 : GEN vo = gel(v,o);
435 : long p;
436 30732287 : for (p = 1; p < m; p++) /* m increases! */
437 : {
438 20943686 : long j = vo[ cy[p] ];
439 20943686 : if (!bit[j]) cy[m++] = j;
440 20943686 : bit[j] = 1;
441 : }
442 : }
443 9774313 : if (m == mold) break;
444 2405294 : k += m - mold;
445 : }
446 7369019 : setlg(cy, m);
447 7369141 : gel(cycle,l++) = gerepileuptoleaf(ltop, cy);
448 : }
449 1411026 : setlg(cycle, l); return cycle;
450 : }
451 : /* memory clean version */
452 : GEN
453 4781 : vecperm_orbits(GEN v, long n)
454 : {
455 4781 : pari_sp av = avma;
456 4781 : return gerepilecopy(av, vecperm_orbits_i(v, n));
457 : }
458 :
459 : static int
460 2667 : isperm(GEN v)
461 : {
462 2667 : pari_sp av = avma;
463 2667 : long i, n = lg(v)-1;
464 : GEN w;
465 2667 : if (typ(v) != t_VECSMALL) return 0;
466 2667 : w = zero_zv(n);
467 26411 : for (i=1; i<=n; i++)
468 : {
469 23779 : long d = v[i];
470 23779 : if (d < 1 || d > n || w[d]) return gc_bool(av,0);
471 23744 : w[d] = 1;
472 : }
473 2632 : return gc_bool(av,1);
474 : }
475 :
476 : /* Compute the cyclic decomposition of a permutation */
477 : GEN
478 13312 : perm_cycles(GEN v)
479 : {
480 13312 : pari_sp av = avma;
481 13312 : return gerepilecopy(av, vecperm_orbits_i(mkvec(v), lg(v)-1));
482 : }
483 :
484 : GEN
485 259 : permcycles(GEN v)
486 : {
487 259 : if (!isperm(v)) pari_err_TYPE("permcycles",v);
488 252 : return perm_cycles(v);
489 : }
490 :
491 : /* Output the order of p */
492 : ulong
493 436353 : perm_orderu(GEN v)
494 : {
495 436353 : pari_sp av = avma;
496 436353 : GEN c = vecperm_orbits_i(mkvec(v), lg(v)-1);
497 : long i, d;
498 3169270 : for(i=1, d=1; i<lg(c); i++) d = ulcm(d, lg(gel(c,i))-1);
499 436389 : return gc_ulong(av,d);
500 : }
501 :
502 : static GEN
503 2002 : _domul(void *data, GEN x, GEN y)
504 : {
505 2002 : GEN (*mul)(GEN,GEN)=(GEN (*)(GEN,GEN)) data;
506 2002 : return mul(x,y);
507 : }
508 :
509 : /* Output the order of p */
510 : GEN
511 427 : perm_order(GEN v)
512 : {
513 427 : pari_sp av = avma;
514 427 : GEN c = vecperm_orbits_i(mkvec(v), lg(v)-1);
515 427 : long i, l = lg(c);
516 427 : GEN V = cgetg(l, t_VEC);
517 2856 : for (i = 1; i < l; i++)
518 2429 : gel(V,i) = utoi(lg(gel(c,i))-1);
519 427 : return gerepileuptoint(av, gen_product(V, (void *)lcmii, _domul));
520 : }
521 :
522 : GEN
523 434 : permorder(GEN v)
524 : {
525 434 : if (!isperm(v)) pari_err_TYPE("permorder",v);
526 427 : return perm_order(v);
527 : }
528 :
529 : /* sign of a permutation */
530 : long
531 956094 : perm_sign(GEN v)
532 : {
533 956094 : pari_sp av = avma;
534 956094 : GEN c = vecperm_orbits_i(mkvec(v), lg(v)-1);
535 956118 : long i, l = lg(c), s = 1;
536 5533831 : for (i = 1; i < l; i++)
537 4577713 : if (odd(lg(gel(c, i)))) s = -s;
538 956118 : return gc_long(av,s);
539 : }
540 :
541 : long
542 273 : permsign(GEN v)
543 : {
544 273 : if (!isperm(v)) pari_err_TYPE("permsign",v);
545 259 : return perm_sign(v);
546 : }
547 :
548 : GEN
549 5915 : Z_to_perm(long n, GEN x)
550 : {
551 : pari_sp av;
552 : ulong i, r;
553 5915 : GEN v = cgetg(n+1, t_VECSMALL);
554 5915 : if (n==0) return v;
555 5908 : uel(v,n) = 1; av = avma;
556 5908 : if (signe(x) <= 0) x = modii(x, mpfact(n));
557 27146 : for (r=n-1; r>=1; r--)
558 : {
559 : ulong a;
560 21238 : x = absdiviu_rem(x, n+1-r, &a);
561 71687 : for (i=r+1; i<=(ulong)n; i++)
562 50449 : if (uel(v,i) > a) uel(v,i)++;
563 21238 : uel(v,r) = a+1;
564 : }
565 5908 : return gc_const(av, v);
566 : }
567 : GEN
568 5915 : numtoperm(long n, GEN x)
569 : {
570 5915 : if (n < 0) pari_err_DOMAIN("numtoperm", "n", "<", gen_0, stoi(n));
571 5915 : if (typ(x) != t_INT) pari_err_TYPE("numtoperm",x);
572 5915 : return Z_to_perm(n, x);
573 : }
574 :
575 : /* destroys v */
576 : static GEN
577 1701 : perm_to_Z_inplace(GEN v)
578 : {
579 1701 : long l = lg(v), i, r;
580 1701 : GEN x = gen_0;
581 1701 : if (!isperm(v)) return NULL;
582 10143 : for (i = 1; i < l; i++)
583 : {
584 8449 : long vi = v[i];
585 8449 : if (vi <= 0) return NULL;
586 8449 : x = i==1 ? utoi(vi-1): addiu(muliu(x,l-i), vi-1);
587 25396 : for (r = i+1; r < l; r++)
588 16947 : if (v[r] > vi) v[r]--;
589 : }
590 1694 : return x;
591 : }
592 : GEN
593 1680 : perm_to_Z(GEN v)
594 : {
595 1680 : pari_sp av = avma;
596 1680 : GEN x = perm_to_Z_inplace(leafcopy(v));
597 1680 : if (!x) pari_err_TYPE("permtonum",v);
598 1680 : return gerepileuptoint(av, x);
599 : }
600 : GEN
601 1708 : permtonum(GEN p)
602 : {
603 1708 : pari_sp av = avma;
604 : GEN v, x;
605 1708 : switch(typ(p))
606 : {
607 1680 : case t_VECSMALL: return perm_to_Z(p);
608 21 : case t_VEC: case t_COL:
609 21 : if (RgV_is_ZV(p)) { v = ZV_to_zv(p); break; }
610 7 : default: pari_err_TYPE("permtonum",p);
611 : return NULL;/*LCOV_EXCL_LINE*/
612 : }
613 21 : x = perm_to_Z_inplace(v);
614 21 : if (!x) pari_err_TYPE("permtonum",p);
615 14 : return gerepileuptoint(av, x);
616 : }
617 :
618 : GEN
619 7300 : cyc_pow(GEN cyc, long exp)
620 : {
621 : long i, j, k, l, r;
622 : GEN c;
623 22167 : for (r = j = 1; j < lg(cyc); j++)
624 : {
625 14867 : long n = lg(gel(cyc,j)) - 1;
626 14867 : r += cgcd(n, exp);
627 : }
628 7300 : c = cgetg(r, t_VEC);
629 22167 : for (r = j = 1; j < lg(cyc); j++)
630 : {
631 14867 : GEN v = gel(cyc,j);
632 14867 : long n = lg(v) - 1, e = umodsu(exp,n), g = (long)ugcd(n, e), m = n / g;
633 31652 : for (i = 0; i < g; i++)
634 : {
635 16785 : GEN p = cgetg(m+1, t_VECSMALL);
636 16785 : gel(c,r++) = p;
637 54716 : for (k = 1, l = i; k <= m; k++)
638 : {
639 37931 : p[k] = v[l+1];
640 37931 : l += e; if (l >= n) l -= n;
641 : }
642 : }
643 : }
644 7300 : return c;
645 : }
646 :
647 : /* Compute the power of a permutation given by product of cycles
648 : * Ouput a perm, not a cyc */
649 : GEN
650 0 : cyc_pow_perm(GEN cyc, long exp)
651 : {
652 : long e, j, k, l, n;
653 : GEN p;
654 0 : for (n = 0, j = 1; j < lg(cyc); j++) n += lg(gel(cyc,j))-1;
655 0 : p = cgetg(n + 1, t_VECSMALL);
656 0 : for (j = 1; j < lg(cyc); j++)
657 : {
658 0 : GEN v = gel(cyc,j);
659 0 : n = lg(v) - 1; e = umodsu(exp, n);
660 0 : for (k = 1, l = e; k <= n; k++)
661 : {
662 0 : p[v[k]] = v[l+1];
663 0 : if (++l == n) l = 0;
664 : }
665 : }
666 0 : return p;
667 : }
668 :
669 : GEN
670 56 : perm_pow(GEN perm, GEN exp)
671 : {
672 56 : long i, r = lg(perm)-1;
673 56 : GEN p = zero_zv(r);
674 56 : pari_sp av = avma;
675 56 : GEN v = cgetg(r+1, t_VECSMALL);
676 196 : for (i=1; i<=r; i++)
677 : {
678 : long e, n, k, l;
679 140 : if (p[i]) continue;
680 56 : v[1] = i;
681 140 : for (n=1, k=perm[i]; k!=i; k=perm[k], n++) v[n+1] = k;
682 56 : e = umodiu(exp, n);
683 196 : for (k = 1, l = e; k <= n; k++)
684 : {
685 140 : p[v[k]] = v[l+1];
686 140 : if (++l == n) l = 0;
687 : }
688 : }
689 56 : return gc_const(av, p);
690 : }
691 :
692 : GEN
693 18690 : perm_powu(GEN perm, ulong exp)
694 : {
695 18690 : ulong i, r = lg(perm)-1;
696 18690 : GEN p = zero_zv(r);
697 18690 : pari_sp av = avma;
698 18690 : GEN v = cgetg(r+1, t_VECSMALL);
699 246540 : for (i=1; i<=r; i++)
700 : {
701 : ulong e, n, k, l;
702 227850 : if (p[i]) continue;
703 84707 : v[1] = i;
704 227850 : for (n=1, k=perm[i]; k!=i; k=perm[k], n++) v[n+1] = k;
705 84707 : e = exp % n;
706 312557 : for (k = 1, l = e; k <= n; k++)
707 : {
708 227850 : p[v[k]] = v[l+1];
709 227850 : if (++l == n) l = 0;
710 : }
711 : }
712 18690 : return gc_const(av, p);
713 : }
714 :
715 : GEN
716 21 : perm_to_GAP(GEN p)
717 : {
718 21 : pari_sp ltop=avma;
719 : GEN gap;
720 : GEN x;
721 : long i;
722 21 : long nb, c=0;
723 : char *s;
724 : long sz;
725 21 : long lp=lg(p)-1;
726 21 : if (typ(p) != t_VECSMALL) pari_err_TYPE("perm_to_GAP",p);
727 21 : x = perm_cycles(p);
728 21 : sz = (long) ((bfffo(lp)+1) * LOG10_2 + 1);
729 : /*Dry run*/
730 133 : for (i = 1, nb = 1; i < lg(x); ++i)
731 : {
732 112 : GEN z = gel(x,i);
733 112 : long lz = lg(z)-1;
734 112 : nb += 1+lz*(sz+2);
735 : }
736 21 : nb++;
737 : /*Real run*/
738 21 : gap = cgetg(nchar2nlong(nb) + 1, t_STR);
739 21 : s = GSTR(gap);
740 133 : for (i = 1; i < lg(x); ++i)
741 : {
742 : long j;
743 112 : GEN z = gel(x,i);
744 112 : if (lg(z) > 2)
745 : {
746 112 : s[c++] = '(';
747 364 : for (j = 1; j < lg(z); ++j)
748 : {
749 252 : if (j > 1)
750 : {
751 140 : s[c++] = ','; s[c++] = ' ';
752 : }
753 252 : sprintf(s+c,"%ld",z[j]);
754 567 : while(s[c++]) /* empty */;
755 252 : c--;
756 : }
757 112 : s[c++] = ')';
758 : }
759 : }
760 21 : if (!c) { s[c++]='('; s[c++]=')'; }
761 21 : s[c] = '\0';
762 21 : return gerepileupto(ltop,gap);
763 : }
764 :
765 : int
766 572495 : perm_commute(GEN s, GEN t)
767 : {
768 572495 : long i, l = lg(t);
769 40373487 : for (i = 1; i < l; i++)
770 39820382 : if (t[ s[i] ] != s[ t[i] ]) return 0;
771 553105 : return 1;
772 : }
773 :
774 : /*************************************************************************/
775 : /** **/
776 : /** Routines for handling groups **/
777 : /** **/
778 : /*************************************************************************/
779 : /* A Group is a t_VEC [gen,orders]
780 : * gen (vecvecsmall): list of generators given by permutations
781 : * orders (vecsmall): relatives orders of generators. */
782 942582 : INLINE GEN grp_get_gen(GEN G) { return gel(G,1); }
783 1597468 : INLINE GEN grp_get_ord(GEN G) { return gel(G,2); }
784 :
785 : /* A Quotient Group is a t_VEC [gen,coset]
786 : * gen (vecvecsmall): coset generators
787 : * coset (vecsmall): gen[coset[p[1]]] generate the p-coset.
788 : */
789 141820 : INLINE GEN quo_get_gen(GEN C) { return gel(C,1); }
790 30058 : INLINE GEN quo_get_coset(GEN C) { return gel(C,2); }
791 :
792 : static GEN
793 52458 : trivialsubgroups(void)
794 52458 : { GEN L = cgetg(2, t_VEC); gel(L,1) = trivialgroup(); return L; }
795 :
796 : /* Compute the order of p modulo the group given by a set */
797 : long
798 220220 : perm_relorder(GEN p, GEN set)
799 : {
800 220220 : pari_sp ltop = avma;
801 220220 : long n = 1, q = p[1];
802 654170 : while (!F2v_coeff(set,q)) { q = p[q]; n++; }
803 220219 : return gc_long(ltop,n);
804 : }
805 :
806 : GEN
807 13076 : perm_generate(GEN S, GEN H, long o)
808 : {
809 13076 : long i, n = lg(H)-1;
810 13076 : GEN L = cgetg(n*o + 1, t_VEC);
811 45885 : for(i=1; i<=n; i++) gel(L,i) = vecsmall_copy(gel(H,i));
812 50673 : for( ; i <= n*o; i++) gel(L,i) = perm_mul(gel(L,i-n), S);
813 13076 : return L;
814 : }
815 :
816 : /*Return the order (cardinality) of a group */
817 : long
818 711039 : group_order(GEN G)
819 : {
820 711039 : return zv_prod(grp_get_ord(G));
821 : }
822 :
823 : /* G being a subgroup of S_n, output n */
824 : long
825 26684 : group_domain(GEN G)
826 : {
827 26684 : GEN gen = grp_get_gen(G);
828 26684 : if (lg(gen) < 2) pari_err_DOMAIN("group_domain", "#G", "=", gen_1,G);
829 26684 : return lg(gel(gen,1)) - 1;
830 : }
831 :
832 : /*Left coset of g mod G: gG*/
833 : GEN
834 304430 : group_leftcoset(GEN G, GEN g)
835 : {
836 304430 : GEN gen = grp_get_gen(G), ord = grp_get_ord(G);
837 304430 : GEN res = cgetg(group_order(G)+1, t_VEC);
838 : long i, j, k;
839 304430 : gel(res,1) = vecsmall_copy(g);
840 304430 : k = 1;
841 560259 : for (i = 1; i < lg(gen); i++)
842 : {
843 255829 : long c = k * (ord[i] - 1);
844 704683 : for (j = 1; j <= c; j++) gel(res,++k) = perm_mul(gel(res,j), gel(gen,i));
845 : }
846 304430 : return res;
847 : }
848 : /*Right coset of g mod G: Gg*/
849 : GEN
850 182245 : group_rightcoset(GEN G, GEN g)
851 : {
852 182245 : GEN gen = grp_get_gen(G), ord = grp_get_ord(G);
853 182245 : GEN res = cgetg(group_order(G)+1, t_VEC);
854 : long i, j, k;
855 182245 : gel(res,1) = vecsmall_copy(g);
856 182245 : k = 1;
857 315553 : for (i = 1; i < lg(gen); i++)
858 : {
859 133308 : long c = k * (ord[i] - 1);
860 419517 : for (j = 1; j <= c; j++) gel(res,++k) = perm_mul(gel(gen,i), gel(res,j));
861 : }
862 182245 : return res;
863 : }
864 : /*Elements of a group from the generators, cf group_leftcoset*/
865 : GEN
866 140812 : group_elts(GEN G, long n)
867 : {
868 140812 : if (lg(G)==3 && typ(gel(G,1))==t_VEC)
869 : {
870 140812 : GEN gen = grp_get_gen(G), ord = grp_get_ord(G);
871 140812 : GEN res = cgetg(group_order(G)+1, t_VEC);
872 : long i, j, k;
873 140812 : gel(res,1) = identity_perm(n);
874 140812 : k = 1;
875 285397 : for (i = 1; i < lg(gen); i++)
876 : {
877 144585 : long c = k * (ord[i] - 1);
878 : /* j = 1, use res[1] = identity */
879 144585 : gel(res,++k) = vecsmall_copy(gel(gen,i));
880 384776 : for (j = 2; j <= c; j++) gel(res,++k) = perm_mul(gel(res,j), gel(gen,i));
881 : }
882 140812 : return res;
883 0 : } else return gcopy(G);
884 : }
885 :
886 : GEN
887 14448 : groupelts_conj_set(GEN elts, GEN p)
888 : {
889 14448 : long i, j, l = lg(elts), n = lg(p)-1;
890 14448 : GEN res = zero_F2v(n);
891 241465 : for(j = 1; j < n; j++)
892 241465 : if (p[j]==1) break;
893 101136 : for(i = 1; i < l; i++)
894 86688 : F2v_set(res, p[mael(elts,i,j)]);
895 14448 : return res;
896 : }
897 :
898 : GEN
899 28098 : groupelts_set(GEN elts, long n)
900 : {
901 28098 : GEN res = zero_F2v(n);
902 28098 : long i, l = lg(elts);
903 137606 : for(i=1; i<l; i++)
904 109508 : F2v_set(res,mael(elts,i,1));
905 28098 : return res;
906 : }
907 :
908 : /*Elements of a group from the generators, returned as a set (bitmap)*/
909 : GEN
910 90699 : group_set(GEN G, long n)
911 : {
912 90699 : GEN res = zero_F2v(n);
913 90699 : pari_sp av = avma;
914 90699 : GEN elts = group_elts(G, n);
915 90699 : long i, l = lg(elts);
916 284270 : for(i=1; i<l; i++)
917 193571 : F2v_set(res,mael(elts,i,1));
918 90699 : return gc_const(av, res);
919 : }
920 :
921 : static int
922 17353 : sgcmp(GEN a, GEN b) { return vecsmall_lexcmp(gel(a,1),gel(b,1)); }
923 :
924 : GEN
925 497 : subgroups_tableset(GEN S, long n)
926 : {
927 497 : long i, l = lg(S);
928 497 : GEN v = cgetg(l, t_VEC);
929 5411 : for(i=1; i<l; i++)
930 4914 : gel(v,i) = mkvec2(group_set(gel(S,i), n), mkvecsmall(i));
931 497 : gen_sort_inplace(v,(void*)sgcmp,cmp_nodata, NULL);
932 497 : return v;
933 : }
934 :
935 : long
936 2002 : tableset_find_index(GEN tbl, GEN set)
937 : {
938 2002 : long i = tablesearch(tbl,mkvec2(set,mkvecsmall(0)),sgcmp);
939 2002 : if (!i) return 0;
940 2002 : return mael3(tbl,i,2,1);
941 : }
942 :
943 : GEN
944 52486 : trivialgroup(void) { retmkvec2(cgetg(1,t_VEC), cgetg(1,t_VECSMALL)); }
945 :
946 : /*Cyclic group generated by g of order s*/
947 : GEN
948 27916 : cyclicgroup(GEN g, long s)
949 27916 : { retmkvec2(mkvec( vecsmall_copy(g) ), mkvecsmall(s)); }
950 :
951 : /*Return the group generated by g1,g2 of relative orders s1,s2*/
952 : GEN
953 1085 : dicyclicgroup(GEN g1, GEN g2, long s1, long s2)
954 1085 : { retmkvec2( mkvec2(vecsmall_copy(g1), vecsmall_copy(g2)),
955 : mkvecsmall2(s1, s2) ); }
956 :
957 : /* return the quotient map G --> G/H */
958 : /*The ouput is [gen,hash]*/
959 : /* gen (vecvecsmall): coset generators
960 : * coset (vecsmall): vecsmall of coset number) */
961 : GEN
962 11725 : groupelts_quotient(GEN elt, GEN H)
963 : {
964 11725 : pari_sp ltop = avma;
965 : GEN p2, p3;
966 11725 : long i, j, a = 1;
967 11725 : long n = lg(gel(elt,1))-1, o = group_order(H);
968 : GEN el;
969 11725 : long le = lg(elt)-1;
970 11725 : GEN used = zero_F2v(le+1);
971 11725 : long l = le/o;
972 11725 : p2 = cgetg(l+1, t_VEC);
973 11725 : p3 = zero_zv(n);
974 11725 : el = zero_zv(n);
975 150073 : for (i = 1; i<=le; i++)
976 138348 : el[mael(elt,i,1)]=i;
977 68656 : for (i = 1; i <= l; ++i)
978 : {
979 : GEN V;
980 150619 : while(F2v_coeff(used,a)) a++;
981 56938 : V = group_leftcoset(H,gel(elt,a));
982 56938 : gel(p2,i) = gel(V,1);
983 195181 : for(j=1;j<lg(V);j++)
984 : {
985 138250 : long b = el[mael(V,j,1)];
986 138250 : if (b==0) pari_err_IMPL("group_quotient for a non-WSS group");
987 138243 : F2v_set(used,b);
988 : }
989 195167 : for (j = 1; j <= o; j++)
990 138236 : p3[mael(V, j, 1)] = i;
991 : }
992 11718 : return gerepilecopy(ltop,mkvec2(p2,p3));
993 : }
994 :
995 : GEN
996 10213 : group_quotient(GEN G, GEN H)
997 : {
998 10213 : return groupelts_quotient(group_elts(G, group_domain(G)), H);
999 : }
1000 :
1001 : /*Compute the image of a permutation by a quotient map.*/
1002 : GEN
1003 30058 : quotient_perm(GEN C, GEN p)
1004 : {
1005 30058 : GEN gen = quo_get_gen(C);
1006 30058 : GEN coset = quo_get_coset(C);
1007 30058 : long j, l = lg(gen);
1008 30058 : GEN p3 = cgetg(l, t_VECSMALL);
1009 283185 : for (j = 1; j < l; ++j)
1010 : {
1011 253127 : p3[j] = coset[p[mael(gen,j,1)]];
1012 253127 : if (p3[j]==0) pari_err_IMPL("quotient_perm for a non-WSS group");
1013 : }
1014 30058 : return p3;
1015 : }
1016 :
1017 : /* H is a subgroup of G, C is the quotient map G --> G/H
1018 : *
1019 : * Lift a subgroup S of G/H to a subgroup of G containing H */
1020 : GEN
1021 50778 : quotient_subgroup_lift(GEN C, GEN H, GEN S)
1022 : {
1023 50778 : GEN genH = grp_get_gen(H);
1024 50778 : GEN genS = grp_get_gen(S);
1025 50778 : GEN genC = quo_get_gen(C);
1026 50778 : long l1 = lg(genH)-1;
1027 50778 : long l2 = lg(genS)-1, j;
1028 50778 : GEN p1 = cgetg(3, t_VEC), L = cgetg(l1+l2+1, t_VEC);
1029 101724 : for (j = 1; j <= l1; ++j) gel(L,j) = gel(genH,j);
1030 118125 : for (j = 1; j <= l2; ++j) gel(L,l1+j) = gel(genC, mael(genS,j,1));
1031 50778 : gel(p1,1) = L;
1032 50778 : gel(p1,2) = vecsmall_concat(grp_get_ord(H), grp_get_ord(S));
1033 50778 : return p1;
1034 : }
1035 :
1036 : /* Let G a group and C a quotient map G --> G/H
1037 : * Assume H is normal, return the group G/H */
1038 : GEN
1039 10206 : quotient_group(GEN C, GEN G)
1040 : {
1041 10206 : pari_sp ltop = avma;
1042 : GEN Qgen, Qord, Qelt, Qset, Q;
1043 10206 : GEN Cgen = quo_get_gen(C);
1044 10206 : GEN Ggen = grp_get_gen(G);
1045 10206 : long i,j, n = lg(Cgen)-1, l = lg(Ggen);
1046 10206 : Qord = cgetg(l, t_VECSMALL);
1047 10206 : Qgen = cgetg(l, t_VEC);
1048 10206 : Qelt = mkvec(identity_perm(n));
1049 10206 : Qset = groupelts_set(Qelt, n);
1050 31164 : for (i = 1, j = 1; i < l; ++i)
1051 : {
1052 20958 : GEN g = quotient_perm(C, gel(Ggen,i));
1053 20958 : long o = perm_relorder(g, Qset);
1054 20958 : gel(Qgen,j) = g;
1055 20958 : Qord[j] = o;
1056 20958 : if (o != 1)
1057 : {
1058 13076 : Qelt = perm_generate(g, Qelt, o);
1059 13076 : Qset = groupelts_set(Qelt, n);
1060 13076 : j++;
1061 : }
1062 : }
1063 10206 : setlg(Qgen,j);
1064 10206 : setlg(Qord,j); Q = mkvec2(Qgen, Qord);
1065 10206 : return gerepilecopy(ltop,Q);
1066 : }
1067 :
1068 : GEN
1069 1512 : quotient_groupelts(GEN C)
1070 : {
1071 1512 : GEN G = quo_get_gen(C);
1072 1512 : long i, l = lg(G);
1073 1512 : GEN Q = cgetg(l, t_VEC);
1074 10612 : for (i = 1; i < l; ++i)
1075 9100 : gel(Q,i) = quotient_perm(C, gel(G,i));
1076 1512 : return Q;
1077 : }
1078 :
1079 : /* Return 1 if g normalizes N, 0 otherwise */
1080 : long
1081 182245 : group_perm_normalize(GEN N, GEN g)
1082 : {
1083 182245 : pari_sp ltop = avma;
1084 182245 : long r = gequal(vecvecsmall_sort_shallow(group_leftcoset(N, g)),
1085 : vecvecsmall_sort_shallow(group_rightcoset(N, g)));
1086 182245 : return gc_long(ltop, r);
1087 : }
1088 :
1089 : /* L is a list of subgroups, C is a coset and r a relative order.*/
1090 : static GEN
1091 65247 : liftlistsubgroups(GEN L, GEN C, long r)
1092 : {
1093 65247 : pari_sp ltop = avma;
1094 65247 : long c = lg(C)-1, l = lg(L)-1, n = lg(gel(C,1))-1, i, k;
1095 : GEN R;
1096 65247 : if (!l) return cgetg(1,t_VEC);
1097 58555 : R = cgetg(l*c+1, t_VEC);
1098 143024 : for (i = 1, k = 1; i <= l; ++i)
1099 : {
1100 84469 : GEN S = gel(L,i), Selt = group_set(S,n);
1101 84469 : GEN gen = grp_get_gen(S);
1102 84469 : GEN ord = grp_get_ord(S);
1103 : long j;
1104 279286 : for (j = 1; j <= c; ++j)
1105 : {
1106 194817 : GEN p = gel(C,j);
1107 194817 : if (perm_relorder(p, Selt) == r && group_perm_normalize(S, p))
1108 108731 : gel(R,k++) = mkvec2(vec_append(gen, p),
1109 : vecsmall_append(ord, r));
1110 : }
1111 : }
1112 58555 : setlg(R, k);
1113 58555 : return gerepilecopy(ltop, R);
1114 : }
1115 :
1116 : /* H is a normal subgroup, C is the quotient map G -->G/H,
1117 : * S is a subgroup of G/H, and G is embedded in Sym(l)
1118 : * Return all the subgroups K of G such that
1119 : * S= K mod H and K inter H={1} */
1120 : static GEN
1121 49266 : liftsubgroup(GEN C, GEN H, GEN S)
1122 : {
1123 49266 : pari_sp ltop = avma;
1124 49266 : GEN V = trivialsubgroups();
1125 49266 : GEN Sgen = grp_get_gen(S);
1126 49266 : GEN Sord = grp_get_ord(S);
1127 49266 : GEN Cgen = quo_get_gen(C);
1128 49266 : long n = lg(Sgen), i;
1129 114513 : for (i = 1; i < n; ++i)
1130 : { /*loop over generators of S*/
1131 65247 : GEN W = group_leftcoset(H, gel(Cgen, mael(Sgen, i, 1)));
1132 65247 : V = liftlistsubgroups(V, W, Sord[i]);
1133 : }
1134 49266 : return gerepilecopy(ltop,V);
1135 : }
1136 :
1137 : /* 1:A4, 2:S4, 3:F36, 0: other */
1138 : long
1139 10038 : group_isA4S4(GEN G)
1140 : {
1141 10038 : GEN elt = grp_get_gen(G);
1142 10038 : GEN ord = grp_get_ord(G);
1143 10038 : long n = lg(ord);
1144 10038 : if (n != 4 && n != 5) return 0;
1145 2219 : if (n==4 && ord[1]==3 && ord[2]==3 && ord[3]==4)
1146 : {
1147 : long i;
1148 7 : GEN p = gel(elt,1), q = gel(elt,2), r = gel(elt,3);
1149 259 : for(i=1; i<=36; i++)
1150 252 : if (p[r[i]]!=r[q[i]]) return 0;
1151 7 : return 3;
1152 : }
1153 2212 : if (ord[1]!=2 || ord[2]!=2 || ord[3]!=3) return 0;
1154 42 : if (perm_commute(gel(elt,1),gel(elt,3))) return 0;
1155 42 : if (n==4) return 1;
1156 21 : if (ord[4]!=2) return 0;
1157 21 : if (perm_commute(gel(elt,3),gel(elt,4))) return 0;
1158 21 : return 2;
1159 : }
1160 : /* compute all the subgroups of a group G */
1161 : GEN
1162 13230 : group_subgroups(GEN G)
1163 : {
1164 13230 : pari_sp ltop = avma;
1165 : GEN p1, H, C, Q, M, sg1, sg2, sg3;
1166 13230 : GEN gen = grp_get_gen(G);
1167 13230 : GEN ord = grp_get_ord(G);
1168 13230 : long lM, i, j, n = lg(gen);
1169 : long t;
1170 13230 : if (n == 1) return trivialsubgroups();
1171 10038 : t = group_isA4S4(G);
1172 10038 : if (t == 3)
1173 : {
1174 7 : GEN H = mkvec2(mkvec3(gel(gen,1), gel(gen,2), perm_sqr(gel(gen,3))),
1175 : mkvecsmall3(3, 3, 2));
1176 7 : GEN S = group_subgroups(H);
1177 7 : GEN V = cgetg(11,t_VEC);
1178 7 : gel(V,1) = cyclicgroup(gel(gen,3),4);
1179 63 : for (i=2; i<10; i++)
1180 56 : gel(V,i) = cyclicgroup(perm_mul(gmael3(V,i-1,1,1),gel(gen,i%3==1 ? 2:1)),4);
1181 7 : gel(V,10) = G;
1182 7 : return gerepilecopy(ltop,shallowconcat(S,V));
1183 : }
1184 10031 : else if (t)
1185 : {
1186 42 : GEN s = gel(gen,1); /*s = (1,2)(3,4) */
1187 42 : GEN t = gel(gen,2); /*t = (1,3)(2,4) */
1188 42 : GEN st = perm_mul(s, t); /*st = (1,4)(2,3) */
1189 42 : H = dicyclicgroup(s, t, 2, 2);
1190 : /* sg3 is the list of subgroups intersecting only partially with H*/
1191 42 : sg3 = cgetg((n==4)?4: 10, t_VEC);
1192 42 : gel(sg3,1) = cyclicgroup(s, 2);
1193 42 : gel(sg3,2) = cyclicgroup(t, 2);
1194 42 : gel(sg3,3) = cyclicgroup(st, 2);
1195 42 : if (n==5)
1196 : {
1197 21 : GEN u = gel(gen,3);
1198 21 : GEN v = gel(gen,4), w, u2;
1199 21 : if (zv_equal(perm_conj(u,s), t)) /*u=(2,3,4)*/
1200 21 : u2 = perm_sqr(u);
1201 : else
1202 : {
1203 0 : u2 = u;
1204 0 : u = perm_sqr(u);
1205 : }
1206 21 : if (perm_orderu(v)==2)
1207 : {
1208 21 : if (!perm_commute(s,v)) /*v=(1,2)*/
1209 : {
1210 0 : v = perm_conj(u,v);
1211 0 : if (!perm_commute(s,v)) v = perm_conj(u,v);
1212 : }
1213 21 : w = perm_mul(v,t); /*w=(1,4,2,3)*/
1214 : }
1215 : else
1216 : {
1217 0 : w = v;
1218 0 : if (!zv_equal(perm_sqr(w), s)) /*w=(1,4,2,3)*/
1219 : {
1220 0 : w = perm_conj(u,w);
1221 0 : if (!zv_equal(perm_sqr(w), s)) w = perm_conj(u,w);
1222 : }
1223 0 : v = perm_mul(w,t); /*v=(1,2)*/
1224 : }
1225 21 : gel(sg3,4) = dicyclicgroup(s,v,2,2);
1226 21 : gel(sg3,5) = dicyclicgroup(t,perm_conj(u,v),2,2);
1227 21 : gel(sg3,6) = dicyclicgroup(st,perm_conj(u2,v),2,2);
1228 21 : gel(sg3,7) = dicyclicgroup(s,w,2,2);
1229 21 : gel(sg3,8) = dicyclicgroup(t,perm_conj(u,w),2,2);
1230 21 : gel(sg3,9) = dicyclicgroup(st,perm_conj(u2,w),2,2);
1231 : }
1232 : }
1233 : else
1234 : {
1235 9989 : ulong osig = mael(factoru(ord[1]), 1, 1);
1236 9989 : GEN sig = perm_powu(gel(gen,1), ord[1]/osig);
1237 9989 : H = cyclicgroup(sig,osig);
1238 9989 : sg3 = NULL;
1239 : }
1240 10031 : C = group_quotient(G,H);
1241 10024 : Q = quotient_group(C,G);
1242 10024 : M = group_subgroups(Q); lM = lg(M);
1243 : /* sg1 is the list of subgroups containing H*/
1244 10017 : sg1 = cgetg(lM, t_VEC);
1245 59283 : for (i = 1; i < lM; ++i) gel(sg1,i) = quotient_subgroup_lift(C,H,gel(M,i));
1246 : /*sg2 is a list of lists of subgroups not intersecting with H*/
1247 10017 : sg2 = cgetg(lM, t_VEC);
1248 : /* Loop over all subgroups of G/H */
1249 59283 : for (j = 1; j < lM; ++j) gel(sg2,j) = liftsubgroup(C, H, gel(M,j));
1250 10017 : p1 = gconcat(sg1, shallowconcat1(sg2));
1251 10017 : if (sg3)
1252 : {
1253 42 : p1 = gconcat(p1, sg3);
1254 42 : if (n==5) /*ensure that the D4 subgroups of S4 are in supersolvable format*/
1255 84 : for(j = 3; j <= 5; j++)
1256 : {
1257 63 : GEN c = gmael(p1,j,1);
1258 63 : if (!perm_commute(gel(c,1),gel(c,3)))
1259 : {
1260 42 : if (perm_commute(gel(c,2),gel(c,3))) { swap(gel(c,1), gel(c,2)); }
1261 : else
1262 21 : perm_mul_inplace2(gel(c,2), gel(c,1));
1263 : }
1264 : }
1265 : }
1266 10017 : return gerepileupto(ltop,p1);
1267 : }
1268 :
1269 : /*return 1 if G is abelian, else 0*/
1270 : long
1271 8932 : group_isabelian(GEN G)
1272 : {
1273 8932 : GEN g = grp_get_gen(G);
1274 8932 : long i, j, n = lg(g);
1275 12852 : for(i=2; i<n; i++)
1276 13034 : for(j=1; j<i; j++)
1277 9114 : if (!perm_commute(gel(g,i), gel(g,j))) return 0;
1278 3990 : return 1;
1279 : }
1280 :
1281 : /*If G is abelian, return its HNF matrix*/
1282 : GEN
1283 385 : group_abelianHNF(GEN G, GEN S)
1284 : {
1285 385 : GEN M, g = grp_get_gen(G), o = grp_get_ord(G);
1286 385 : long i, j, k, n = lg(g);
1287 385 : if (!group_isabelian(G)) return NULL;
1288 315 : if (n==1) return cgetg(1,t_MAT);
1289 301 : if (!S) S = group_elts(G, group_domain(G));
1290 301 : M = cgetg(n,t_MAT);
1291 980 : for(i=1; i<n; i++)
1292 : {
1293 679 : GEN P, C = cgetg(n,t_COL);
1294 679 : pari_sp av = avma;
1295 679 : gel(M,i) = C;
1296 679 : P = perm_inv(perm_powu(gel(g,i), o[i]));
1297 959 : for(j=1; j<lg(S); j++)
1298 959 : if (zv_equal(P, gel(S,j))) break;
1299 679 : set_avma(av);
1300 679 : if (j==lg(S)) pari_err_BUG("galoisisabelian [inconsistent group]");
1301 679 : j--;
1302 1218 : for(k=1; k<i; k++)
1303 : {
1304 539 : long q = j / o[k];
1305 539 : gel(C,k) = stoi(j - q*o[k]);
1306 539 : j = q;
1307 : }
1308 679 : gel(C,k) = stoi(o[i]);
1309 1218 : for (k++; k<n; k++) gel(C,k) = gen_0;
1310 : }
1311 301 : return M;
1312 : }
1313 :
1314 : /*If G is abelian, return its abstract SNF matrix*/
1315 : GEN
1316 336 : group_abelianSNF(GEN G, GEN L)
1317 : {
1318 336 : pari_sp ltop = avma;
1319 336 : GEN H = group_abelianHNF(G,L);
1320 336 : if (!H) return NULL;
1321 266 : return gerepileupto(ltop, smithclean( ZM_snf(H) ));
1322 : }
1323 :
1324 : GEN
1325 434 : abelian_group(GEN v)
1326 : {
1327 434 : long card = zv_prod(v), i, d = 1, l = lg(v);
1328 434 : GEN G = cgetg(3,t_VEC), gen = cgetg(l,t_VEC);
1329 434 : gel(G,1) = gen;
1330 434 : gel(G,2) = vecsmall_copy(v);
1331 882 : for(i=1; i<l; i++)
1332 : {
1333 448 : GEN p = cgetg(card+1, t_VECSMALL);
1334 448 : long o = v[i], u = d*(o-1), j, k, l;
1335 448 : gel(gen, i) = p;
1336 : /* The following loop is over-optimized. Remember that I wrote it for
1337 : * testpermutation. Something has survived... BA */
1338 1036 : for(j=1;j<=card;)
1339 : {
1340 2296 : for(k=1;k<o;k++)
1341 4543 : for(l=1;l<=d; l++,j++) p[j] = j+d;
1342 1995 : for (l=1; l<=d; l++,j++) p[j] = j-u;
1343 : }
1344 448 : d += u;
1345 : }
1346 434 : return G;
1347 : }
1348 :
1349 : static long
1350 14609 : groupelts_subgroup_isnormal(GEN G, GEN H)
1351 : {
1352 14609 : long i, n = lg(G);
1353 64470 : for(i = 1; i < n; i++)
1354 62895 : if (!group_perm_normalize(H, gel(G,i))) return 0;
1355 1575 : return 1;
1356 : }
1357 :
1358 : /*return 1 if H is a normal subgroup of G*/
1359 : long
1360 336 : group_subgroup_isnormal(GEN G, GEN H)
1361 : {
1362 336 : if (lg(grp_get_gen(H)) > 1 && group_domain(G) != group_domain(H))
1363 0 : pari_err_DOMAIN("group_subgroup_isnormal","domain(H)","!=",
1364 : strtoGENstr("domain(G)"), H);
1365 336 : return groupelts_subgroup_isnormal(grp_get_gen(G), H);
1366 : }
1367 :
1368 : static GEN
1369 4816 : group_subgroup_kernel_set(GEN G, GEN H)
1370 : {
1371 : pari_sp av;
1372 4816 : GEN g = grp_get_gen(G);
1373 4816 : long i, n = lg(g);
1374 : GEN S, elts;
1375 4816 : long d = group_domain(G);
1376 4816 : if (lg(grp_get_gen(H)) > 1 && group_domain(G) != group_domain(H))
1377 0 : pari_err_DOMAIN("group_subgroup_isnormal","domain(H)","!=",
1378 : strtoGENstr("domain(G)"), H);
1379 4816 : elts = group_elts(H,d);
1380 4816 : S = groupelts_set(elts, d);
1381 4816 : av = avma;
1382 19264 : for(i=1; i<n; i++)
1383 : {
1384 14448 : F2v_and_inplace(S, groupelts_conj_set(elts,gel(g,i)));
1385 14448 : set_avma(av);
1386 : }
1387 4816 : return S;
1388 : }
1389 :
1390 : int
1391 4816 : group_subgroup_is_faithful(GEN G, GEN H)
1392 : {
1393 4816 : pari_sp av = avma;
1394 4816 : GEN K = group_subgroup_kernel_set(G,H);
1395 4816 : F2v_clear(K,1);
1396 4816 : return gc_long(av, F2v_equal0(K));
1397 : }
1398 :
1399 : long
1400 0 : groupelts_exponent(GEN elts)
1401 : {
1402 0 : long i, n = lg(elts)-1, expo = 1;
1403 0 : for(i=1; i<=n; i++) expo = ulcm(expo, perm_orderu(gel(elts,i)));
1404 0 : return expo;
1405 : }
1406 :
1407 : GEN
1408 700 : groupelts_center(GEN S)
1409 : {
1410 700 : pari_sp ltop = avma;
1411 700 : long i, j, n = lg(S)-1, l = n;
1412 700 : GEN V, elts = zero_F2v(n+1);
1413 25732 : for(i=1; i<=n; i++)
1414 : {
1415 25032 : if (F2v_coeff(elts,i)) { l--; continue; }
1416 573384 : for(j=1; j<=n; j++)
1417 563192 : if (!perm_commute(gel(S,i),gel(S,j)))
1418 : {
1419 14322 : F2v_set(elts,i);
1420 14322 : F2v_set(elts,j); l--; break;
1421 : }
1422 : }
1423 700 : V = cgetg(l+1,t_VEC);
1424 25732 : for (i=1, j=1; i<=n ;i++)
1425 25032 : if (!F2v_coeff(elts,i)) gel(V,j++) = vecsmall_copy(gel(S,i));
1426 700 : return gerepileupto(ltop,V);
1427 : }
1428 :
1429 : GEN
1430 4270 : groupelts_conjclasses(GEN elts, long *pnbcl)
1431 : {
1432 4270 : long i, j, cl = 0, n = lg(elts)-1;
1433 4270 : GEN c = const_vecsmall(n,0);
1434 4270 : pari_sp av = avma;
1435 52850 : for (i=1; i<=n; i++)
1436 : {
1437 48580 : GEN g = gel(elts,i);
1438 48580 : if (c[i]) continue;
1439 34965 : c[i] = ++cl;
1440 486871 : for(j=1; j<=n; j++)
1441 451906 : if (j != i)
1442 : {
1443 416941 : GEN h = perm_conj(gel(elts,j), g);
1444 416941 : long i2 = gen_search(elts,h,(void*)&vecsmall_lexcmp,&cmp_nodata);
1445 416941 : c[i2] = cl; set_avma(av);
1446 : }
1447 : }
1448 4270 : if (pnbcl) *pnbcl = cl;
1449 4270 : return c;
1450 : }
1451 :
1452 : GEN
1453 4270 : conjclasses_repr(GEN conj, long nb)
1454 : {
1455 4270 : long i, l = lg(conj);
1456 4270 : GEN e = const_vecsmall(nb, 0);
1457 52850 : for(i=1; i<l; i++)
1458 : {
1459 48580 : long ci = conj[i];
1460 48580 : if (!e[ci]) e[ci] = i;
1461 : }
1462 4270 : return e;
1463 : }
1464 :
1465 : /* elts of G sorted wrt vecsmall_lexcmp order: g in G is determined by g[1]
1466 : * so sort by increasing g[1] */
1467 : static GEN
1468 3885 : galois_elts_sorted(GEN gal)
1469 : {
1470 : long i, l;
1471 3885 : GEN elts = gal_get_group(gal), v = cgetg_copy(elts, &l);
1472 43141 : for (i = 1; i < l; i++) { GEN g = gel(elts,i); gel(v, g[1]) = g; }
1473 3885 : return v;
1474 : }
1475 : GEN
1476 4291 : group_to_cc(GEN G)
1477 : {
1478 4291 : GEN elts = checkgroupelts(G), z = cgetg(5,t_VEC);
1479 4270 : long n, flag = 1;
1480 4270 : if (typ(gel(G,1)) == t_POL)
1481 3885 : elts = galois_elts_sorted(G); /* galoisinit */
1482 : else
1483 : {
1484 385 : long i, l = lg(elts);
1485 385 : elts = gen_sort_shallow(elts,(void*)vecsmall_lexcmp,cmp_nodata);
1486 5824 : for (i = 1; i < l; i++)
1487 5586 : if (gel(elts,i)[1] != i) { flag = 0; break; }
1488 : }
1489 4270 : gel(z,1) = elts;
1490 4270 : gel(z,2) = groupelts_conjclasses(elts,&n);
1491 4270 : gel(z,3) = conjclasses_repr(gel(z,2),n);
1492 4270 : gel(z,4) = utoi(flag); return z;
1493 : }
1494 :
1495 : /* S a list of generators */
1496 : GEN
1497 0 : groupelts_abelian_group(GEN S)
1498 : {
1499 0 : pari_sp ltop = avma;
1500 : GEN Qgen, Qord, Qelt;
1501 0 : long i, j, n = lg(gel(S,1))-1, l = lg(S);
1502 0 : Qord = cgetg(l, t_VECSMALL);
1503 0 : Qgen = cgetg(l, t_VEC);
1504 0 : Qelt = mkvec(identity_perm(n));
1505 0 : for (i = 1, j = 1; i < l; ++i)
1506 : {
1507 0 : GEN g = gel(S,i);
1508 0 : long o = perm_relorder(g, groupelts_set(Qelt, n));
1509 0 : gel(Qgen,j) = g;
1510 0 : Qord[j] = o;
1511 0 : if (o != 1) { Qelt = perm_generate(g, Qelt, o); j++; }
1512 : }
1513 0 : setlg(Qgen,j);
1514 0 : setlg(Qord,j);
1515 0 : return gerepilecopy(ltop, mkvec2(Qgen, Qord));
1516 : }
1517 :
1518 : GEN
1519 14 : group_export_GAP(GEN G)
1520 : {
1521 14 : pari_sp av = avma;
1522 14 : GEN s, comma, g = grp_get_gen(G);
1523 14 : long i, k, l = lg(g);
1524 14 : if (l == 1) return strtoGENstr("Group(())");
1525 7 : s = cgetg(2*l, t_VEC);
1526 7 : comma = strtoGENstr(", ");
1527 7 : gel(s,1) = strtoGENstr("Group(");
1528 28 : for (i=1, k=2; i < l; ++i)
1529 : {
1530 21 : if (i > 1) gel(s,k++) = comma;
1531 21 : gel(s,k++) = perm_to_GAP(gel(g,i));
1532 : }
1533 7 : gel(s,k++) = strtoGENstr(")");
1534 7 : return gerepilecopy(av, shallowconcat1(s));
1535 : }
1536 :
1537 : GEN
1538 14 : group_export_MAGMA(GEN G)
1539 : {
1540 14 : pari_sp av = avma;
1541 14 : GEN s, comma, g = grp_get_gen(G);
1542 14 : long i, k, l = lg(g);
1543 14 : if (l == 1) return strtoGENstr("PermutationGroup<1|>");
1544 7 : s = cgetg(2*l, t_VEC);
1545 7 : comma = strtoGENstr(", ");
1546 7 : gel(s,1) = gsprintf("PermutationGroup<%ld|",group_domain(G));
1547 28 : for (i=1, k=2; i < l; ++i)
1548 : {
1549 21 : if (i > 1) gel(s,k++) = comma;
1550 21 : gel(s,k++) = GENtoGENstr( vecsmall_to_vec(gel(g,i)) );
1551 : }
1552 7 : gel(s,k++) = strtoGENstr(">");
1553 7 : return gerepilecopy(av, shallowconcat1(s));
1554 : }
1555 :
1556 : GEN
1557 28 : group_export(GEN G, long format)
1558 : {
1559 28 : switch(format)
1560 : {
1561 14 : case 0: return group_export_GAP(G);
1562 14 : case 1: return group_export_MAGMA(G);
1563 : }
1564 0 : pari_err_FLAG("galoisexport");
1565 0 : return NULL; /*-Wall*/
1566 : }
1567 :
1568 : static GEN
1569 3577 : groupelts_cyclic_subgroups(GEN G)
1570 : {
1571 3577 : pari_sp av = avma;
1572 3577 : long i, j, n = lg(G)-1;
1573 : GEN elts, f, gen, ord;
1574 3577 : if (n==1) return cgetg(1,t_VEC);
1575 3577 : elts = zero_F2v(lg(gel(G,1))-1);
1576 3577 : gen = cgetg(n+1, t_VECSMALL);
1577 3577 : ord = cgetg(n+1, t_VECSMALL);
1578 53109 : for (i=1, j=1; i<=n; i++)
1579 : {
1580 49532 : long k = 1, o, c = 0;
1581 49532 : GEN p = gel(G, i);
1582 49532 : if (F2v_coeff(elts, p[1])) continue;
1583 35903 : o = perm_orderu(p);
1584 35903 : gen[j] = i; ord[j] = o; j++;
1585 : do
1586 : {
1587 96229 : if (cgcd(o, ++c)==1) F2v_set(elts, p[k]);
1588 96229 : k = p[k];
1589 96229 : } while (k!=1);
1590 : }
1591 3577 : setlg(gen, j);
1592 3577 : setlg(ord, j);
1593 3577 : f = vecsmall_indexsort(ord);
1594 3577 : return gerepilecopy(av, mkvec2(vecsmallpermute(gen, f),
1595 : vecsmallpermute(ord, f)));
1596 : }
1597 :
1598 : GEN
1599 3584 : groupelts_to_group(GEN G)
1600 : {
1601 3584 : pari_sp av = avma;
1602 : GEN L, cyc, ord;
1603 3584 : long i, l, n = lg(G)-1;
1604 3584 : if (n==1) return trivialgroup();
1605 3563 : L = groupelts_cyclic_subgroups(G);
1606 3563 : cyc = gel(L,1); ord = gel(L,2);
1607 3563 : l = lg(cyc);
1608 16324 : for (i = l-1; i >= 2; i--)
1609 : {
1610 15645 : GEN p = gel(G,cyc[i]);
1611 15645 : long o = ord[i];
1612 15645 : GEN H = cyclicgroup(p, o);
1613 15645 : if (o == n) return gerepileupto(av, H);
1614 14273 : if (groupelts_subgroup_isnormal(G, H))
1615 : {
1616 1512 : GEN C = groupelts_quotient(G, H);
1617 1512 : GEN Q = quotient_groupelts(C);
1618 1512 : GEN R = groupelts_to_group(Q);
1619 1512 : if (!R) return gc_NULL(av);
1620 1512 : return gerepilecopy(av, quotient_subgroup_lift(C, H, R));
1621 : }
1622 : }
1623 679 : if (n==12 && l==9 && ord[2]==2 && ord[3]==2 && ord[5]==3)
1624 602 : return gerepilecopy(av,
1625 301 : mkvec2(mkvec3(gel(G,cyc[2]), gel(G,cyc[3]), gel(G,cyc[5])), mkvecsmall3(2,2,3)));
1626 378 : if (n==24 && l==18 && ord[11]==3 && ord[15]==4 && ord[16]==4)
1627 : {
1628 350 : GEN t21 = perm_sqr(gel(G,cyc[15]));
1629 350 : GEN t22 = perm_sqr(gel(G,cyc[16]));
1630 350 : GEN s = perm_mul(t22, gel(G,cyc[15]));
1631 700 : return gerepilecopy(av,
1632 350 : mkvec2(mkvec4(t21,t22, gel(G,cyc[11]), s), mkvecsmall4(2,2,3,2)));
1633 : }
1634 28 : if (n==36 && l==24 && ord[11]==3 && ord[15]==4)
1635 : {
1636 7 : GEN t1 = gel(G,cyc[11]), t3 = gel(G,cyc[15]);
1637 7 : return gerepilecopy(av,
1638 : mkvec2(mkvec3(perm_conj(t3, t1), t1, t3), mkvecsmall3(3,3,4)));
1639 : }
1640 21 : return gc_NULL(av);
1641 : }
1642 :
1643 : static GEN
1644 1176 : subg_get_gen(GEN subg) { return gel(subg, 1); }
1645 :
1646 : static GEN
1647 8694 : subg_get_set(GEN subg) { return gel(subg, 2); }
1648 :
1649 : static GEN
1650 812 : groupelt_subg_normalize(GEN elt, GEN subg, GEN cyc)
1651 : {
1652 812 : GEN gen = subg_get_gen(subg), set = subg_get_set(subg);
1653 812 : long i, j, u, n = lg(elt)-1, lgen = lg(gen);
1654 812 : GEN b = F2v_copy(cyc), res = zero_F2v(n);
1655 49532 : for(i = 1; i <= n; i++)
1656 : {
1657 : GEN g;
1658 48720 : if (!F2v_coeff(b, i)) continue;
1659 22386 : g = gel(elt,i);
1660 763532 : for(u=1; u<=n; u++)
1661 763532 : if (g[u]==1) break;
1662 25172 : for(j=1; j<lgen; j++)
1663 : {
1664 23786 : GEN h = gel(elt,gen[j]);
1665 23786 : if (!F2v_coeff(set,g[h[u]])) break;
1666 : }
1667 22386 : if (j < lgen) continue;
1668 1386 : F2v_set(res,i);
1669 84546 : for(j=1; j <= n; j++)
1670 83160 : if (F2v_coeff(set, j))
1671 6720 : F2v_clear(b,g[gel(elt,j)[1]]);
1672 : }
1673 812 : return res;
1674 : }
1675 :
1676 : static GEN
1677 14 : triv_subg(GEN elt)
1678 : {
1679 14 : GEN v = cgetg(3, t_VEC);
1680 14 : gel(v,1) = cgetg(1,t_VECSMALL);
1681 14 : gel(v,2) = zero_F2v(lg(elt)-1);
1682 14 : F2v_set(gel(v,2),1);
1683 14 : return v;
1684 : }
1685 :
1686 : static GEN
1687 364 : subg_extend(GEN U, long e, long o, GEN elt)
1688 : {
1689 364 : long i, j, n = lg(elt)-1;
1690 364 : GEN g = gel(elt, e);
1691 364 : GEN gen = vecsmall_append(subg_get_gen(U), e);
1692 364 : GEN set = subg_get_set(U);
1693 364 : GEN Vset = zv_copy(set);
1694 22204 : for(i = 1; i <= n; i++)
1695 21840 : if (F2v_coeff(set, i))
1696 : {
1697 1260 : long h = gel(elt, i)[1];
1698 2800 : for(j = 1; j < o; j++)
1699 : {
1700 1540 : h = g[h];
1701 1540 : F2v_set(Vset, h);
1702 : }
1703 : }
1704 364 : return mkvec2(gen, Vset);
1705 : }
1706 :
1707 : static GEN
1708 434 : cyclic_subg(long e, long o, GEN elt)
1709 : {
1710 434 : long j, n = lg(elt)-1, h = 1;
1711 434 : GEN g = gel(elt, e);
1712 434 : GEN gen = mkvecsmall(e);
1713 434 : GEN set = zero_F2v(n);
1714 434 : F2v_set(set,1);
1715 1260 : for(j = 1; j < o; j++)
1716 : {
1717 826 : h = g[h];
1718 826 : F2v_set(set, h);
1719 : }
1720 434 : return mkvec2(gen, set);
1721 : }
1722 :
1723 : static GEN
1724 14 : groupelts_to_regular(GEN elt)
1725 : {
1726 14 : long i, j, n = lg(elt)-1;
1727 14 : GEN V = cgetg(n+1,t_VEC);
1728 854 : for (i=1; i<=n; i++)
1729 : {
1730 840 : pari_sp av = avma;
1731 840 : GEN g = gel(elt, i);
1732 840 : GEN W = cgetg(n+1,t_VEC);
1733 51240 : for(j=1; j<=n; j++)
1734 50400 : gel(W,j) = perm_mul(g, gel(elt,j));
1735 840 : gel(V, i) = gerepileuptoleaf(av,vecvecsmall_indexsort(W));
1736 : }
1737 14 : vecvecsmall_sort_inplace(V, NULL);
1738 14 : return V;
1739 : }
1740 :
1741 : static long
1742 434 : groupelts_pow(GEN elt, long j, long n)
1743 : {
1744 434 : GEN g = gel(elt,j);
1745 434 : long i, h = 1;
1746 1694 : for (i=1; i<=n; i++)
1747 1260 : h = g[h];
1748 434 : return h;
1749 : }
1750 :
1751 : static GEN
1752 14 : groupelts_cyclic_primepow(GEN elt, GEN *pt_pr, GEN *pt_po)
1753 : {
1754 14 : GEN R = groupelts_cyclic_subgroups(elt);
1755 14 : GEN gen = gel(R,1), ord = gel(R,2);
1756 14 : long i, n = lg(elt)-1, l = lg(gen);
1757 14 : GEN set = zero_F2v(n);
1758 14 : GEN pr = zero_Flv(n);
1759 14 : GEN po = zero_Flv(n);
1760 462 : for (i = 1; i < l; i++)
1761 : {
1762 448 : long h = gen[i];
1763 : ulong p;
1764 448 : if (uisprimepower(ord[i], &p))
1765 : {
1766 434 : F2v_set(set, h);
1767 434 : uel(pr,h) = p;
1768 434 : po[h] = groupelts_pow(elt, h, p);
1769 : }
1770 : }
1771 14 : *pt_pr = pr; *pt_po = po;
1772 14 : return set;
1773 : }
1774 :
1775 : static GEN
1776 50400 : perm_bracket(GEN p, GEN q)
1777 : {
1778 50400 : return perm_mul(perm_mul(p,q), perm_inv(perm_mul(q,p)));
1779 : }
1780 :
1781 : static GEN
1782 826 : set_groupelts(GEN S, GEN x)
1783 : {
1784 826 : long i, n = F2v_hamming(x), k=1, m = x[1];
1785 826 : GEN v = cgetg(n+1, t_VEC);
1786 50386 : for (i=1; i<=m; i++)
1787 49560 : if (F2v_coeff(x,i))
1788 4914 : gel(v,k++) = gel(S,i);
1789 826 : return v;
1790 : }
1791 :
1792 : static GEN
1793 14 : set_idx(GEN x)
1794 : {
1795 14 : long i, n = F2v_hamming(x), k=1, m = x[1];
1796 14 : GEN v = cgetg(n+1, t_VECSMALL);
1797 854 : for (i=1; i<=m; i++)
1798 840 : if (F2v_coeff(x,i))
1799 840 : uel(v,k++) = i;
1800 14 : return v;
1801 : }
1802 :
1803 : static GEN
1804 14 : set_derived(GEN set, GEN elts)
1805 : {
1806 14 : long i, j, l = lg(elts);
1807 14 : GEN V = zero_F2v(l-1);
1808 854 : for(i = 1; i < l; i++)
1809 840 : if (F2v_coeff(set, i))
1810 51240 : for(j = 1; j < l; j++)
1811 50400 : if (F2v_coeff(set, j))
1812 50400 : F2v_set(V, perm_bracket(gel(elts,i),gel(elts,j))[1]);
1813 14 : return V;
1814 : }
1815 :
1816 : static GEN
1817 14 : groupelts_residuum(GEN elts)
1818 : {
1819 14 : pari_sp av = avma;
1820 14 : long o = lg(elts)-1, oo;
1821 14 : GEN set = const_F2v(o);
1822 : do
1823 : {
1824 14 : oo = o;
1825 14 : set = set_derived(set, elts);
1826 14 : o = F2v_hamming(set);
1827 14 : } while (o > 1 && o < oo);
1828 14 : if (o==1) return NULL;
1829 14 : return gerepilecopy(av,mkvec2(set_idx(set), set));
1830 : }
1831 :
1832 : static GEN
1833 14 : all_cyclic_subg(GEN pr, GEN po, GEN elt)
1834 : {
1835 14 : long i, n = lg(pr)-1, m = 0, k = 1;
1836 : GEN W;
1837 854 : for (i=1; i <= n; i++)
1838 840 : m += po[i]==1;
1839 14 : W = cgetg(m+1, t_VEC);
1840 854 : for (i=1; i <= n; i++)
1841 840 : if (po[i]==1)
1842 434 : gel(W, k++) = cyclic_subg(i, pr[i], elt);
1843 14 : return W;
1844 : }
1845 :
1846 : static GEN
1847 14 : groupelts_subgroups_raw(GEN elts)
1848 : {
1849 14 : pari_sp av = avma;
1850 14 : GEN elt = groupelts_to_regular(elts);
1851 14 : GEN pr, po, cyc = groupelts_cyclic_primepow(elt, &pr, &po);
1852 14 : long n = lg(elt)-1;
1853 14 : long i, j, nS = 1;
1854 14 : GEN S, L, R = NULL;
1855 14 : S = cgetg(1+bigomegau(n)+1, t_VEC);
1856 14 : gel(S, nS++) = mkvec(triv_subg(elt));
1857 14 : gel(S, nS++) = L = all_cyclic_subg(pr, po, elt);
1858 14 : if (DEBUGLEVEL) err_printf("subgroups: level %ld: %ld\n",nS-1,lg(L)-1);
1859 70 : while (lg(L) > 1)
1860 : {
1861 56 : pari_sp av2 = avma;
1862 56 : long nW = 1, lL = lg(L);
1863 56 : long ng = n;
1864 56 : GEN W = cgetg(1+ng, t_VEC);
1865 868 : for (i=1; i<lL; i++)
1866 : {
1867 812 : GEN U = gel(L, i), set = subg_get_set(U);
1868 812 : GEN G = groupelt_subg_normalize(elt, U, cyc);
1869 7154 : for (j=1; j<nW; j++)
1870 : {
1871 6342 : GEN Wj = subg_get_set(gel(W, j));
1872 6342 : if (F2v_subset(set, Wj))
1873 728 : F2v_negimply_inplace(G, Wj);
1874 : }
1875 49532 : for (j=1; j<=n; j++)
1876 48720 : if(F2v_coeff(G,j))
1877 : {
1878 770 : long p = pr[j];
1879 770 : if (F2v_coeff(set, j)) continue;
1880 364 : if (F2v_coeff(set, po[j]))
1881 : {
1882 364 : GEN U2 = subg_extend(U, j, p, elt);
1883 364 : F2v_negimply_inplace(G, subg_get_set(U2));
1884 364 : if (nW > ng) { ng<<=1; W = vec_lengthen(W, ng); }
1885 364 : gel(W, nW++) = U2;
1886 : }
1887 : }
1888 : }
1889 56 : setlg(W, nW);
1890 56 : L = W;
1891 56 : if (nW > 1) gel(S, nS++) = L = gerepilecopy(av2, W);
1892 56 : if (DEBUGLEVEL) err_printf("subgroups: level %ld: %ld\n",nS-1,nW-1);
1893 56 : if (lg(L)==1 && !R)
1894 : {
1895 14 : R = groupelts_residuum(elt);
1896 14 : if (!R) break;
1897 14 : gel(S, nS++) = L = mkvec(R);
1898 : }
1899 : }
1900 14 : setlg(S, nS);
1901 14 : return gerepilecopy(av, shallowconcat1(S));
1902 : }
1903 :
1904 : static GEN
1905 14 : subg_to_elts(GEN S, GEN x)
1906 840 : { pari_APPLY_type(t_VEC, set_groupelts(S, gmael(x,i,2))); }
1907 :
1908 : GEN
1909 14 : groupelts_solvablesubgroups(GEN G)
1910 : {
1911 14 : pari_sp av = avma;
1912 14 : GEN S = vecvecsmall_sort(checkgroupelts(G));
1913 14 : GEN L = groupelts_subgroups_raw(S);
1914 14 : return gerepilecopy(av, subg_to_elts(S, L));
1915 : }
|
