OpenXM_contrib/pari-2.2/src/basemath/subgroup.c - diff

Return to subgroup.c CVS log

Up to [local] / OpenXM_contrib / pari-2.2 / src / basemath

Diff for /OpenXM_contrib/pari-2.2/src/basemath/Attic/subgroup.c between version 1.1 and 1.2

-version 1.1, 2001/10/02 11:17:05
+version 1.2, 2002/09/11 07:26:52
 Line 18  extern GEN hnf0(GEN x, long r);
 Line 18  extern GEN hnf0(GEN x, long r);
 Line 18  extern GEN hnf0(GEN x, long r);
  void push_val(entree *ep, GEN a);
  void pop_val(entree *ep);
+ typedef struct {
+   entree *ep;
+   char *s;
+ } expr_t;
+ typedef struct slist {
+   struct slist *next;
+   long *data;
+ } slist;
+ typedef struct {
+   GEN hnfgroup;
+   GEN listKer;
+   ulong count;
+   slist *list;
+ } sublist_t;
+ /* MAX: [G:H] <= bound, EXACT: [G:H] = bound, TYPE: type(H) = bound */
+ enum { b_NONE, b_MAX, b_EXACT, b_TYPE };
  /* SUBGROUPS
-  * Assume: G = Gp x I, with Gp a p-group and (|I|,p)=1, and I small.
+  * G = Gp x I, with Gp a p-Sylow (I assumed small).
   * Compute subgroups of I by recursive calls
   * Loop through subgroups Hp of Gp using Birkhoff's algorithm.
   * If (I is non trivial)
   *   lift Hp to G (mul by exponent of I)
   *   for each subgp of I, lift it to G (mult by exponent of Gp)
   *   consider the group generated by the two subgroups (concat)
-  */
+  *
- static long *powerlist, *mmu, *lam, *c, *maxc, *a, *maxa, **g, **maxg;
+  * type(H) = mu --> H = Z/p^mu[1] x ... x Z/p^mu[len(mu)] */
- static GEN **H, subq, subqpart, hnfgroup;
+ typedef struct subgp_iter {
- static GEN BINMAT;
+   long *M, *L; /* mu = p-subgroup type, lambda = p-group type */
- static long countsub, expoI;
+   long *powlist; /* [i] = p^i, i = 0.. */
- static long *available, indexbound, lsubq, lsubqpart;
+   long *c, *maxc, *a, *maxa, **g, **maxg;
- static char *gpch;
+   long *available;
- static entree *ep;
+   GEN **H; /* p-subgroup of type mu, in matrix form */
- static void(*treatsub_fun)(GEN);
+   GEN cyc; /* cyclic factors of G */
- typedef struct slist {
+   GEN subq;/* subgrouplist(I) */
-   struct slist *next;
+   GEN subqpart; /* J in subq s.t [I:J][Gp:Hp] <= indexbound */
-   long *data;
+   GEN bound; /* if != NULL, impose a "bound" on [G:H] (see boundtype) */
- } slist;
+   int boundtype;
+   long countsub; /* number of subgroups of type M (so far) */
+   long count; /* number of p-subgroups so far [updated when M completed] */
+   GEN expoI; /* exponent of I */
+   void(*fun)(struct subgp_iter *, GEN); /* applied to each subgroup */
+   void *fundata; /* for fun */
+ } subgp_iter;
- static slist *sublist;
+ #define len(x)      (x)[0]
+ #define setlen(x,l) len(x)=(l)
  void
- printtyp(long *typ)
+ printtyp(const long *typ)
  {
-   long i;
+   long i, l = len(typ);
-   for (i=1; i<=typ[0]; i++) fprintferr(" %ld ",typ[i]);
+   for (i=1; i<=l; i++) fprintferr(" %ld ",typ[i]);
    fprintferr("\n");
  }
-Line 54  printtyp(long *typ)
+Line 81  printtyp(long *typ)
 Line 54  printtyp(long *typ)
 Line 81  printtyp(long *typ)
  static long*
  conjugate(long *typ)
  {
-   long *t, i, k = typ[0], l, last;
+   long *t, i, k = len(typ), l, last;
-   if (!k) { t = new_chunk(1); t[0]=0; return t; }
+   if (!k) { t = new_chunk(1); setlen(t,0); return t; }
    l = typ[1]; t = new_chunk(l+2);
    t[1] = k; last = k;
    for (i=2; i<=l; i++)
-Line 64  conjugate(long *typ)
+Line 91  conjugate(long *typ)
 Line 64  conjugate(long *typ)
 Line 91  conjugate(long *typ)
      while (typ[last] < i) last--;
      t[i] = last;
    }
