Line data Source code
1 : /* Copyright (C) 2000 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 : #define DEBUGLEVEL DEBUGLEVEL_subgrouplist
19 :
20 : typedef struct slist {
21 : struct slist *next;
22 : long *data;
23 : long prec;
24 : } slist;
25 :
26 : typedef struct {
27 : GEN cyc, gen;
28 : ulong count;
29 : slist *list;
30 : } sublist_t;
31 :
32 : /* SUBGROUPS
33 : * G = Gp x I, with Gp a p-Sylow (I assumed small).
34 : * Compute subgroups of I by recursive calls
35 : * Loop through subgroups Hp of Gp using Birkhoff's algorithm.
36 : * If (I is non trivial)
37 : * lift Hp to G (mul by exponent of I)
38 : * for each subgp of I, lift it to G (mult by exponent of Gp)
39 : * consider the group generated by the two subgroups (concat)
40 : *
41 : * type(H) = mu --> H = Z/p^mu[1] x ... x Z/p^mu[len(mu)] */
42 : typedef struct subgp_iter {
43 : long *M, *L; /* mu = p-subgroup type, lambda = p-group type */
44 : GEN *powlist; /* [i] = p^i, i = 0.. */
45 : long *c, *maxc;
46 : GEN *a, *maxa, **g, **maxg;
47 : long *available;
48 : GEN **H; /* p-subgroup of type mu, in matrix form */
49 : GEN cyc; /* cyclic factors of G */
50 : GEN subq;/* subgrouplist(I) */
51 : GEN subqpart; /* J in subq s.t [I:J][Gp:Hp] <= indexbound */
52 : GEN bound; /* if != NULL, impose a "bound" on [G:H] (see boundtype) */
53 : long boundtype;
54 : long countsub; /* number of subgroups of type M (so far) */
55 : long count; /* number of p-subgroups so far [updated when M completed] */
56 : GEN expoI; /* exponent of I */
57 :
58 : long(*fun)(void*, GEN); /* callback applied to each subgroup */
59 : void *fundata; /* data for fun */
60 : long stop;
61 : } subgp_iter;
62 :
63 : /* MAX: [G:H] <= bound, EXACT: [G:H] = bound, TYPE: type(H) = bound */
64 : enum { b_NONE, b_MAX, b_EXACT, b_TYPE };
65 :
66 : #define len(x) (x)[0]
67 : #define setlen(x,l) len(x)=(l)
68 :
69 : static void
70 35 : printtyp(const long *typ) /*Used only for ddebugging */
71 : {
72 35 : long i, l = len(typ);
73 133 : for (i=1; i<=l; i++) err_printf(" %ld ",typ[i]);
74 35 : err_printf("\n");
75 35 : }
76 :
77 : /* compute conjugate partition of typ */
78 : static long*
79 70 : conjugate(long *typ)
80 : {
81 70 : long *t, i, k = len(typ), l, last;
82 :
83 70 : if (!k) { t = new_chunk(1); setlen(t,0); return t; }
84 63 : l = typ[1]; t = new_chunk(l+2);
85 63 : t[1] = k; last = k;
86 63 : for (i=2; i<=l; i++)
87 : {
88 0 : while (typ[last] < i) last--;
89 0 : t[i] = last;
90 : }
91 63 : t[i] = 0; setlen(t,l);
92 63 : return t;
93 : }
94 : /* ----subgp_iter 'fun' attached to subgrouplist ------------- */
95 : static void
96 48200 : addcell(sublist_t *S, GEN H)
97 : {
98 48200 : long *pt,i,j,L, n = lg(H)-1;
99 : slist *cell;
100 :
101 48200 : L = 3;
102 118575 : for (j=1; j<=n; j++)
103 : { /* H in HNF, largest entries are on diagonal */
104 70375 : long l = lgefint(gcoeff(H,j,j));
105 70375 : if (l > L) L = l;
106 : }
107 48200 : L -= 2;
108 48200 : cell = (slist*) pari_malloc(sizeof(slist)
109 48200 : + ((n*(n+1)) >> 1) * sizeof(long) * L);
110 48202 : S->list->next = cell; cell->data = pt = (long*) (cell + 1);
111 48202 : cell->prec = L;
112 118580 : for (j=1; j<=n; j++)
113 171717 : for(i=1; i<=j; i++) {
114 101339 : GEN z = gcoeff(H,i,j);
115 101339 : long h, lz = lgefint(z) - 2;
116 119750 : for (h = 0; h < L - lz; h++) *pt++ = 0;
117 184267 : for (h = 0; h < lz; h++) *pt++ = z[h+2];
118 : }
119 48202 : S->list = cell;
120 48202 : S->count++;
121 48202 : }
122 :
123 : static long
124 101523 : list_fun(void *E, GEN x)
125 : {
126 101523 : sublist_t *S = (sublist_t*)E;
127 101523 : GEN H = ZM_hnfmodid(x, S->cyc);
128 101526 : if (!S->gen || subgroup_conductor_ok(H, S->gen)) addcell(S, H);
129 101523 : return 0;
130 : }
131 : /* -------------------------------------------------------------- */
132 :
133 : /* treat subgroup Hp (not in HNF, T->fun should do it if desired) */
134 : static void
135 99857 : treatsub(subgp_iter *T, GEN H)
136 : {
137 : long i;
138 99857 : if (!T->subq) {T->stop = T->fun(T->fundata, H); T->countsub++; }
139 : else
140 : { /* not a p group, add the trivial part */
141 756 : GEN Hp = gmul(T->expoI, H); /* lift H to G */
142 756 : long n = lg(T->subqpart)-1;
143 1834 : for (i=1; i<=n; i++)
144 1078 : if (T->fun(T->fundata, shallowconcat(Hp, gel(T->subqpart,i))))
145 0 : { T->stop = 1; break; }
146 756 : T->countsub += n;
147 : }
148 99857 : }
149 :
150 : /* x a t_INT, x++. Could be optimized... */
151 : static void
152 118603 : inc(GEN x) { affii(addiu(x,1), x); }
153 :
154 : /* assume t>0 and l>1 */
155 : static void
156 39186 : dogroup(subgp_iter *T)
157 : {
158 39186 : const GEN *powlist = T->powlist;
159 39186 : long *M = T->M;
160 39186 : long *L = T->L;
161 39186 : long *c = T->c;
162 39186 : GEN *a = T->a, *maxa = T->maxa;
163 39186 : GEN **g = T->g, **maxg = T->maxg;
164 39186 : GEN **H = T->H;
165 : pari_sp av;
166 39186 : long i,j,k,r,n,t2,ind, t = len(M), l = len(L);
167 :
168 39186 : t2 = (l==t)? t-1: t;
169 39186 : n = t2 * l - (t2*(t2+1))/2; /* number of gamma_ij */
170 64848 : for (i=1, r=t+1; ; i++)
171 : {
172 64848 : if (T->available[i]) c[r++] = i;
173 64848 : if (r > l) break;
174 : }
175 39186 : if (DEBUGLEVEL>6) { err_printf(" column selection:"); printtyp(c); }
176 : /* a/g and maxa/maxg access the same data indexed differently */
177 90713 : for (ind=0,i=1; i<=t; ind+=(l-i),i++)
178 : {
179 51527 : maxg[i] = maxa + (ind - (i+1)); /* only access maxg[i][i+1..l] */
180 51527 : g[i] = a + (ind - (i+1));
181 116767 : for (r=i+1; r<=l; r++)
182 65240 : if (c[r] < c[i]) maxg[i][r] = powlist[M[i]-M[r]-1];
183 39753 : else if (L[c[r]] < M[i]) maxg[i][r] = powlist[L[c[r]]-M[r]];
184 39438 : else maxg[i][r] = powlist[M[i]-M[r]];
185 : }
186 : /* allocate correct lg */
187 104426 : for (i = 0; i<= n-1; i++) a[i] = icopy(maxa[i]);
188 65240 : affui(0, a[n-1]); for (i=0; i<n-1; i++) affui(1, a[i]);
189 39186 : av = avma;
190 109635 : for(;!T->stop;)
191 : {
192 109635 : inc(a[n-1]);
193 109625 : if (cmpii(a[n-1], maxa[n-1]) > 0)
194 : {
195 76970 : j=n-2; while (j>=0 && equalii(a[j], maxa[j])) j--;
196 48151 : if (j < 0) return;
197 :
198 8967 : inc(a[j]);
199 20699 : for (k=j+1; k<n; k++) affsi(1, a[k]);
200 : }
201 174818 : for (i=1; i<=t; i++)
202 : {
203 145154 : for (r=1; r<i; r++) H[i][c[r]] = gen_0;
204 104372 : H[i][c[r]] = powlist[L[c[r]] - M[r]];
205 249553 : for (r=i+1; r<=l; r++)
206 : {
207 145185 : GEN e = g[i][r];
208 145185 : long d = L[c[r]] - M[i];
209 145185 : if (c[r] < c[i])
210 30828 : e = mulii(e, powlist[d+1]);
211 114357 : else if (d > 0)
212 91 : e = mulii(e, powlist[d]);
213 145181 : H[i][c[r]] = e;
214 : }
215 : }
216 70446 : treatsub(T, (GEN)H); set_avma(av);
217 : }
218 : }
219 :
220 : /* T->c[1],...,T->c[r-1] filled */
221 : static void
222 70693 : loop(subgp_iter *T, long r)
223 : {
224 : long j;
225 :
226 70693 : if (r > len(T->M)) {
227 39186 : pari_sp av = avma; dogroup(T); set_avma(av);
228 39185 : return;
229 : }
230 :
231 31507 : if (r!=1 && (T->M[r-1] == T->M[r])) j = T->c[r-1]+1; else j = 1;
232 84531 : for ( ; j<=T->maxc[r]; j++)
233 53024 : if (T->available[j])
234 : {
235 52807 : T->c[r] = j; T->available[j] = 0;
236 52807 : loop(T, r+1); T->available[j] = 1;
237 : }
238 : }
239 :
240 : static void
241 47292 : dopsubtyp(subgp_iter *T)
242 : {
243 47292 : pari_sp av = avma;
244 47292 : long i,r, l = len(T->L), t = len(T->M);
245 :
246 47292 : if (!t) { treatsub(T, mkmat( zerocol(l) )); set_avma(av); return; }
247 17934 : if (l==1) /* imply t = 1 */
248 : {
249 49 : GEN p1 = gtomat(T->powlist[T->L[1]-T->M[1]]);
250 49 : treatsub(T, p1); set_avma(av); return;
251 : }
252 17885 : T->c = new_chunk(l+1); setlen(T->c, l);
253 17885 : T->maxc = new_chunk(l+1);
254 17885 : T->available = new_chunk(l+1);
255 17885 : T->a = (GEN*)new_chunk(l*(t+1));
256 17885 : T->maxa= (GEN*)new_chunk(l*(t+1));
257 17885 : T->g = (GEN**)new_chunk(t+1);
258 17885 : T->maxg = (GEN**)new_chunk(t+1);
259 :
260 17885 : if (DEBUGLEVEL>4) { err_printf(" subgroup:"); printtyp(T->M); }
261 39837 : for (i=1; i<=t; i++)
262 : {
263 74249 : for (r=1; r<=l; r++)
264 52493 : if (T->M[i] > T->L[r]) break;
265 21952 : T->maxc[i] = r-1;
266 : }
267 17885 : T->H = (GEN**)cgetg(t+1, t_MAT);
268 39837 : for (i=1; i<=t; i++) T->H[i] = (GEN*)cgetg(l+1, t_COL);
269 57505 : for (i=1; i<=l; i++) T->available[i]=1;
270 39837 : for (i=1; i<=t; i++) T->c[i]=0;
271 : /* go through all column selections */
272 17885 : loop(T, 1); set_avma(av); return;
273 : }
274 :
275 : static long
276 163934 : weight(long *typ)
277 : {
278 163934 : long i, l = len(typ), w = 0;
279 327524 : for (i=1; i<=l; i++) w += typ[i];
280 163934 : return w;
281 : }
282 :
283 : static void
284 47270 : dopsub(subgp_iter *T, GEN p, GEN indexsubq)
285 : {
286 47270 : long *M, *L = T->L;
287 47270 : long w,i,j,k,lsubq, wG = weight(L), wmin = 0, wmax = wG, n = len(L);
288 :
289 47270 : if (DEBUGLEVEL>4) { err_printf("\ngroup:"); printtyp(L); }
290 47270 : T->count = 0;
291 47270 : switch(T->boundtype)
292 : {
293 42 : case b_MAX: /* upper bound */
294 42 : wmin = (long) (wG - (log(gtodouble(T->bound)) / log(gtodouble(p))));
295 42 : if (cmpii(powiu(p, wG - wmin), T->bound) > 0) wmin++;
296 42 : break;
297 47088 : case b_EXACT: /* exact value */
298 47088 : wmin = wmax = wG - Z_pval(T->bound, p);
299 47089 : break;
300 : }
301 47271 : T->M = M = new_chunk(n+1);
302 47271 : if (T->subq)
303 : {
304 238 : lsubq = lg(T->subq);
305 238 : T->subqpart = T->bound? cgetg(lsubq, t_VEC): T->subq;
306 : }
307 : else
308 47033 : lsubq = 0; /* -Wall */
309 68509 : M[1] = -1; for (i=2; i<=n; i++) M[i]=0;
310 163939 : for(;!T->stop;) /* go through all vectors mu_{i+1} <= mu_i <= lam_i */
311 : {
312 163939 : M[1]++;
313 163939 : if (M[1] > L[1])
314 : {
315 93909 : for (j=2; j<=n; j++)
316 46639 : if (M[j] < L[j] && M[j] < M[j-1]) break;
317 68556 : if (j > n) return;
318 :
319 46687 : M[j]++; for (k=1; k<j; k++) M[k]=M[j];
320 : }
321 211754 : for (j=1; j<=n; j++)
322 163805 : if (!M[j]) break;
323 116668 : setlen(M, j-1); w = weight(M);
324 116668 : if (w >= wmin && w <= wmax)
325 : {
326 47292 : GEN p1 = gen_1;
327 :
328 47292 : if (T->subq && T->bound) /* G not a p-group */
329 : {
330 217 : pari_sp av = avma;
331 217 : GEN indexH = powiu(p, wG - w);
332 217 : GEN B = divii(T->bound, indexH);
333 217 : k = 1;
334 546 : for (i=1; i<lsubq; i++)
335 329 : if (cmpii(gel(indexsubq,i), B) <= 0)
336 273 : T->subqpart[k++] = T->subq[i];
337 217 : setlg(T->subqpart, k); set_avma(av);
338 : }
339 47292 : if (DEBUGLEVEL>4)
340 : {
341 35 : long *Lp = conjugate(L);
342 35 : long *Mp = conjugate(M);
343 35 : GEN BINMAT = matqpascal(len(L)+1, p);
344 :
345 35 : if (DEBUGLEVEL>7)
346 : {
347 0 : err_printf(" lambda = "); printtyp(L);
348 0 : err_printf(" lambda'= "); printtyp(Lp);
349 0 : err_printf(" mu = "); printtyp(M);
350 0 : err_printf(" mu'= "); printtyp(Mp);
351 : }
352 63 : for (j=1; j<=len(Mp); j++)
353 : {
354 28 : p1 = mulii(p1, powiu(p, Mp[j+1]*(Lp[j]-Mp[j])));
355 28 : p1 = mulii(p1, gcoeff(BINMAT, Lp[j]-Mp[j+1]+1, Mp[j]-Mp[j+1]+1));
356 : }
357 35 : err_printf(" alpha_lambda(mu,p) = %Ps\n",p1);
358 : }
359 :
360 47292 : T->countsub = 0; dopsubtyp(T);
361 47291 : T->count += T->countsub;
362 :
363 47291 : if (DEBUGLEVEL>4)
364 : {
365 35 : err_printf(" countsub = %ld\n", T->countsub);
366 35 : if (T->subq) p1 = muliu(p1,lg(T->subqpart)-1);
367 35 : if (!absequaliu(p1,T->countsub))
368 : {
369 0 : err_printf(" alpha = %Ps\n",p1);
370 0 : pari_err_BUG("forsubgroup (alpha != countsub)");
371 : }
372 : }
373 : }
374 : }
375 : }
376 :
377 : static void
378 62984 : set_bound(subgp_iter *T, GEN B)
379 : {
380 : GEN b;
381 62984 : T->bound = B;
382 62984 : if (!B) { T->boundtype = b_NONE; return; }
383 62844 : switch(typ(B))
384 : {
385 63 : case t_INT: /* upper bound */
386 63 : T->boundtype = b_MAX;
387 63 : break;
388 62774 : case t_VEC: /* exact value */
389 62774 : b = gel(B,1);
390 62774 : if (lg(B) != 2 || typ(b) != t_INT) pari_err_TYPE("subgroup", B);
391 62774 : T->boundtype = b_EXACT;
392 62774 : T->bound = b;
393 62774 : break;
394 7 : case t_COL: /* exact type */
395 7 : pari_err_IMPL("exact type in subgrouplist");
396 0 : if (lg(B) > len(T->L)+1) pari_err_TYPE("subgroup",B);
397 0 : T->boundtype = b_TYPE;
398 0 : break;
399 0 : default: pari_err_TYPE("subgroup",B);
400 : }
401 62837 : if (signe(T->bound) <= 0)
402 7 : pari_err_DOMAIN("subgroup", "index bound", "<=", gen_0, T->bound);
403 : }
404 :
405 : static GEN
406 329 : expand_sub(GEN x, long n)
407 : {
408 329 : long i,j, m = lg(x);
409 329 : GEN p = matid(n-1), q,c;
410 :
411 770 : for (i=1; i<m; i++)
412 : {
413 441 : q = gel(p,i); c = gel(x,i);
414 1218 : for (j=1; j<m; j++) gel(q,j) = gel(c,j);
415 693 : for ( ; j<n; j++) gel(q,j) = gen_0;
416 : }
417 329 : return p;
418 : }
419 :
420 : static GEN
421 47271 : init_powlist(long k, GEN p)
422 : {
423 47271 : GEN z = new_chunk(k+1);
424 : long i;
425 47271 : gel(z,0) = gen_1; gel(z,1) = p;
426 47901 : for (i=2; i<=k; i++) gel(z,i) = mulii(p, gel(z,i-1));
427 47271 : return z;
428 : }
429 :
430 : static GEN subgrouplist_i(GEN cyc, GEN bound, GEN expoI, GEN gen);
431 : static void
432 62972 : subgroup_engine(subgp_iter *T)
433 : {
434 62972 : pari_sp av = avma;
435 62972 : GEN B, L, p, listL, P = NULL, indexsubq = NULL, cyc = T->cyc;
436 62972 : long i, j, k, imax, lprim, n = lg(cyc);
437 :
438 62972 : if (n == 1) {
439 15701 : switch(T->boundtype)
440 : {
441 15687 : case b_EXACT: if (!is_pm1(T->bound)) break;
442 1435 : default: T->fun(T->fundata, cyc);
443 : }
444 15701 : set_avma(av); return;
445 : }
446 47271 : P = gel(Z_factor(cyc_get_expo(cyc)), 1);
447 47271 : listL = cgetg_copy(P, &lprim);
448 47271 : imax = k = 0;
449 94850 : for (i=1; i<lprim; i++)
450 : {
451 47579 : L = new_chunk(n); p = gel(P,i);
452 116410 : for (j=1; j<n; j++)
453 : {
454 68957 : L[j] = Z_pval(gel(cyc,j), p);
455 68957 : if (!L[j]) break;
456 : }
457 47579 : j--; setlen(L, j);
458 47579 : if (j > k) { k = j; imax = i; }
459 47579 : gel(listL,i) = L;
460 : }
461 47271 : L = gel(listL,imax); p = gel(P,imax);
462 47271 : k = L[1];
463 47271 : T->L = L;
464 47271 : T->powlist = (GEN*)init_powlist(k, p);
465 47271 : B = T->bound;
466 47271 : if (lprim == 2)
467 : {
468 47033 : T->subq = NULL;
469 47033 : if (T->boundtype == b_EXACT)
470 : {
471 46900 : (void)Z_pvalrem(T->bound,p,&B);
472 46899 : if (!is_pm1(B)) { set_avma(av); return; }
473 : }
474 : }
475 : else
476 : { /* not a p-group */
477 238 : GEN cycI = leafcopy(cyc);
478 : long lsubq;
479 490 : for (i=1; i<n; i++)
480 : {
481 308 : gel(cycI,i) = divii(gel(cycI,i), T->powlist[L[i]]);
482 308 : if (is_pm1(gel(cycI,i))) break;
483 : }
484 238 : setlg(cycI, i); /* cyclic factors of I */
485 238 : if (T->boundtype == b_EXACT)
486 : {
487 189 : (void)Z_pvalrem(T->bound,p,&B);
488 189 : B = mkvec(B);
489 : }
490 238 : T->expoI = gel(cycI,1);
491 238 : T->subq = subgrouplist_i(cycI, B, T->expoI, NULL);
492 :
493 238 : lsubq = lg(T->subq);
494 567 : for (i=1; i<lsubq; i++)
495 329 : gel(T->subq,i) = expand_sub(gel(T->subq,i), n);
496 238 : if (T->bound)
497 : {
498 203 : indexsubq = cgetg(lsubq,t_VEC);
499 462 : for (i=1; i<lsubq; i++)
500 259 : gel(indexsubq,i) = ZM_det_triangular(gel(T->subq,i));
501 : }
502 : /* lift subgroups of I to G */
503 567 : for (i=1; i<lsubq; i++)
504 329 : gel(T->subq,i) = gmul(T->powlist[k],gel(T->subq,i));
505 238 : if (DEBUGLEVEL>6)
506 0 : err_printf("(lifted) subgp of prime to %Ps part:\n%Ps\n",p, T->subq);
507 : }
508 47270 : dopsub(T, p,indexsubq);
509 47271 : if (DEBUGLEVEL>4) err_printf("nb subgroup = %ld\n",T->count);
510 47271 : set_avma(av);
511 : }
512 :
513 : static GEN
514 62992 : get_snf(GEN x, long *N)
515 : {
516 : GEN cyc;
517 : long n;
518 62992 : switch(typ(x))
519 : {
520 7 : case t_MAT:
521 7 : if (!RgM_isdiagonal(x)) return NULL;
522 7 : cyc = RgM_diagonal_shallow(x); break;
523 62985 : case t_VEC: if (lg(x) == 4 && typ(gel(x,2)) == t_VEC) x = gel(x,2);
524 62985 : case t_COL: cyc = leafcopy(x); break;
525 0 : default: return NULL;
526 : }
527 62992 : *N = lg(cyc)-1;
528 63167 : for (n = *N; n > 0; n--) /* take care of trailing 1s */
529 : {
530 49147 : GEN c = gel(cyc,n);
531 49147 : if (typ(c) != t_INT || signe(c) <= 0) return NULL;
532 49147 : if (!is_pm1(c)) break;
533 : }
534 62992 : setlg(cyc, n+1);
535 133229 : for ( ; n > 0; n--)
536 : {
537 70244 : GEN c = gel(cyc,n);
538 70244 : if (typ(c) != t_INT || signe(c) <= 0 ||
539 70244 : (n != *N && !dvdii(c, gel(cyc,n+1)))) return NULL;
540 : }
541 62985 : return cyc;
542 : }
543 :
544 : void
545 56 : forsubgroup(void *E, long call(void*, GEN), GEN cyc, GEN bound)
546 : {
547 : subgp_iter T;
548 : long N;
549 :
550 56 : T.fun = call;
551 56 : T.cyc = get_snf(cyc,&N);
552 56 : if (!T.cyc) pari_err_TYPE("forsubgroup [not a finite group]",cyc);
553 49 : set_bound(&T, bound);
554 49 : T.fundata = E;
555 49 : T.stop = 0;
556 49 : subgroup_engine(&T);
557 49 : }
558 :
559 : void
560 56 : forsubgroup0(GEN cyc, GEN bound, GEN code)
561 56 : { EXPRVOID_WRAP(code, forsubgroup(EXPR_ARGVOID, cyc, bound)) }
562 :
563 : static GEN
564 101339 : packtoi(long *pt, long L)
565 : {
566 : long i, l;
567 : GEN z;
568 119750 : for (i=0; i<L; i++, pt++)
569 101339 : if (*pt) break;
570 101339 : L -= i;
571 101339 : if (!L) return gen_0;
572 82928 : l = L+2; z = cgeti(l);
573 82928 : z[1] = evalsigne(1) | evallgefint(l);
574 165856 : for (i = 2; i<l; i++) z[i] = *pt++;
575 82928 : return z;
576 : }
577 :
578 : /* in place, remove trailing 1s */
579 : static void
580 62789 : snf_clean(GEN c)
581 : {
582 : long n;
583 64483 : for (n = lg(c)-1; n > 0; n--)
584 48811 : if (!is_pm1(gel(c,n))) break;
585 62789 : setlg(c, n+1);
586 62789 : }
587 : static GEN
588 62921 : update_cyc(subgp_iter *T, GEN cyc)
589 : {
590 : ulong k;
591 62921 : switch(T->boundtype)
592 : {
593 62746 : case b_EXACT:
594 62746 : cyc = ZV_snf_gcd(cyc, T->bound);
595 62747 : snf_clean(cyc); break;
596 42 : case b_MAX:
597 42 : if ((k = itos_or_0(T->bound)))
598 : {
599 42 : GEN fa = absZ_factor_limit_strict(cyc_get_expo(cyc), k + 1, NULL);
600 42 : cyc = T->cyc = ZV_snf_gcd(cyc, factorback(fa));
601 42 : snf_clean(cyc); break;
602 : }
603 : }
604 62922 : return cyc;
605 : }
606 :
607 : static GEN
608 62936 : subgrouplist_i(GEN CYC, GEN bound, GEN expoI, GEN gen)
609 : {
610 62936 : pari_sp av = avma;
611 : subgp_iter T;
612 : sublist_t S;
613 : slist *list, *sublist;
614 62936 : long ii,i,j,nbsub,n,N = 0; /* -Wall */
615 : GEN z, H, cyc;
616 :
617 62936 : cyc = get_snf(CYC, &N);
618 62935 : if (!cyc) pari_err_TYPE("subgrouplist [not a finite group]",CYC);
619 62935 : set_bound(&T, bound);
620 62921 : cyc = update_cyc(&T, cyc);
621 62922 : n = lg(cyc)-1; /* not necessarily = N */
622 62922 : S.list = sublist = (slist*) pari_malloc(sizeof(slist));
623 62923 : S.cyc = cyc;
624 62923 : S.gen = gen;
625 62923 : S.count = 0;
626 62923 : T.fun = &list_fun;
627 62923 : T.fundata = (void*)&S;
628 62923 : T.stop = 0;
629 62923 : T.cyc = cyc;
630 62923 : T.expoI = expoI;
631 62923 : subgroup_engine(&T);
632 :
633 62923 : nbsub = (long)S.count;
634 62923 : set_avma(av);
635 62923 : z = cgetg(nbsub+1,t_VEC);
636 111125 : for (ii=1; ii<=nbsub; ii++)
637 : {
638 : long *pt, L;
639 48202 : list = sublist; sublist = list->next; pari_free(list);
640 48202 : pt = sublist->data;
641 48202 : L = sublist->prec;
642 48202 : H = cgetg(N+1,t_MAT); gel(z,ii) = H;
643 118580 : for (j=1; j<=n; j++)
644 : {
645 70378 : gel(H,j) = cgetg(N+1, t_COL);
646 171717 : for (i=1; i<=j; i++) { gcoeff(H,i,j) = packtoi(pt, L); pt += L; }
647 102235 : for ( ; i<=N; i++) gcoeff(H,i,j) = gen_0;
648 : }
649 50029 : for ( ; j<=N; j++) gel(H,j) = col_ei(N, j);
650 : }
651 62923 : pari_free(sublist); return z;
652 : }
653 :
654 : GEN
655 6825 : subgrouplist(GEN cyc, GEN bound)
656 6825 : { return subgrouplist_i(cyc,bound,NULL,NULL); }
657 : GEN
658 55872 : subgroupcondlist(GEN cyc, GEN bound, GEN L)
659 55872 : { return subgrouplist_i(cyc,bound,NULL,L); }
|
