[BACK]Return to perm.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib / pari-2.2 / src / basemath

Annotation of OpenXM_contrib/pari-2.2/src/basemath/perm.c, Revision 1.1

1.1     ! noro        1: /* $Id: perm.c,v 1.9 2002/05/29 17:58:23 bill Exp $
        !             2:
        !             3: Copyright (C) 2000  The PARI group.
        !             4:
        !             5: This file is part of the PARI/GP package.
        !             6:
        !             7: PARI/GP is free software; you can redistribute it and/or modify it under the
        !             8: terms of the GNU General Public License as published by the Free Software
        !             9: Foundation. It is distributed in the hope that it will be useful, but WITHOUT
        !            10: ANY WARRANTY WHATSOEVER.
        !            11:
        !            12: Check the License for details. You should have received a copy of it, along
        !            13: with the package; see the file 'COPYING'. If not, write to the Free Software
        !            14: Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */
        !            15:
        !            16: #include "pari.h"
        !            17:
        !            18: /*************************************************************************/
        !            19: /**                                                                     **/
        !            20: /**                   Routine for handling VECSMALL                     **/
        !            21: /**                                                                     **/
        !            22: /*************************************************************************/
        !            23:
        !            24: GEN vecsmall_const(long n, long c)
        !            25: {
        !            26:   long i;
        !            27:   GEN V=cgetg(n+1,t_VECSMALL);
        !            28:   for(i=1;i<=n;i++) V[i]=c;
        !            29:   return V;
        !            30: }
        !            31:
        !            32: GEN vecsmall_shorten(GEN v, long n)
        !            33: {
        !            34:   long i;
        !            35:   GEN V=cgetg(n+1,t_VECSMALL);
        !            36:   for(i=1;i<=n;i++) V[i]=v[i];
        !            37:   return V;
        !            38:
        !            39: }
        !            40:
        !            41: /*in place sort.*/
        !            42: void vecsmall_sort(GEN V)
        !            43: {
        !            44:   long i,j,k,l,m;
        !            45:   for(l=1;l<lg(V);l<<=1)
        !            46:     for(j=1;(j>>1)<lg(V);j+=l<<1)
        !            47:       for(k=j+l,i=j; i<k && k<j+(l<<1) && k<lg(V);)
        !            48:        if (V[i]>V[k])
        !            49:        {
        !            50:          long z=V[k];
        !            51:          for (m=k;m>i;m--)
        !            52:            V[m]=V[m-1];
        !            53:          V[i]=z;
        !            54:          k++;
        !            55:        }
        !            56:        else
        !            57:          i++;
        !            58: }
        !            59:
        !            60: GEN vecsmall_uniq(GEN V)
        !            61: {
        !            62:   gpmem_t ltop=avma;
        !            63:   GEN W;
        !            64:   long i,j;
        !            65:   if ( lg(V) == 1 ) return gcopy(V);
        !            66:   W=cgetg(lg(V),t_VECSMALL);
        !            67:   W[1]=V[1];
        !            68:   for(i=2,j=1;i<lg(V);i++)
        !            69:     if (V[i]!=W[j])
        !            70:       W[++j]=V[i];
        !            71:   setlg(W,j+1);
        !            72:   return gerepileupto(ltop,W);
        !            73: }
        !            74:
        !            75: int
        !            76: vecsmall_lexcmp(GEN x, GEN y)
        !            77: {
        !            78:   long lx,ly,l,i;
        !            79:   lx = lg(x);
        !            80:   ly = lg(y); l = min(lx,ly);
        !            81:   for (i=1; i<l; i++)
        !            82:     if (x[i]!=y[i])
        !            83:       return x[i]<y[i]?-1:1;
        !            84:   if (lx == ly) return 0;
        !            85:   return (lx < ly)? -1 : 1;
        !            86: }
        !            87:
        !            88: int
        !            89: vecsmall_prefixcmp(GEN x, GEN y)
        !            90: {
        !            91:   long lx,ly,l,i;
        !            92:   lx = lg(x);
        !            93:   ly = lg(y); l = min(lx,ly);
        !            94:   for (i=1; i<l; i++)
        !            95:     if (x[i]!=y[i])
        !            96:       return x[i]<y[i]?-1:1;
        !            97:   return 0;
        !            98: }
        !            99:
        !           100: /*Can be used on vector but with no copy*/
        !           101: GEN vecsmall_prepend(GEN V, long s)
        !           102: {
        !           103:   GEN res;
        !           104:   long l2 = lg(V);
        !           105:   long i;
        !           106:   res = cgetg(l2+1, typ(V));
        !           107:   res[1]=s;
        !           108:   for (i = 2; i <= l2; ++i)
        !           109:     res[i] = V[i - 1];
        !           110:   return res;
        !           111: }
        !           112:
        !           113: /*Can be used on vector but with no copy*/
        !           114: GEN vecsmall_append(GEN V, long s)
        !           115: {
        !           116:   GEN res;
        !           117:   long l2 = lg(V);
        !           118:   long i;
        !           119:   res = cgetg(l2+1, typ(V));
        !           120:   for (i = 1; i < l2; ++i)
        !           121:     res[i] = V[i];
        !           122:   res[l2]=s;
        !           123:   return res;
        !           124: }
        !           125:
        !           126: GEN vecsmall_concat(GEN u, GEN v)
        !           127: {
        !           128:   GEN res;
        !           129:   long l1 = lg(u)-1;
        !           130:   long l2 = lg(v)-1;
        !           131:   long i;
        !           132:   res = cgetg(l1+l2+1, t_VECSMALL);
        !           133:   for (i = 1; i <= l1; ++i)
        !           134:     res[i] = u[i];
        !           135:   for (i = 1; i <= l2; ++i)
        !           136:     res[i+l1] = v[i];
        !           137:   return res;
        !           138: }
        !           139:
        !           140:
        !           141: /*************************************************************************/
        !           142: /**                                                                     **/
        !           143: /**               Routine for handling bit vector                       **/
        !           144: /**                                                                     **/
        !           145: /*************************************************************************/
        !           146:
        !           147: GEN
        !           148: bitvec_alloc(long n)
        !           149: {
        !           150:   long l=1+(n>>TWOPOTBITS_IN_LONG);
        !           151:   return vecsmall_const(l,0);
        !           152: }
        !           153:
        !           154:
        !           155: GEN
        !           156: bitvec_shorten(GEN bitvec, long n)
        !           157: {
        !           158:   long l=1+(n>>TWOPOTBITS_IN_LONG);
        !           159:   return vecsmall_shorten(bitvec,l);
        !           160: }
        !           161:
        !           162: long
        !           163: bitvec_test(GEN bitvec, long b)
        !           164: {
        !           165:   long q=b>>TWOPOTBITS_IN_LONG;
        !           166:   long r=b&(BITS_IN_LONG-1);
        !           167:   return (bitvec[1+q]>>r)&1L;
        !           168: }
        !           169:
        !           170: void
        !           171: bitvec_set(GEN bitvec, long b)
        !           172: {
        !           173:   long q=b>>TWOPOTBITS_IN_LONG;
        !           174:   long r=b&(BITS_IN_LONG-1);
        !           175:   bitvec[1+q]|=1L<<r;
        !           176: }
        !           177:
        !           178: void
        !           179: bitvec_clear(GEN bitvec, long b)
        !           180: {
        !           181:   long q=b>>TWOPOTBITS_IN_LONG;
        !           182:   long r=b&(BITS_IN_LONG-1);
        !           183:   bitvec[1+q]&=~(1L<<r);
        !           184: }
        !           185:
        !           186: /*************************************************************************/
        !           187: /**                                                                     **/
        !           188: /**               Routine for handling vector of VECSMALL               **/
        !           189: /**                                                                     **/
        !           190: /*************************************************************************/
        !           191:
        !           192: GEN
        !           193: vecvecsmall_sort(GEN x)
        !           194: {
        !           195:   return gen_sort(x, 0, vecsmall_lexcmp);
        !           196: }
        !           197:
        !           198: GEN
        !           199: vecvecsmall_indexsort(GEN x)
        !           200: {
        !           201:   GEN p=gen_sort(x, cmp_IND | cmp_C, vecsmall_lexcmp);
        !           202:   return p;
        !           203: }
        !           204:
        !           205: long
        !           206: vecvecsmall_search(GEN x, GEN y, long flag)
        !           207: {
        !           208:   return gen_search(x,y,flag,vecsmall_prefixcmp);
        !           209: }
        !           210:
        !           211: /*************************************************************************/
        !           212: /**                                                                     **/
        !           213: /**                   Routine for handling permutations                 **/
        !           214: /**                                                                     **/
        !           215: /*************************************************************************/
        !           216:
        !           217: /* Permutations may be given by
        !           218:  * perm (VECSMALL): a bijection from 1...n to 1...n i-->perm[i]
        !           219:  * cyc (VEC of VECSMALL): a product of disjoint cycles.
        !           220:  */
        !           221:
        !           222: /* indentity permutation */
        !           223: /* Not a good name since l is not a perm...*/
        !           224: GEN
        !           225: perm_identity(long l)
        !           226: {
        !           227:   GEN     perm;
        !           228:   int     i;
        !           229:   perm = cgetg(l + 1, t_VECSMALL);
        !           230:   for (i = 1; i <= l; i++)
        !           231:     perm[i] = i;
        !           232:   return perm;
        !           233: }
        !           234:
        !           235: GEN
        !           236: cyclicperm(long l, long d)
        !           237: {
        !           238:   GEN     perm;
        !           239:   int     i;
        !           240:   perm = cgetg(l + 1, t_VECSMALL);
        !           241:   for (i = 1; i <= l-d; i++)
        !           242:     perm[i] = i+d;
        !           243:   for (i = l-d+1; i <= l; i++)
        !           244:     perm[i] = i-l+d;
        !           245:   return perm;
        !           246: }
        !           247:
        !           248: /*
        !           249:  * Multiply (compose) two permutations.
        !           250:  * Can be used if s is a vector but with no copy
        !           251:  */
        !           252: GEN
        !           253: perm_mul(GEN s, GEN t)
        !           254: {
        !           255:   GEN     u;
        !           256:   int     i;
        !           257:   if (lg(s) < lg(t))
        !           258:     err(talker, "First permutation shorter than second in perm_mul");
        !           259:   u = cgetg(lg(s), typ(s));
        !           260:   for (i = 1; i < lg(s); i++)
        !           261:     u[i] = s[t[i]];
        !           262:   return u;
        !           263: }
        !           264: /* Compute the inverse (reciprocal) of a permutation.
        !           265:  */
        !           266: GEN
        !           267: perm_inv(GEN x)
        !           268: {
        !           269:       long i,lx = lg(x);
        !           270:       GEN y = cgetg(lx,t_VECSMALL);
        !           271:       for (i=1; i<lx; i++) y[x[i]] = i;
        !           272:       return y;
        !           273: }
        !           274:
        !           275: /* Orbits of the subgroup generated by v on {1,..,n}
        !           276:  */
        !           277: GEN
        !           278: vecperm_orbits(GEN v, long n)
        !           279: {
        !           280:   gpmem_t ltop = avma;
        !           281:   int     j, k, l, m, o, p, flag;
        !           282:   GEN     bit, cycle, cy;
        !           283:   long    mj=1;
        !           284:   cycle = cgetg(n+1, t_VEC);
        !           285:   bit = bitvec_alloc(n);
        !           286:   for (k = 1, l = 1; k <= n;)
        !           287:   {
        !           288:     for (  ; bitvec_test(bit,mj); mj++);
        !           289:     cy = cgetg(n+1, t_VECSMALL);
        !           290:     m = 1;
        !           291:     cy[m++] = mj;
        !           292:     bitvec_set(bit,mj++);
        !           293:     k++;
        !           294:     do
        !           295:     {
        !           296:       flag = 0;
        !           297:       for (o = 1; o < lg(v); o++)
        !           298:       {
        !           299:        for (p = 1; p < m; p++) /* m varie! */
        !           300:        {
        !           301:          j = mael(v,o,cy[p]);
        !           302:          if (!bitvec_test(bit,j))
        !           303:          {
        !           304:            flag = 1;
        !           305:            bitvec_set(bit,j);
        !           306:            k++;
        !           307:            cy[m++] = j;
        !           308:          }
        !           309:        }
        !           310:       }
        !           311:     }
        !           312:     while (flag);
        !           313:     setlg(cy, m);
        !           314:     cycle[l++] = (long) cy;
        !           315:   }
        !           316:   setlg(cycle, l);
        !           317:   return gerepilecopy(ltop, cycle);
        !           318: }
        !           319:
        !           320: /* Compute the cyclic decomposition of a permutation
        !           321:  */
        !           322:
        !           323: GEN
        !           324: perm_cycles(GEN v)
        !           325: {
        !           326:   gpmem_t ltop = avma;
        !           327:   GEN  u = cgetg(2, t_VEC);
        !           328:   u[1] = (long) v;
        !           329:   return gerepileupto(ltop, vecperm_orbits(u, lg(v)-1));
        !           330: }
        !           331: /* Compute the power of a permutation given by product of cycles
        !           332:  * Ouput a perm, not a cyc.
        !           333:  * */
        !           334: GEN
        !           335: cyc_powtoperm(GEN cyc, long exp)
        !           336: {
        !           337:   int     j, k, n;
        !           338:   GEN     p;
        !           339:   for (n = 0, j = 1; j < lg(cyc); j++)
        !           340:     n += lg(cyc[j]) - 1;
        !           341:   p = cgetg(n + 1, t_VECSMALL);
        !           342:   for (j = 1; j < lg(cyc); j++)
        !           343:   {
        !           344:     n = lg(cyc[j]) - 1;
        !           345:     for (k = 1; k <= n; k++)
        !           346:       p[mael(cyc,j,k)] = mael(cyc,j,1 + (k + exp - 1) % n);
        !           347:   }
        !           348:   return p;
        !           349: }
        !           350:
        !           351: /*
        !           352:  * Compute the power of a permutation.
        !           353:  * TODO: make it more clever with small exp.
        !           354:  */
        !           355: GEN
        !           356: perm_pow(GEN perm, long exp)
        !           357: {
        !           358:   return cyc_powtoperm(perm_cycles(perm),exp);
        !           359: }
        !           360:
        !           361: /*************************************************************************/
        !           362: /**                                                                     **/
        !           363: /**                   Routine for handling groups                       **/
        !           364: /**                                                                     **/
        !           365: /*************************************************************************/
        !           366:
        !           367: /* Groups are vector [gen,orders]
        !           368:  * gen (vecvecsmall): list of generators given by permutations
        !           369:  * orders (vecsmall): relatives orders of generators.
        !           370:  */
        !           371:
        !           372: GEN trivialsubgroups(void)
        !           373: {
        !           374:   GEN p2,p3;     /* vec */
        !           375:   p2 = cgetg(2, t_VEC);
        !           376:   p3 = cgetg(3, t_VEC);
        !           377:   p3[1] = (long) cgetg(1,t_VEC);
        !           378:   p3[2] = (long) cgetg(1,t_VECSMALL);
        !           379:   p2[1] = (long) p3;
        !           380:   return p2;
        !           381: }
        !           382:
        !           383:
        !           384:
        !           385: /* Compute the orders of p modulo the group S given by a set.
        !           386:  */
        !           387: long
        !           388: perm_relorder(GEN p, GEN S)
        !           389: {
        !           390:   gpmem_t ltop=avma;
        !           391:   long n = 1;
        !           392:   GEN  q = p;
        !           393:   while (!vecvecsmall_search(S, q, 0))
        !           394:   {
        !           395:     q = perm_mul(q, p);
        !           396:     ++n;
        !           397:   }
        !           398:   ltop=avma;
        !           399:   return n;
        !           400: }
        !           401:
        !           402: GEN perm_generate(GEN S, GEN H, long o)
        !           403: {
        !           404:   long i,k;
        !           405:   long n = lg(H)-1;
        !           406:   GEN L = cgetg(1+n*o, t_VEC);
        !           407:   for(i=1; i<=n; i++)
        !           408:     L[i]=lcopy((GEN)H[i]);
        !           409:   for(k=n+1; k <= n*o; ++k)
        !           410:     L[k] = (long) perm_mul((GEN) L[k-n], S);
        !           411:   return L;
        !           412: }
        !           413:
        !           414: /*Return the order (cardinal) of a group */
        !           415:
        !           416: long group_order(GEN G)
        !           417: {
        !           418:   GEN ord=(GEN) G[2];
        !           419:   long i;
        !           420:   long card=1;
        !           421:   for (i = 1; i < lg(ord); i++)
        !           422:     card *= ord[i];
        !           423:   return card;
        !           424: }
        !           425:
        !           426: /*Compute the left coset of g mod G: gG*/
        !           427:
        !           428: GEN group_leftcoset(GEN G, GEN g)
        !           429: {
        !           430:   GEN res;
        !           431:   long i,j,k;
        !           432:   GEN gen=(GEN) G[1];
        !           433:   GEN ord=(GEN) G[2];
        !           434:   long card=group_order(G);
        !           435:   res = cgetg(card + 1, t_VEC);
        !           436:   res[1] = lcopy(g);
        !           437:   k = 1;
        !           438:   for (i = 1; i < lg(gen); i++)
        !           439:   {
        !           440:     int     c = k * (ord[i] - 1);
        !           441:     for (j = 1; j <= c; j++)   /* I like it */
        !           442:       res[++k] = (long) perm_mul((GEN) res[j], (GEN) gen[i]);
        !           443:   }
        !           444:   return res;
        !           445: }
        !           446:
        !           447: /*Compute the right coset of g mod G: Gg*/
        !           448:
        !           449: GEN group_rightcoset(GEN G, GEN g)
        !           450: {
        !           451:   GEN res;
        !           452:   long i,j,k;
        !           453:   GEN gen=(GEN) G[1];
        !           454:   GEN ord=(GEN) G[2];
        !           455:   long card=group_order(G);
        !           456:   res = cgetg(card + 1, t_VEC);
        !           457:   res[1] = lcopy(g);
        !           458:   k = 1;
        !           459:   for (i = 1; i < lg(gen); i++)
        !           460:   {
        !           461:     int     c = k * (ord[i] - 1);
        !           462:     for (j = 1; j <= c; j++)   /* I like it */
        !           463:       res[++k] = (long) perm_mul((GEN) gen[i], (GEN) res[j]);
        !           464:   }
        !           465:   return res;
        !           466: }
        !           467:
        !           468: /*Compute the elements of a group from the generators*/
        !           469: /*Not stack clean!*/
        !           470:
        !           471: GEN group_elts(GEN G, long n)
        !           472: {
        !           473:   return group_leftcoset(G,perm_identity(n));
        !           474: }
        !           475:
        !           476: /*Return the cyclic group generated by g of order s*/
        !           477:
        !           478: GEN cyclicgroup(GEN g, long s)
        !           479: {
        !           480:   GEN p2,p3,p4;
        !           481:   p2 = cgetg(3, t_VEC);
        !           482:   p3 = cgetg(2, t_VEC);
        !           483:   p3[1] = lcopy(g);
        !           484:   p4 = cgetg(2,t_VECSMALL);
        !           485:   p4[1] = s;
        !           486:   p2[1] = (long) p3;
        !           487:   p2[2] = (long) p4;
        !           488:   return p2;
        !           489: }
        !           490:
        !           491: /*Return the group generated by g1,g2 of rel orders s1,s2*/
        !           492:
        !           493: GEN dicyclicgroup(GEN g1, GEN g2, long s1, long s2)
        !           494: {
        !           495:   GEN H = cgetg(3, t_VEC);
        !           496:   GEN p3,p4;
        !           497:   p3 = cgetg(3, t_VEC);
        !           498:   p3[1] = lcopy((GEN)g1);
        !           499:   p3[2] = lcopy((GEN)g2);
        !           500:   p4 = cgetg(3,t_VECSMALL);
        !           501:   p4[1] = s1;
        !           502:   p4[2] = s2;
        !           503:   H[1] = (long) p3;
        !           504:   H[2] = (long) p4;
        !           505:   return H;
        !           506: }
        !           507:
        !           508: /* return the quotient map G --> G/H */
        !           509: /*The ouput is [gen,hash]*/
        !           510: /* gen (vecvecsmall): coset generators
        !           511:  * hash (vecvecsmall): sorted vecsmall of concat(perm,coset number)
        !           512:  */
        !           513: GEN group_quotient(GEN G, GEN H)
        !           514: {
        !           515:   gpmem_t ltop=avma;
        !           516:   GEN p1,p2,p3;
        !           517:   long i,j,k;
        !           518:   long a=1;
        !           519:   long n=lg(mael(G,1,1))-1;
        !           520:   long o=group_order(H);
        !           521:   GEN elt = vecvecsmall_sort(group_elts(G,n));
        !           522:   GEN used = bitvec_alloc(lg(elt));
        !           523:   long l = (lg(elt)-1)/o;
        !           524:   p2 = cgetg(l+1, t_VEC);
        !           525:   p3 = cgetg(lg(elt), t_VEC);
        !           526:   for (i = 1, k = 1; i <= l; ++i)
        !           527:   {
        !           528:     GEN V;
        !           529:     while(bitvec_test(used,a)) a++;
        !           530:     V = group_leftcoset(H,(GEN)elt[a]);
        !           531:     p2[i] = V[1];
        !           532:     for(j=1;j<lg(V);j++)
        !           533:     {
        !           534:       long b=vecvecsmall_search(elt,(GEN)V[j],0);
        !           535:       bitvec_set(used,b);
        !           536:     }
        !           537:     for (j = 1; j <= o; j++)
        !           538:       p3[k++] = (long) vecsmall_append((GEN) V[j],i);
        !           539:   }
        !           540:   setlg(p3,k);
        !           541:   p1 = cgetg(3,t_VEC);
        !           542:   p1[1] = lcopy(p2);
        !           543:   p1[2]= (long) vecvecsmall_sort(p3);
        !           544:   return gerepileupto(ltop,p1);
        !           545: }
        !           546:
        !           547: /*Find in which coset a perm lie.*/
        !           548:
        !           549: long
        !           550: cosets_perm_search(GEN C, GEN p)
        !           551: {
        !           552:   long n=gen_search((GEN) C[2],p,0,vecsmall_prefixcmp);
        !           553:   if (!n)
        !           554:     err(talker, "coset not found in cosets_perm_search");
        !           555:   return mael3(C,2,n,lg(p));
        !           556: }
        !           557:
        !           558: /*Compute the image of a permutation by a quotient map.*/
        !           559:
        !           560: GEN quotient_perm(GEN C, GEN p)
        !           561: {
        !           562:   GEN p3;
        !           563:   long j;
        !           564:   long l2 = lg(C[1]);
        !           565:   p3 = cgetg(l2, t_VECSMALL);
        !           566:   for (j = 1; j < l2; ++j)
        !           567:     p3[j] = cosets_perm_search(C, perm_mul(p, gmael(C,1,j)));
        !           568:   return p3;
        !           569: }
        !           570:
        !           571: /* H is a subgroup of G, C is the quotient map G --> G/H
        !           572:  *
        !           573:  * Lift a subgroup S of G/H to a subgroup of G containing H
        !           574:  */
        !           575:
        !           576: GEN quotient_subgroup_lift(GEN C, GEN H, GEN S)
        !           577: {
        !           578:   GEN p1,L;
        !           579:   long l1=lg(H[1])-1;
        !           580:   long l2=lg(S[1])-1;
        !           581:   long j;
        !           582:   p1 = cgetg(3, t_VEC);
        !           583:   L = cgetg(l1+l2+1, t_VEC);
        !           584:   for (j = 1; j <= l1; ++j)
        !           585:     L[j] = mael(H,1,j);
        !           586:   for (j = 1; j <= l2; ++j)
        !           587:     L[l1+j] = (long) gmael(C,1,mael3(S,1,j,1));
        !           588:   p1[1] = (long) L;
        !           589:   p1[2] = (long) vecsmall_concat((GEN)H[2],(GEN)S[2]);
        !           590:   return p1;
        !           591: }
        !           592:
        !           593: /* Let G a group and H a quotient map G --> G/H
        !           594:  * Assume H is normal, return the group G/H.
        !           595:  */
        !           596:
        !           597: GEN quotient_group(GEN C, GEN G)
        !           598: {
        !           599:   gpmem_t ltop=avma;
        !           600:   GEN Qgen,Qord,Qelt;
        !           601:   GEN Q;
        !           602:   long i,j;
        !           603:   long n=lg(C[1])-1;
        !           604:   long l=lg(G[1]);
        !           605:   Qord = cgetg(l, t_VECSMALL);
        !           606:   Qgen = cgetg(l, t_VEC);
        !           607:   Qelt = cgetg(2, t_VEC);
        !           608:   Qelt[1] = (long) perm_identity(n);
        !           609:   for (i = 1, j = 1; i < l; ++i)
        !           610:   {
        !           611:     Qgen[j] = (long) quotient_perm(C, gmael(G,1,i));
        !           612:     Qord[j] = (long) perm_relorder((GEN) Qgen[j], vecvecsmall_sort(Qelt));
        !           613:     if (Qord[j]!=1)
        !           614:     {
        !           615:       Qelt=perm_generate((GEN) Qgen[j], Qelt, Qord[j]);
        !           616:       j++;
        !           617:     }
        !           618:   }
        !           619:   setlg(Qgen,j); setlg(Qord,j);
        !           620:   Q=cgetg(3,t_VEC);
        !           621:   Q[1]=(long)Qgen;
        !           622:   Q[2]=(long)Qord;
        !           623:   return gerepilecopy(ltop,Q);
        !           624: }
        !           625:
        !           626: /* Test if g normalize N*/
        !           627: long group_perm_normalize(GEN N, GEN g)
        !           628: {
        !           629:   gpmem_t ltop=avma;
        !           630:   long l1 = gegal(vecvecsmall_sort(group_leftcoset(N, g)),
        !           631:                   vecvecsmall_sort(group_rightcoset(N, g)));
        !           632:   avma=ltop;
        !           633:   return l1;
        !           634: }
        !           635:
        !           636: /* L is a list of subgroups, C is a coset and r a rel. order.*/
        !           637: static
        !           638: GEN liftlistsubgroups(GEN L, GEN C, long r)
        !           639: {
        !           640:   gpmem_t ltop=avma;
        !           641:   GEN p4;
        !           642:   long i, k;
        !           643:   long c=lg(C)-1;
        !           644:   long l=lg(L)-1;
        !           645:   long n=lg(C[1])-1;
        !           646:   if ( !l )
        !           647:     return cgetg(1,t_VEC);
        !           648:   p4 = cgetg(l*c+1, t_VEC);
        !           649:   for (i = 1, k = 1; i <= l; ++i)
        !           650:   {
        !           651:     long j;
        !           652:     GEN S = (GEN) L[i];
        !           653:     GEN Selt = vecvecsmall_sort(group_elts(S,n));
        !           654:     for (j = 1; j <= c; ++j)
        !           655:       if ((perm_relorder((GEN) C[j], Selt) == r) && group_perm_normalize(S, (GEN) C[j]))
        !           656:       {
        !           657:         GEN p7 = cgetg(3, t_VEC);
        !           658:         p7[1] = (long) vecsmall_append((GEN) S[1], C[j]);
        !           659:         p7[2] = (long) vecsmall_append((GEN) S[2], r);
        !           660:         p4[k++] = (long) p7;
        !           661:       }
        !           662:   }
        !           663:   setlg(p4,k);
        !           664:   return gerepilecopy(ltop,p4);
        !           665: }
        !           666:
        !           667: /* H is a normal subgroup, C is the quotient map G -->G/H,
        !           668:  * S is a subgroup of G/H, and G is embedded in Sym(l)
        !           669:  * Return all the subgroups K of G such that
        !           670:  * S= K mod H and K inter H={1}.
        !           671:  */
        !           672: static GEN liftsubgroup(GEN C, GEN H, GEN S)
        !           673: {
        !           674:   gpmem_t ltop=avma;
        !           675:   GEN V = trivialsubgroups();
        !           676:   long n = lg(S[1]);
        !           677:   long i;
        !           678:   /*Loop over generators of S*/
        !           679:   for (i = 1; i < n; ++i)
        !           680:   {
        !           681:     GEN W = group_leftcoset(H, gmael(C, 1, mael3(S, 1, i, 1)));
        !           682:     V = liftlistsubgroups(V, W, mael(S, 2, i));
        !           683:   }
        !           684:   return gerepilecopy(ltop,V);
        !           685: }
        !           686: /* compute all the subgroups of a group G
        !           687:  */
        !           688: GEN group_subgroups(GEN G)
        !           689: {
        !           690:   gpmem_t ltop=avma;
        !           691:   GEN p1;
        !           692:   GEN C,Q,M;
        !           693:   long lM;
        !           694:   GEN sg1,sg2,sg3;
        !           695:   long i, j;
        !           696:   GEN gen=(GEN)G[1], ord=(GEN)G[2];
        !           697:   GEN H;
        !           698:   long l, n = lg(gen);
        !           699:   if (n == 1)
        !           700:     return trivialsubgroups();
        !           701:   l = lg(gen[1]);/*now lg(gen)>1*/
        !           702:   if ( ( n == 4 || n == 5) && ord[1]==2 && ord[2]==2 && ord[3]==3
        !           703:       && (n == 4 || ord[4]==2) )
        !           704:   {
        !           705:     GEN u=perm_mul((GEN) gen[1],(GEN) gen[2]);
        !           706:     H = dicyclicgroup((GEN) gen[1],(GEN) gen[2],2,2);
        !           707:   /* sg3 is the list of subgroups intersecting only partially with H*/
        !           708:     sg3 = cgetg((n==4)?4:10, t_VEC);
        !           709:     sg3[1] = (long) cyclicgroup((GEN) gen[1], 2);
        !           710:     sg3[2] = (long) cyclicgroup((GEN) gen[2], 2);
        !           711:     sg3[3] = (long) cyclicgroup(u, 2);
        !           712:     if (n==5)
        !           713:     {
        !           714:       GEN s=(GEN) gen[1];
        !           715:       GEN t=(GEN) gen[2];
        !           716:       GEN u=(GEN) gen[3],u2=perm_mul(u,u);
        !           717:       GEN v=(GEN) gen[4];
        !           718:       GEN st=perm_mul(s,t);
        !           719:       GEN w=perm_mul(perm_mul(u2,perm_mul(s,v)),u2);
        !           720:       sg3[4] = (long) dicyclicgroup(s,v,2,2);
        !           721:       sg3[5] = (long) dicyclicgroup(t,perm_mul(u,perm_mul(v,u2)),2,2);
        !           722:       sg3[6] = (long) dicyclicgroup(st,perm_mul(u2,perm_mul(v,u)),2,2);
        !           723:       sg3[7] = (long) dicyclicgroup(s,w,2,2);
        !           724:       sg3[8] = (long) dicyclicgroup(t,perm_mul(u,perm_mul(w,u2)),2,2);
        !           725:       sg3[9] = (long) dicyclicgroup(st,perm_mul(u2,perm_mul(w,u)),2,2);
        !           726:     }
        !           727:   }
        !           728:   else
        !           729:   {
        !           730:     long osig = itos((GEN) coeff(decomp(stoi(ord[1])), 1, 1));
        !           731:     GEN sig = perm_pow((GEN) gen[1], ord[1]/osig);
        !           732:     H = cyclicgroup(sig,osig);
        !           733:     sg3=NULL;
        !           734:   }
        !           735:   C = group_quotient(G,H);
        !           736:   Q = quotient_group(C,G);
        !           737:   M = group_subgroups(Q);
        !           738:   lM = lg(M);
        !           739:   /* sg1 is the list of subgroups containing H*/
        !           740:   sg1 = cgetg(lM, t_VEC);
        !           741:   for (i = 1; i < lM; ++i)
        !           742:     sg1[i] = (long) quotient_subgroup_lift(C,H,(GEN)M[i]);
        !           743:   /*sg2 is a list of lists of subgroups not intersecting with H*/
        !           744:   sg2 = cgetg(lM, t_VEC);
        !           745:   /* Loop over all subgroups of G/H */
        !           746:   for (j = 1; j < lM; ++j)
        !           747:     sg2[j] = (long) liftsubgroup(C, H, (GEN) M[j]);
        !           748:   p1 = concat(sg1, concat(sg2, NULL));
        !           749:   if (sg3)
        !           750:     p1 = concat(p1, sg3);
        !           751:   return gerepileupto(ltop,p1);
        !           752: }
        !           753:
        !           754: /*return 1 if G is abelian, else 0*/
        !           755: long
        !           756: group_isabelian(GEN G)
        !           757: {
        !           758:   gpmem_t ltop=avma;
        !           759:   long i,j;
        !           760:   for(i=2;i<lg(G[1]);i++)
        !           761:     for(j=1;j<i;j++)
        !           762:     {
        !           763:       long test=gegal(perm_mul(gmael(G,1,i),gmael(G,1,j)),
        !           764:          perm_mul(gmael(G,1,j),gmael(G,1,i)));
        !           765:       avma=ltop;
        !           766:       if (!test) return 0;
        !           767:     }
        !           768:   return 1;
        !           769: }
        !           770:
        !           771: /*If G is abelian, return its HNF matrix*/
        !           772:
        !           773: GEN group_abelianHNF(GEN G)
        !           774: {
        !           775:   long i, j;
        !           776:   long n=lg(G[1]);
        !           777:   GEN M,S;
        !           778:   if (!group_isabelian(G)) return NULL;
        !           779:   S=group_elts(G,lg(mael(G,1,1)));
        !           780:   M=cgetg(n,t_MAT);
        !           781:   for(i=1;i<n;i++)
        !           782:   {
        !           783:     gpmem_t btop;
        !           784:     GEN P;
        !           785:     long k;
        !           786:     M[i]=lgetg(n,t_COL);
        !           787:     btop=avma;
        !           788:     P=perm_pow(gmael(G,1,i),mael(G,2,i));
        !           789:     for(j=1;j<lg(S);j++)
        !           790:       if (gegal(P,(GEN) S[j]))
        !           791:          break;
        !           792:     avma=btop;
        !           793:     if (j==lg(S)) err(talker,"wrong argument in galoisisabelian");
        !           794:     j--;
        !           795:     for(k=1;k<i;k++)
        !           796:     {
        !           797:       mael(M,i,k)=lstoi(j%mael(G,2,k));
        !           798:       j/=mael(G,2,k);
        !           799:     }
        !           800:     mael(M,i,k++)=lstoi(mael(G,2,i));
        !           801:     for(  ;k<n;k++)
        !           802:       mael(M,i,k)=zero;
        !           803:   }
        !           804:   return M;
        !           805: }
        !           806:
        !           807: #if 0
        !           808: /* Compute generators for the subgroup of (Z/nZ)* given in HNF.
        !           809:  * I apologize for the following spec:
        !           810:  * If zns=znstar(n) then
        !           811:  * zn2=gtovecsmall((GEN)zns[2]);
        !           812:  * zn3=lift((GEN)zns[3]);
        !           813:  * gen and ord : VECSMALL of length lg(zn3).
        !           814:  * the result is in gen.
        !           815:  * ord contains the relatives orders of the generators.
        !           816:  */
        !           817:
        !           818: GEN
        !           819: znstar_group(long n, GEN ZN, GEN H)
        !           820: {
        !           821:   gpmem_t ltop=avma;
        !           822:   int j,h;
        !           823:   GEN m=stoi(n);
        !           824:   GEN gen;
        !           825:   for (j = 1; j < lg(gen); j++)
        !           826:   {
        !           827:     gen[j] = 1;
        !           828:     for (h = 1; h < lg(lss); h++)
        !           829:       gen[j] = mulssmod(gen[j], itos(powmodulo((GEN)zn3[h],gmael(lss,j,h),m)),n);
        !           830:     ord[j] = zn2[j] / itos(gmael(lss,j,j));
        !           831:   }
        !           832:   avma=ltop;
        !           833:   return gen;
        !           834: }
        !           835: #endif
        !           836:
        !           837:
        !           838: GEN
        !           839: abelian_group(GEN v)
        !           840: {
        !           841:   GEN G=cgetg(3,t_VEC);
        !           842:   long card;
        !           843:   long i;
        !           844:   long d=1;
        !           845:   G[1]=lgetg(lg(v),t_VEC);
        !           846:   G[2]=lcopy(v);
        !           847:   card=group_order(G);
        !           848:   for(i=1;i<lg(v);i++)
        !           849:   {
        !           850:     GEN p=cgetg(1+card,t_VECSMALL);
        !           851:     long o=v[i],u=d*(o-1);
        !           852:     long j,k,l;
        !           853:     mael(G,1,i) = (long) p;
        !           854:     /*The following loop is a bit over-optimised. Oh well.
        !           855:      *Remember that I wrote the loop in testpermutation.
        !           856:      *Something have survived... BA*/
        !           857:     for(j=1;j<=card;)
        !           858:     {
        !           859:       for(k=1;k<o;k++)
        !           860:         for(l=1;l<=d;l++,j++)
        !           861:           p[j]=j+d;
        !           862:       for(l=1;l<=d;l++,j++)
        !           863:         p[j]=j-u;
        !           864:     }
        !           865:     d+=u;
        !           866:   }
        !           867:   return G;
        !           868: }

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