-   t[i] = 0; t[0] = l;
+   t[i] = 0; setlen(t,l);
    return t;
  }
+ /* --- subgp_iter 'fun' associated to forsubgroup -------------- */
  static void
- std_fun(GEN x)
+ std_fun(subgp_iter *T, GEN x)
  {
-   ep->value = (void*)x;
+   expr_t *E = (expr_t *)T->fundata;
-   lisseq(gpch); countsub++;
+   E->ep->value = (void*)x;
+   (void)lisseq(E->s); T->countsub++;
  }
+ /* ----subgp_iter 'fun' associated to subgrouplist ------------- */
- void
+ static void
- addcell(GEN H)
+ addcell(sublist_t *S, GEN H)
  {
    long *pt,i,j, k = 0, n = lg(H)-1;
    slist *cell = (slist*) gpmalloc(sizeof(slist) + n*(n+1)/2 * sizeof(long));
-   sublist->next = cell; cell->data = pt = (long*) (cell + 1);
+   S->list->next = cell; cell->data = pt = (long*) (cell + 1);
    for (j=1; j<=n; j++)
      for(i=1; i<=j; i++) pt[k++] = itos(gcoeff(H,i,j));
-   sublist = cell;
+   S->list = cell;
+   S->count++;
  }
+ extern int hnfdivide(GEN A, GEN B);
+ /* 1 if h^(-1)*list[i] integral for some i, 0 otherwise */
+ static int
+ hnflistdivise(GEN list,GEN h)
+ {
+   long i, n = lg(list);
+   for (i=1; i<n; i++)
+     if (hnfdivide(h, (GEN)list[i])) break;
+   return i < n;
+ }
  static void
- list_fun(GEN x)
+ list_fun(subgp_iter *T, GEN x)
  {
-   addcell(hnf(concatsp(hnfgroup,x))); countsub++;
+   sublist_t *S = (sublist_t*)T->fundata;
+   GEN H = hnf(concatsp(S->hnfgroup,x));
+   if (S->listKer && hnflistdivise(S->listKer, H)) return; /* test conductor */
+   addcell(S, H); T->countsub++;
  }
+ /* -------------------------------------------------------------- */
- /* treat subgroup Hp (not in HNF, treatsub_fun should do it if desired) */
+ /* treat subgroup Hp (not in HNF, T->fun should do it if desired) */
  static void
- treatsub(GEN Hp)
+ treatsub(subgp_iter *T, GEN H)
  {
    long i;
-   if (!subq) treatsub_fun(Hp);
+   if (!T->subq) T->fun(T, H);
    else
    { /* not a p group, add the trivial part */
-     Hp = gmulsg(expoI,Hp); /* lift Hp to G */
+     GEN Hp = gmul(T->expoI, H); /* lift H to G */
-     for (i=1; i<lsubqpart; i++)
+     long l = lg(T->subqpart);
-       treatsub_fun(concatsp(Hp, (GEN)subqpart[i]));
+     for (i=1; i<l; i++)
+       T->fun(T, concatsp(Hp, (GEN)T->subqpart[i]));
    }
  }
  /* assume t>0 and l>1 */
  static void
- dogroup(void)
+ dogroup(subgp_iter *T)
  {
-   long av = avma,av1, e,i,j,k,r,n,t2,ind, t = mmu[0], l = lam[0];
+   const long *powlist = T->powlist;
+   long *M = T->M;
+   long *L = T->L;
+   long *c = T->c;
+   long  *a = T->a,  *maxa = T->maxa;
+   long **g = T->g, **maxg = T->maxg;
+   GEN **H = T->H;
+   gpmem_t av = avma;
+   long e,i,j,k,r,n,t2,ind, t = len(M), l = len(L);
    t2 = (l==t)? t-1: t;
    n = t2 * l - (t2*(t2+1))/2; /* number of gamma_ij */
    for (i=1, r=t+1; ; i++)
    {
-     if (available[i]) c[r++] = i;
+     if (T->available[i]) c[r++] = i;
      if (r > l) break;
    }
    if (DEBUGLEVEL>2) { fprintferr("    column selection:"); printtyp(c); }
-Line 127  dogroup(void)
+Line 181  dogroup(void)
 Line 127  dogroup(void)
 Line 181  dogroup(void)
      maxg[i] = maxa + (ind - (i+1)); /* only access maxg[i][i+1..l] */
      g[i] = a + (ind - (i+1));
      for (r=i+1; r<=l; r++)
-       if (c[r] < c[i])             maxg[i][r] = powerlist[mmu[i]-mmu[r]-1];
+       if (c[r] < c[i])         maxg[i][r] = powlist[M[i]-M[r]-1];
-       else if (lam[c[r]] < mmu[i]) maxg[i][r] = powerlist[lam[c[r]]-mmu[r]];
+       else if (L[c[r]] < M[i]) maxg[i][r] = powlist[L[c[r]]-M[r]];
-       else                         maxg[i][r] = powerlist[mmu[i]-mmu[r]];
+       else                     maxg[i][r] = powlist[M[i]-M[r]];
    }
