1 #include "defs.h"
2 #include "permfns.h"
3
4 extern char wrd,nt,isbase,inf[],outf1[],outf2[],fixed[];
5 extern int perm[],sv[],cp[],actgen[],orb[],
6 base[],lorb[],order[],pno[], *pptr[],*svptr[],
7 mp,mb,mnpt;
8 extern int psp,svsp;
9 int npt;
10 FILE *ip,*op;
11
12 /* The data structures in this program are similar to most permutation group
13 programs. Permutations are numbered 0,1,2,3,... (where 2x+1 is usually the
14 inverse of 2x) and stored in the array perm. Perm no x is pointed to by
15 pptr[x]. npt=no. of points. Usually pptr[i][npt+1] gives the number of the
16 group in the stabilizer chain in which the perm lies. So this no is i if the
17 perm fixes the first i-1 base points. nb=no of base points.
18 The base and the lengths of the orbits in the stab chain are stored in base
19 and lorb.
20 Schreier vectors are stored in sv, and pointed to by svptr[i], i=1,..,nb.
21 pno is a list of *pno (=pno[0]) perm nos.
22 cp is a similar list of length *cp, but it is used to represent the product
23 of the perms cp[1]cp[2]..cp[*cp]. This product can be evaluated by the
24 procedure image in permfns.c
25 */
26 int
gpprog(void)27 gpprog (void)
28 { int i,j,k,l,m,n,nperms,nb,jobt,np2,seek,cord,ocord,given,ordknown,trivrel,
29 lsv,u,v,w,x,y,z,mxp,mnb,bpt,bno,id;
30 int quot;
31 float grporder;
32
33 if ((ip=fopen(inf,"r"))== 0)
34 { fprintf(stderr,"Cannot open file %s\n",inf); return(-1); }
35 fscanf(ip,"%d%d%d%d",&npt,&nperms,&nb,&jobt);
36 if (npt>mnpt) { fprintf(stderr,"npt too big. Increase NPT.\n"); return(-1); }
37 if (jobt>0) { fprintf(stderr,"Wrong input format.\n"); return(-1); }
38 quot=psp/(npt+1); if (quot>mp) quot=mp; mxp=quot;
39 quot=svsp/npt; if (quot>mb) quot=mb; mnb=quot; if (mnb>=npt) mnb=npt-1;
40 if (nb>mnb)
41 { fprintf(stderr,"nb to big. Increase SVSP (or MB).\n"); return(-1); }
42 /* pptr[i] is the i th permutation, svptr[i] is the i the Schreier vector. */
43 for (i=0;i<mxp;i++) pptr[i]=perm+i*(npt+1)-1;
44 for (i=1;i<=mnb;i++) svptr[i]=sv+(i-1)*npt-1;
45
46 /* Next we read any base points and orbit lengths */
47 if (nb!=0)
48 { for (i=1;i<=npt;i++) orb[i]=0;
49 for (i=1;i<=nb;i++)
50 { fscanf(ip,"%d",base+i); if (orb[base[i]]!=0)
51 { fprintf(stderr,"Repeated base element.\n"); return(-1); }
52 orb[base[i]]=1;
53 }
54 }
55 if (jobt<0)
56 { jobt= -jobt; if (jobt>nb) {fprintf(stderr,"jobt too big.\n"); return(-1); }
57 seek=jobt; given=jobt; ordknown=1;
58 for (i=1;i<=jobt;i++) fscanf(ip,"%d",order+i);
59 }
60 else ordknown=0;
61
62 /* Now we read the generating permutations */
63 np2=2*nperms-2;
64 for (i=0;i<=np2;i+=2)
65 { k=i/2 +1; j=readperm(pptr[i]);
66 if (j==2)
67 { fprintf(stderr,"Generator no %d is not a permutation.\n",k); return(-1); }
68 if (j==1)
69 { fprintf(stderr,"Generator no %d is the identity.\n",k);
70 if (nt) return(-1);
71 nperms-=1; np2-=2; i-=2;
72 }
73 else
74 { invert(pptr[i],pptr[i+1]); actgen[i]=1; x=1;
75 for (z=base[x];x<=nb && pptr[i][z]==z;z=base[x]) x++;
76 pptr[i][npt+1]=x;
77 if (x>nb)
78 { if (isbase)
79 { fprintf(stderr,"Given base is not a base!\n"); return(-1);}
80 nb++; for (z=1;pptr[i][z]==z;z++); base[nb]=z;
81 printf("New base element no %d is %d.\n",nb,z);
82 }
83 if (x==1) printf("Generator no %d moves first base point.\n",k);
84 else printf("Generator no %d fixes first %d base point(s).\n",k,x-1);
85 }
86 if (nperms==0) { fprintf(stderr,"Trivial group!\n"); return(-1); }
87 }
88 fclose(ip);
89
90 if (wrd) {op=fopen(outf2,"w"); fprintf(op,"%4d\n",nperms); }
91 bno=nb; for (i=0;i<=mnb;i++) lorb[i]=0;
92
93 /* Now the main algorithm begins */
94 loop:
95 *pno=0; if (ordknown) ocord= (bno==nb) ? 1 : order[bno+1];
96 /* We make a list of the perm nos that fix the first bno-1 base pts */
97 for (i=0;i<=np2;i+=2)
98 { if (pptr[i][npt+1]>=bno && actgen[i]<=bno)
99 { (*pno)++; pno[*pno]=i; }
100 }
101 lorb[bno]=orbitsv(base[bno],svptr[bno],0);
102 if (ordknown)
103 { cord= (bno>=given) ? ocord*lorb[bno] : lorb[bno];
104 if (bno==seek && cord==order[bno])
105 { seek--; printf("Order is now as given for bno = %d.\n",bno);
106 bno--; if (bno==0) goto foundorder; goto loop;
107 }
108 }
109 /* Now we start generating the Schreier generators that fix the firat bno
110 base points, test them for membership, and add them as new generators
111 if necessary.
112 */
113 if (*pno!=0)
114 { nperms++; np2+=2; y=np2+1;
115 if (y>=mxp)
116 {fprintf(stderr,"Out of space. Increase PSP (or MP).\n"); return(-1); }
117 for (i=1;i<=npt;i++) fixed[i]=0;
118 for (i=1;i<bno;i++)
119 { u=base[i]; fixed[u]=1; pptr[np2][u]=u; pptr[y][u]=u; }
120 for (i=1;i<=lorb[bno];i++)
121 { *cp=0; addsv(orb[i],svptr[bno]);
122 for (w=1,x= *cp;w<=x;w++,x--)
123 { if (w==x) cp[w]-=1; else {u=cp[w]; cp[w]=cp[x]-1; cp[x]=u-1;}}
124 lsv= *cp;
125 for (j=1;j<=*pno;j++)
126 { *cp=lsv; u=pno[j]+1;
127 trivrel = (*cp>0) ? cp[*cp]==(u-1) : 0;
128 if (trivrel==0)
129 { (*cp)++; cp[*cp]=u; id=0;
130 for (l=bno;l<=nb;l++) fixed[base[l]]=0;
131 for (l=bno;l<=nb;l++)
132 { v=base[l]; u=image(v);
133 if (svptr[l][u]==0) goto newgen;
134 addsv(u,svptr[l]); pptr[np2][v]=v; pptr[y][v]=v; fixed[v]=1;
135 }
136 l=nb+1; id=1;
137 newgen: if (isbase==0)
138 for (k=1;k<=npt;k++)
139 { if (fixed[k]==0)
140 { u=image(k); pptr[np2][k]=u; pptr[y][u]=k;
141 if (id && k!=u)
142 { id=0; nb++;
143 if (nb>mnb) { fprintf(stderr,"nb to big. Increase SVSP (or MB).\n"); return(-1); }
144 base[nb]=k;
145 printf("New base point no %d is %d.\n",nb,k);
146 }
147 }
148 }
149 if (id==0)
150 { pptr[np2][npt+1]=l; actgen[np2]=bno+1;
151 if (isbase) for (k=1;k<=npt;k++)
152 { u=image(k); pptr[np2][k]=u; pptr[y][u]=k;}
153 printf("New generator no %d, fixing first %d base pts is:\n",
154 nperms,l-1);
155 for (m=1;m<=*cp;m++)
156 { z=cp[m]; y=z/2; x= (z==2*y) ? y+1 : -(y+1);
157 printf("%3d",x); if (m== *cp) printf("\n"); else printf("*");
158 }
159 if (wrd)
160 {for (m=0;m<=*cp;m++) fprintf(op,"%4d",cp[m]);fprintf(op,"\n");}
161 bno=l; goto loop;
162 }
163 }
164 }
165 }
166 nperms--; np2-=2;
167 }
168 if (ordknown && bno>given) order[bno]=cord;
169 bno--;
170
171 foundorder:
172 if (bno==0)
173 { printf("The order of the group is:\n");
174 for (i=1;i<=nb;i++)
175 { printf("%3d",lorb[i]);
176 if (i==nb) printf(" =\n"); else printf("*");
177 }
178 grporder=1; for (i=1;i<=nb;i++) grporder *= lorb[i];
179 printf("%g\n",grporder);
180 }
181 else goto loop;
182 if (wrd) { fprintf(op,"%d\n",-1); fclose(op); }
183
184 /* Now we output the generating perms and Schreier vectors */
185 op=fopen(outf1,"w");
186 fprintf(op,"%4d%4d%4d%4d\n",npt,nperms,nb,1);
187 printbaselo(nb,base,lorb); *pno=0;
188 for (i=1;i<=nperms;i++) {(*pno)++; pno[*pno]=2*(i-1);}
189 printpsv(nb,pno,svptr);
190 fclose(op); return(0);
191 }
192