-   av1=avma; a[n-1]=0; for (i=0; i<n-1; i++) a[i]=1;
+   av = avma; a[n-1]=0; for (i=0; i<n-1; i++) a[i]=1;
    for(;;)
    {
      a[n-1]++;
      if (a[n-1] > maxa[n-1])
      {
        j=n-2; while (j>=0 && a[j]==maxa[j]) j--;
-       if (j<0) { avma=av; return; }
+       if (j < 0) { avma = av; return; }
        a[j]++; for (k=j+1; k<n; k++) a[k]=1;
      }
      for (i=1; i<=t; i++)
      {
        for (r=1; r<i; r++) affsi(0, H[i][c[r]]);
-       affsi(powerlist[lam[c[r]]-mmu[r]], H[r][c[r]]);
+       affsi(powlist[L[c[r]] - M[r]], H[r][c[r]]);
        for (r=i+1; r<=l; r++)
        {
          if (c[r] < c[i])
-           e = g[i][r]*powerlist[lam[c[r]]-mmu[i]+1];
+           e = g[i][r] * powlist[L[c[r]] - M[i]+1];
          else
-           if (lam[c[r]] < mmu[i]) e = g[i][r];
+           if (L[c[r]] < M[i]) e = g[i][r];
-           else e = g[i][r]*powerlist[lam[c[r]]-mmu[i]];
+           else e = g[i][r] * powlist[L[c[r]] - M[i]];
          affsi(e, H[i][c[r]]);
        }
      }
-     treatsub((GEN)H); avma=av1;
+     treatsub(T, (GEN)H); avma = av;
    }
  }
- /* c[1],...,c[r-1] filled */
+ /* T->c[1],...,T->c[r-1] filled */
  static void
- loop(long r)
+ loop(subgp_iter *T, long r)
  {
    long j;
-   if (r > mmu[0]) { dogroup(); return; }
+   if (r > len(T->M)) { dogroup(T); return; }
-   if (r!=1 && (mmu[r-1] == mmu[r])) j = c[r-1]+1; else j = 1;
+   if (r!=1 && (T->M[r-1] == T->M[r])) j = T->c[r-1]+1; else j = 1;
-   for (  ; j<=maxc[r]; j++)
+   for (  ; j<=T->maxc[r]; j++)
-     if (available[j])
+     if (T->available[j])
      {
-       c[r] = j;  available[j] = 0;
+       T->c[r] = j;  T->available[j] = 0;
-       loop(r+1); available[j] = 1;
+       loop(T, r+1); T->available[j] = 1;
      }
  }
  static void
- dopsubtyp(void)
+ dopsubtyp(subgp_iter *T)
  {
-   long av = avma, i,r, l = lam[0], t = mmu[0];
+   gpmem_t av = avma;
+   long i,r, l = len(T->L), t = len(T->M);
    if (!t)
    {
      GEN p1 = cgetg(2,t_MAT);
      p1[1] = (long)zerocol(l);
-     treatsub(p1); avma=av; return;
+     treatsub(T, p1); avma = av; return;
    }
    if (l==1) /* imply t = 1 */
    {
-     GEN p1 = gtomat(stoi(powerlist[lam[1]-mmu[1]]));
+     GEN p1 = gtomat(stoi(T->powlist[T->L[1]-T->M[1]]));
-     treatsub(p1); avma=av; return;
+     treatsub(T, p1); avma = av; return;
    }
-   c = new_chunk(l+1); c[0] = l;
+   T->c         = new_chunk(l+1); setlen(T->c, l);
-   maxc = new_chunk(l+1);
+   T->maxc      = new_chunk(l+1);
-   available = new_chunk(l+1);
+   T->available = new_chunk(l+1);
-   a  = new_chunk(l*(t+1));
+   T->a   = new_chunk(l*(t+1));
-   maxa=new_chunk(l*(t+1));
+   T->maxa= new_chunk(l*(t+1));
-   g = (long**)new_chunk(t+1);
+   T->g    = (long**)new_chunk(t+1);
-   maxg = (long**)new_chunk(t+1);
+   T->maxg = (long**)new_chunk(t+1);
-   if (DEBUGLEVEL) { fprintferr("  subgroup:"); printtyp(mmu); }
+   if (DEBUGLEVEL) { fprintferr("  subgroup:"); printtyp(T->M); }
    for (i=1; i<=t; i++)
    {
      for (r=1; r<=l; r++)
-       if (mmu[i] > lam[r]) break;
+       if (T->M[i] > T->L[r]) break;
-     maxc[i] = r-1;
+     T->maxc[i] = r-1;
    }
-   H = (GEN**)cgetg(t+1, t_MAT);
+   T->H = (GEN**)cgetg(t+1, t_MAT);
    for (i=1; i<=t; i++)
    {
-     H[i] = (GEN*)cgetg(l+1, t_COL);
+     T->H[i] = (GEN*)cgetg(l+1, t_COL);
-     for (r=1; r<=l; r++) H[i][r] = cgeti(3);
+     for (r=1; r<=l; r++) T->H[i][r] = cgeti(3);
    }
-   for (i=1; i<=l; i++) available[i]=1;
+   for (i=1; i<=l; i++) T->available[i]=1;
-   for (i=1; i<=t; i++) c[i]=0;
+   for (i=1; i<=t; i++) T->c[i]=0;
    /* go through all column selections */
-   loop(1); avma=av; return;
+   loop(T, 1); avma = av; return;
  }
  static long
  weight(long *typ)
  {
-   long i,w = 0;
+   long i, l = len(typ), w = 0;
-   for (i=1; i<=typ[0]; i++) w += typ[i];
+   for (i=1; i<=l; i++) w += typ[i];
    return w;
  }
- static long
+ static void
- dopsub(long p, long *gtyp, long *indexsubq)
+ dopsub(subgp_iter *T, GEN p, GEN indexsubq)
  {
-   long w,i,j,k,n, wg = 0, wmin = 0, count = 0;
+   long *M, *L = T->L;
+   long w,i,j,k, wG = weight(L), wmin = 0, wmax = wG, n = len(L);
-   if (DEBUGLEVEL) { fprintferr("\ngroup:"); printtyp(gtyp); }
+   if (DEBUGLEVEL) { fprintferr("\ngroup:"); printtyp(L); }
-   if (indexbound)
+   T->count = 0;
+   switch(T->boundtype)
    {
-     wg = weight(gtyp);
+   case b_MAX: /* upper bound */
-     wmin = (long) (wg - (log((double)indexbound) / log((double)p)));
+     wmin = (long) (wG - (log(gtodouble(T->bound)) / log(gtodouble(p))));
-     if (cmpis(gpuigs(stoi(p), wg - wmin), indexbound) > 0) wmin++;
+     if (cmpii(gpowgs(p, wG - wmin), T->bound) > 0) wmin++;
+     break;
+   case b_EXACT: /* exact value */
+     wmin = wmax = wG - ggval(T->bound, p);
+     break;
    }
-   lam = gtyp; n = lam[0]; mmu = new_chunk(n+1);
+   T->M = M = new_chunk(n+1);
-   mmu[1] = -1; for (i=2; i<=n; i++) mmu[i]=0;
+   M[1] = -1; for (i=2; i<=n; i++) M[i]=0;
    for(;;) /* go through all vectors mu_{i+1} <= mu_i <= lam_i */
    {
-     mmu[1]++;
+     M[1]++;
-     if (mmu[1] > lam[1])
+     if (M[1] > L[1])
      {
        for (j=2; j<=n; j++)
-         if (mmu[j] < lam[j] && mmu[j] < mmu[j-1]) break;
+         if (M[j] < L[j] && M[j] < M[j-1]) break;
-       if (j>n) return count;
+       if (j > n) return;
-       mmu[j]++; for (k=1; k<j; k++) mmu[k]=mmu[j];
+       M[j]++; for (k=1; k<j; k++) M[k]=M[j];
      }
      for (j=1; j<=n; j++)
-       if (!mmu[j]) break;
+       if (!M[j]) break;
-     mmu[0] = j-1; w = weight(mmu);
+     setlen(M, j-1); w = weight(M);
-     if (w >= wmin)
+     if (w >= wmin && w <= wmax)
      {
        GEN p1 = gun;
-       if (subq) /* G not a p-group */
+       if (T->subq) /* G not a p-group */
        {
-         if (indexbound)
+         if (T->bound)
          {
-           long indexH = itos(gpuigs(stoi(p), wg - w));
+           GEN indexH = gpowgs(p, wG - w);
-           long bound = indexbound / indexH;
+           GEN B = divii(T->bound, indexH);
-           subqpart = cgetg(lsubq, t_VEC);
+           long k, l = lg(T->subq);
-           lsubqpart = 1;
+           T->subqpart = cgetg(l, t_VEC);
-           for (i=1; i<lsubq; i++)
+           k = 1;
-             if (indexsubq[i] <= bound) subqpart[lsubqpart++] = subq[i];
+           for (i=1; i<l; i++)
+             if (cmpii((GEN)indexsubq[i], B) <= 0)
+               T->subqpart[k++] = T->subq[i];
+           setlg(T->subqpart, k);
          }
-         else { subqpart = subq; lsubqpart = lsubq; }
+         else T->subqpart = T->subq;
        }
        if (DEBUGLEVEL)
        {
-         long *lp = conjugate(lam);
+         long *Lp = conjugate(L);
-         long *mp = conjugate(mmu);
+         long *Mp = conjugate(M);
+         GEN BINMAT = matqpascal(len(L)+1, p);
          if (DEBUGLEVEL > 3)
          {
-           fprintferr("    lambda = "); printtyp(lam);
+           fprintferr("    lambda = "); printtyp(L);
-           fprintferr("    lambda'= "); printtyp(lp);
+           fprintferr("    lambda'= "); printtyp(Lp);
-           fprintferr("    mu = "); printtyp(mmu);
+           fprintferr("    mu = "); printtyp(M);
-           fprintferr("    mu'= "); printtyp(mp);
+           fprintferr("    mu'= "); printtyp(Mp);
          }
-         for (j=1; j<=mp[0]; j++)
+         for (j=1; j<=len(Mp); j++)
          {
-           p1 = mulii(p1, gpuigs(stoi(p), mp[j+1]*(lp[j]-mp[j])));
+           p1 = mulii(p1, gpuigs(p, Mp[j+1]*(Lp[j]-Mp[j])));
-           p1 = mulii(p1, gcoeff(BINMAT, lp[j]-mp[j+1]+1, mp[j]-mp[j+1]+1));
+           p1 = mulii(p1, gcoeff(BINMAT, Lp[j]-Mp[j+1]+1, Mp[j]-Mp[j+1]+1));
          }
          fprintferr("  alpha_lambda(mu,p) = %Z\n",p1);
        }
-       countsub = 0;
+       T->countsub = 0;
-       dopsubtyp();
-       count += countsub;
+       dopsubtyp(T);
+       T->count += T->countsub;
        if (DEBUGLEVEL)
        {
-         fprintferr("  countsub = %ld\n", countsub);
+         fprintferr("  countsub = %ld\n", T->countsub);
          msgtimer("for this type");
-         if (subq) p1 = mulis(p1,lsubqpart-1);
+         if (T->fun != list_fun || !((sublist_t*)(T->fundata))->listKer)
-         if (cmpis(p1,countsub))
          {
-           fprintferr("  alpha = %Z\n",p1);
+           if (T->subq) p1 = mulis(p1,lg(T->subqpart)-1);
-           err(bugparier,"forsubgroup (alpha != countsub)");
+           if (cmpis(p1,T->countsub))
+           {
+             fprintferr("  alpha = %Z\n",p1);
+             err(bugparier,"forsubgroup (alpha != countsub)");
+           }
          }
        }
      }
    }
  }
+ void
+ parse_bound(subgp_iter *T)
+ {
+   GEN b, B = T->bound;
+   if (!B) { T->boundtype = b_NONE; return; }
+   switch(typ(B))
+   {
+   case t_INT: /* upper bound */
+     T->boundtype = b_MAX;
+     break;
+   case t_VEC: /* exact value */
+     b = (GEN)B[1];
+     if (lg(B) != 2 || typ(b) != t_INT) err(typeer,"subgroup");
+     T->boundtype = b_EXACT;
+     T->bound = b;
+     break;
+   case t_COL: /* exact type */
+     if (lg(B) > len(T->L)+1) err(typeer,"subgroup");
+     err(impl,"exact type in subgrouplist");
+     T->boundtype = b_TYPE;
+     break;
+   default: err(typeer,"subgroup");
+   }
+   if (signe(T->bound) <= 0)
+     err(talker,"subgroup: index bound must be positive");
+ }
  static GEN
  expand_sub(GEN x, long n)
  {
-Line 324  expand_sub(GEN x, long n)
+Line 423  expand_sub(GEN x, long n)
 Line 324  expand_sub(GEN x, long n)
 Line 423  expand_sub(GEN x, long n)
  }
  extern GEN matqpascal(long n, GEN q);
+ static GEN subgrouplist_i(GEN cyc, GEN bound, GEN expoI, GEN listKer);
- static long
+ static GEN
- subgroup_engine(GEN cyc, long bound)
+ init_powlist(long k, long p)
  {
-   long av=avma,i,j,k,imax,nbprim,count, n = lg(cyc);
+   GEN z = new_chunk(k+1);
-   GEN gtyp,fa,junk,primlist,p,listgtyp,indexsubq = NULL;
+   long i;
-   long oindexbound = indexbound;
+   z[0] = 1; z[1] = (long)p;
-   long oexpoI      = expoI;
+   for (i=1; i<=k; i++)
-   long *opowerlist = powerlist;
+   {
-   GEN osubq        = subq;
+     GEN t = muluu(p, z[i-1]);
-   GEN oBINMAT      = BINMAT;
+     z[i] = itos(t);
-   long olsubq      = lsubq;
+   }
+   return z;
+ }
-   long *ommu=mmu, *olam=lam, *oc=c, *omaxc=maxc, *oa=a, *omaxa=maxa, **og=g, **omaxg=maxg;
+ static void
-   GEN **oH=H, osubqpart=subqpart;
+ subgroup_engine(subgp_iter *T)
-   long ocountsub=countsub;
+ {
-   long *oavailable=available, olsubqpart=lsubqpart;
+   gpmem_t av = avma;
+   GEN B,L,fa,junk,primlist,p,listL,indexsubq = NULL;
-   if (typ(cyc) != t_VEC)
+   GEN cyc = T->cyc;
+   long i,j,k,imax,nbprim, n = lg(cyc);
+   if (typ(cyc) != t_VEC)
    {
      if (typ(cyc) != t_MAT) err(typeer,"forsubgroup");
      cyc = mattodiagonal(cyc);
-Line 350  subgroup_engine(GEN cyc, long bound)
+Line 455  subgroup_engine(GEN cyc, long bound)
 Line 350  subgroup_engine(GEN cyc, long bound)
 Line 455  subgroup_engine(GEN cyc, long bound)
    for (i=1; i<n-1; i++)
      if (!divise((GEN)cyc[i], (GEN)cyc[i+1]))
        err(talker,"not a group in forsubgroup");
-   if (n == 1 || gcmp1((GEN)cyc[1])) { treatsub(cyc); return 1; }
+   if (n == 1) { T->fun(T, cyc); return; }
    if (!signe(cyc[1]))
      err(talker,"infinite group in forsubgroup");
-   if (DEBUGLEVEL) timer2();
+   if (DEBUGLEVEL) (void)timer2();
-   indexbound = bound;
    fa = factor((GEN)cyc[1]); primlist = (GEN)fa[1];
    nbprim = lg(primlist);
-   listgtyp = new_chunk(n); imax = k = 0;
+   listL = new_chunk(n); imax = k = 0;
    for (i=1; i<nbprim; i++)
    {
-     gtyp = new_chunk(n); p = (GEN)primlist[i];
+     L = new_chunk(n); p = (GEN)primlist[i];
      for (j=1; j<n; j++)
      {
-       gtyp[j] = pvaluation((GEN)cyc[j], p, &junk);
+       L[j] = pvaluation((GEN)cyc[j], p, &junk);
-       if (!gtyp[j]) break;
+       if (!L[j]) break;
      }
-     j--; gtyp[0] = j;
+     j--; setlen(L, j);
      if (j > k) { k = j; imax = i; }
-     listgtyp[i] = (long)gtyp;
+     listL[i] = (long)L;
    }
-   gtyp = (GEN)listgtyp[imax]; p = (GEN)primlist[imax];
+   L = (GEN)listL[imax]; p = (GEN)primlist[imax];
-   k = gtyp[1];
+   k = L[1];
-   powerlist = new_chunk(k+1); powerlist[0] = 1;
+   T->L = L;
-   powerlist[1] = itos(p);
+   T->powlist = init_powlist(k, itos(p));
-   for (j=1; j<=k; j++) powerlist[j] = powerlist[1] * powerlist[j-1];
+   B = T->bound;
+   parse_bound(T);
-   if (DEBUGLEVEL) BINMAT = matqpascal(gtyp[0]+1, p);
+   if (nbprim == 2)
-   if (nbprim == 2) subq = NULL;
+   {
+     T->subq = NULL;
+     if (T->boundtype == b_EXACT)
+     {
+       (void)pvaluation(T->bound,p,&B);
+       if (!gcmp1(B)) { avma = av; return; }
+     }
+   }
    else
    { /* not a p-group */
-     GEN cyc2 = dummycopy(cyc);
+     GEN cycI = dummycopy(cyc);
-     GEN ohnfgroup = hnfgroup;
+     long lsubq;
      for (i=1; i<n; i++)
      {
-       cyc2[i] = ldivis((GEN)cyc2[i], powerlist[gtyp[i]]);
+       cycI[i] = ldivis((GEN)cycI[i], T->powlist[L[i]]);
-       if (gcmp1((GEN)cyc2[i])) break;
+       if (gcmp1((GEN)cycI[i])) break;
      }
-     setlg(cyc2, i);
+     setlg(cycI, i); /* cyclic factors of I */
-     if (is_bigint(cyc[1]))
+     if (T->boundtype == b_EXACT)
-       err(impl,"subgrouplist for large cyclic factors");
-     expoI = itos((GEN)cyc2[1]);
-     hnfgroup = diagonal(cyc2);
-     subq = subgrouplist(cyc2, bound);
-     hnfgroup = ohnfgroup;
-     lsubq = lg(subq);
-     for (i=1; i<lsubq; i++)
-       subq[i] = (long)expand_sub((GEN)subq[i], n);
-     if (indexbound)
      {
-       indexsubq = new_chunk(lsubq);
+       (void)pvaluation(T->bound,p,&B);
+       B = _vec(B);
+     }
+     T->expoI = (GEN)cycI[1];
+     T->subq = subgrouplist_i(cycI, B, T->expoI, NULL);
+     lsubq = lg(T->subq);
+     for (i=1; i<lsubq; i++)
+       T->subq[i] = (long)expand_sub((GEN)T->subq[i], n);
+     if (T->bound)
+     {
+       indexsubq = cgetg(lsubq,t_VEC);
        for (i=1; i<lsubq; i++)
-         indexsubq[i] = itos(dethnf((GEN)subq[i]));
+         indexsubq[i] = (long)dethnf_i((GEN)T->subq[i]);
      }
      /* lift subgroups of I to G */
      for (i=1; i<lsubq; i++)
-       subq[i] = lmulsg(powerlist[k],(GEN)subq[i]);
+       T->subq[i] = lmulsg(T->powlist[k],(GEN)T->subq[i]);
      if (DEBUGLEVEL>2)
      {
        fprintferr("(lifted) subgp of prime to %Z part:\n",p);
-       outbeaut(subq);
+       outbeaut(T->subq);
      }
    }
-   count = dopsub(powerlist[1],gtyp,indexsubq);
+   dopsub(T, p,indexsubq);
-   if (DEBUGLEVEL) fprintferr("nb subgroup = %ld\n",count);
+   if (DEBUGLEVEL) fprintferr("nb subgroup = %ld\n",T->count);
+   avma = av;
-   mmu=ommu; lam=olam; c=oc; maxc=omaxc; a=oa; maxa=omaxa; g=og; maxg=omaxg;
-   H=oH; subqpart=osubqpart,
-   countsub=ocountsub;
-   available=oavailable; lsubqpart=olsubqpart;
-   indexbound = oindexbound;
-   expoI    = oexpoI;
-   powerlist = opowerlist;
-   subq      = osubq;
-   BINMAT    = oBINMAT;
-   lsubq     = olsubq;
-   avma=av; return count;
  }
  static GEN
- get_snf(GEN x)
+ get_snf(GEN x, long *N)
  {
    GEN cyc;
    long n;
-Line 439  get_snf(GEN x)
+Line 541  get_snf(GEN x)
 Line 439  get_snf(GEN x)
 Line 541  get_snf(GEN x)
      case t_MAT:
        if (!isdiagonal(x)) return NULL;
        cyc = mattodiagonal_i(x); break;
-     case t_VEC:
+     case t_VEC: if (lg(x) == 4 && typ(x[2]) == t_VEC) x = (GEN)x[2];
      case t_COL: cyc = dummycopy(x); break;
      default: return NULL;
    }
-   for (n = lg(cyc)-1; n > 1; n--) /* take care of trailing 1s */
+   *N = lg(cyc)-1;
+   for (n = *N; n > 0; n--) /* take care of trailing 1s */
    {
      GEN c = (GEN)cyc[n];
      if (typ(c) != t_INT) return NULL;
-Line 459  get_snf(GEN x)
+Line 562  get_snf(GEN x)
 Line 459  get_snf(GEN x)
 Line 562  get_snf(GEN x)
  }
  void
- forsubgroup(entree *oep, GEN cyc, long bound, char *och)
+ forsubgroup(entree *ep, GEN cyc, GEN bound, char *ch)
  {
-   entree *saveep = ep;
+   subgp_iter T;
-   char *savech = gpch;
+   expr_t E;
-   void(*savefun)(GEN) = treatsub_fun;
+   long N;
-   treatsub_fun = std_fun;
+   T.fun = &std_fun;
-   cyc = get_snf(cyc);
+   cyc = get_snf(cyc,&N);
    if (!cyc) err(typeer,"forsubgroup");
-   gpch = och;
+   T.bound = bound;
-   ep = oep;
+   T.cyc = cyc;
+   E.s = ch;
+   E.ep= ep; T.fundata = (void*)&E;
    push_val(ep, gzero);
-   (void)subgroup_engine(cyc,bound);
+   subgroup_engine(&T);
    pop_val(ep);
-   treatsub_fun = savefun;
-   gpch = savech;
-   ep = saveep;
  }
- GEN
+ static GEN
- subgrouplist(GEN cyc, long bound)
+ subgrouplist_i(GEN cyc, GEN bound, GEN expoI, GEN listKer)
  {
-   void(*savefun)(GEN) = treatsub_fun;
+   gpmem_t av = avma;
-   slist *olist = sublist, *list;
+   subgp_iter T;
-   long ii,i,j,k,nbsub,n, N = lg(cyc)-1, av = avma;
+   sublist_t S;
+   slist *list, *sublist;
+   long ii,i,j,k,nbsub,n,N;
    GEN z,H;
-   GEN ohnfgroup = hnfgroup;
-   sublist = list = (slist*) gpmalloc(sizeof(slist));
+   cyc = get_snf(cyc, &N);
-   treatsub_fun = list_fun;
-   cyc = get_snf(cyc);
    if (!cyc) err(typeer,"subgrouplist");
-   n = lg(cyc)-1;
+   n = lg(cyc)-1; /* not necessarily = N */
-   hnfgroup = diagonal(cyc);
-   nbsub = subgroup_engine(cyc,bound);
+   S.list = sublist = (slist*) gpmalloc(sizeof(slist));
-   hnfgroup = ohnfgroup; avma = av;
+   S.hnfgroup = diagonal(cyc);
-   z = cgetg(nbsub+1,t_VEC); sublist = list;
+   S.listKer = listKer;
+   S.count = 0;
+   T.fun = &list_fun;
+   T.fundata = (void*)&S;
+   T.cyc = cyc;
+   T.bound = bound;
+   T.expoI = expoI;
+   subgroup_engine(&T);
+   nbsub = S.count;
+   avma = av;
+   z = cgetg(nbsub+1,t_VEC);
    for (ii=1; ii<=nbsub; ii++)
    {
      list = sublist; sublist = list->next; free(list);
-Line 512  subgrouplist(GEN cyc, long bound)
+Line 628  subgrouplist(GEN cyc, long bound)
 Line 512  subgrouplist(GEN cyc, long bound)
 Line 628  subgrouplist(GEN cyc, long bound)
        for (i=1; i<=N; i++) coeff(H,i,j) = (i==j)? un: zero;
      }
    }
-   free(sublist); sublist = olist;
+   free(sublist); return z;
-   treatsub_fun = savefun; return z;
+ }
+ GEN
+ subgrouplist(GEN cyc, GEN bound)
+ {
+   return subgrouplist_i(cyc,bound,NULL,NULL);
+ }
+ GEN
+ subgroupcondlist(GEN cyc, GEN bound, GEN listKer)
+ {
+   return subgrouplist_i(cyc,bound,NULL,listKer);
  }

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>