1 #include "defs.h"
2 
3 /* This version contains the routine subact() called at end of each main
4    cycle in calcfm.
5 */
6 extern char act,inf[],inf1[],inf3[],outf1[],safilech;
7 extern short facexp,prime,exp,depth,*rpf,*rpb,*eexpnt,*enexpnt,**pcb,
8        dim,nng,onng,expnt[],nexpnt[],*pcptr[],**comptr[],*vec[],**mat[],cp[],
9        dpth[],pinv[],wt[],*wf,*wc,**extno,**subno,chsdim,chpdim,marg;
10 extern int ptrsp,wsp;
11 FILE *ip,*ips,*op;
12 
13 int
comprels(void)14 comprels (void)
15 /* Used when crel is true to compute values of relators of P in the chosen
16    stable extension of M by P.
17 */
18 { short i,j,k,l,m,n,k1,*p,*q,*r,stabdim,*stabhom,**stabno,*orpf,
19         nb,np,*rno,*covrel,**pgen,*ps,*pf,*v1,*v2,**cbmat;
20   char sgn;
21   stabdim=chpdim-chsdim; stabno=subno-stabdim;
22   stabhom=rpf-1; rpf+=onng; orpf=rpf;
23   for (i=1;i<=stabdim;i++) {stabno[i]=rpf-1; rpf+=onng; }
24   k=0;
25   for (i=1;i<=onng;i++) if (subno[i]==0 && extno[i]!=0)
26   { k++; p=stabno[k];
27     for (j=1;j<=onng;j++) p[j]=0; p[i]=1;
28     for (j=1;j<=onng;j++) if ((q=subno[j])!=0)
29     { expand(q,rpf,onng); if ((l=rpf[i])!=0) p[j]=prime-l; }
30   }
31   if (k!=stabdim) {fprintf(stderr,"stabdim error.\n"); return(-1);}
32   if (stabdim>1) fprintf(stderr,"Basis of stable homomorphisms:\n");
33   if (stabdim>1) for (i=1;i<=stabdim;i++)
34   { for (j=1;j<=onng;j++) fprintf(stderr,"%3d",stabno[i][j]); printf("\n"); }
35   for (i=1;i<=onng;i++) stabhom[i]=0;
36   if (stabdim>1)
37   fprintf(stderr,"Choose required stable hom as a vector in this basis!\n");
38   for (i=1;i<=stabdim;i++)
39   { if (stabdim>1) scanf("%hd",rpf);
40     else *rpf=1;
41     if (*rpf<0)
42     { for (j=1;j<=onng;j++) scanf("%hd",stabhom+j); break;}
43     if (*rpf!=0)
44     { p=stabhom; r=p+onng; q=stabno[i]+1;
45       while (++p<=r) { *p+= *q*(*rpf); *p%=prime; q++; }
46     }
47   }
48   printf("Chosen hom is:\n");
49   for (i=1;i<=onng;i++) printf("%3d",stabhom[i]); printf("\n");
50   rpf=orpf;
51   strcpy(inf1,inf); strcat(inf1,"psgwds"); ip=fopen(inf1,"r");
52   if (ip==0)
53   { fprintf(stderr,"Cannot open file %s.\n",inf1); return(-1); }
54   pgen=pcb; fscanf(ip,"%hd",&np);
55   for (i=1;i<=np;i++)
56   { pgen[i]=rpf; fscanf(ip,"%hd",rpf); p=rpf; rpf+=(1+ *p);
57     while (++p<rpf) fscanf(ip,"%hd",p);
58   }
59   fclose(ip);
60 /* Compute matrix to change basis of module back to the original, by
61    using the base change matrices output by matcalc.
62 */
63   strcpy(inf1,inf); strcat(inf1,"cbmats");
64   if ((ip=fopen(inf1,"r"))==0)
65   { fprintf(stderr,"Cannot open file %s.\n",inf1); return(-1); }
66   fscanf(ip,"%hd%hd%hd",&i,&j,&k);
67   if (i!=prime || j!=dim || k!=2)
68   { fprintf(stderr,"Error in line 1 of %s.\n",inf1); return(-1);}
69   readmat(mat[1]); readmat(mat[2]);
70   *cp=2; cp[1]=2; cp[2]=1; prod(cp,mat[3]);
71   inv(mat[3],mat[2]); trans(mat[2],mat[1]);
72   cbmat=mat[1]; v1=mat[2][1]; v2=mat[3][1];
73   fclose(ip);
74   strcpy(inf1,inf); strcat(inf1,"psg.rel"); ip=fopen(inf1,"r");
75   if (ip==0)
76   { fprintf(stderr,"Cannot open file %s.\n",inf1); return(-1); }
77   strcpy(outf1,inf); strcat(outf1,"psg.er"); op=fopen(outf1,"w");
78   fscanf(ip,"%hd",&nb); rno=rpf-1; rpf+=nb;
79   for (i=1;i<=nb;i++) fscanf(ip,"%hd",rno+i);
80   fprintf(op,"%4d%4d\n",nb,dim);
81   for (i=1;i<=nb;i++) fprintf(op,"%4d",rno[i]); fprintf(op,"\n");
82   wf=rpf;
83   for (i=1;i<=rno[1];i++)
84   { fscanf(ip,"%hd",&l); covrel=rpb-l; p=covrel;
85     while (++p<=rpb) fscanf(ip,"%hd",p);
86     zero(expnt,eexpnt);
87     for (j=l;j>=1;j--)
88     { wc=wf-2; k=covrel[j]; k1=k/2+1; m= *pgen[k1];
89       if (k%2==0) { sgn=1; ps=pgen[k1]+1; pf=ps+m-2; }
90       else { sgn= -1; pf=pgen[k1]+1; ps=pf+m-2; }
91       while(1)
92       { wc+=2; *wc= *ps; *(wc+1)= *(ps+1)*sgn;
93         if (ps==pf) break; ps+=(2*sgn);
94       }
95       collect(wc,wf,1);
96     }
97     fprintf(op,"%4d  ",l);
98     for (j=1;j<=l;j++) fprintf(op,"%4d",covrel[j]); fprintf(op,"\n");
99     zero(v1,v1+dim);
100     for (n=1;n<=exp;n++) if ((l=expnt[n])!=0)
101     { if (n<=facexp) {fprintf(stderr,"relation error.\n"); return(-1); }
102       for (j=1;j<=dim;j++)
103       { p= *(comptr[exp+j]+n);
104         if (p!=0)
105         { r=p+ *p;
106           while (++p<r) if ((k=stabhom[*p])!=0)
107           { v1[j]+=(k*l* *(++p)); v1[j]%=prime; }
108           else ++p;
109         }
110       }
111     }
112     im(v1,v2,cbmat);
113     l=0; p=v2;
114     while (++p<=v2+dim) if (*p !=0) l+=2;
115     fprintf(op,"%4d  ",l); p=v2;
116     while (++p<=v2+dim) if (*p!=0) fprintf(op,"%4d%4d",p-v2,*p);
117     fprintf(op,"\n");
118     if ((i==rno[1]) && (fscanf(ip,"%hd",&j)>0))
119     { fprintf(op,"%4d\n",j); rno[1]+=j;}
120   }
121   return(0);
122 }
123 
124 int
subact(void)125 subact (void)
126 { short substeps,subexp,*orpf,*commno,**commer,**comm,**subg,**pcp,subdp,
127   subdp0,ncommer,ncomm,commst,commlim,i,j,k,l,m,n,*ptr,*ptrlim,*ptre,*ptr2;
128   char subfile[3],infc[80],str[2];
129   bgc();
130   strcpy(subfile,"xa");
131   orpf=rpf;
132 restart:
133   pcp=pcb+1; rpf=orpf;
134   if (act)
135   { strcpy(infc,inf3); str[0]=safilech; str[1]='\0'; strcat(infc,str);}
136   else strcpy(infc,inf);
137   strcat(infc,subfile);
138   if ((ips=fopen(infc,"r"))==0) return(0);
139   fscanf(ips,"%hd%hd",&subexp,&substeps);
140   if (depth<=substeps) {fclose(ips); subfile[1]++; goto restart;}
141   subg=pcp-1; pcp+=facexp;
142   if (pcp-pcptr>ptrsp)
143   { printf("Not enough pointer space. Increase PTRSP.\n"); return(-1);}
144   subdp0=0; commst=facexp+1;
145   while (subexp==facexp)
146   { subdp0++; fscanf(ips,"%hd",&subexp);
147     if (subexp==0) subexp=facexp;
148     while (commst<=exp && dpth[commst]<=subdp0) commst++;
149   }
150 
151   for (i=1;i<=subexp;i++)
152   { subg[i]=rpf; fscanf(ips,"%hd",rpf); l= *rpf; rpf++;
153     for (j=1;j<=l;j++) {fscanf(ips,"%hd",rpf); rpf++;}
154   }
155   if (rpb+1-rpf<marg)
156   { printf("Running out of space in subact.\n");
157     if (wsp>marg) bgc(); else return(-1);
158   }
159   commer=pcp-1; commlim=commst; i=0;
160   while (commlim<=exp && dpth[commlim]<=depth-substeps+subdp0)
161   { *(pcp++)=rpf; *(rpf++)=2; *(rpf++)=commlim; *(rpf++)=1;
162     if (pcp-pcptr>ptrsp)
163     { printf("Not enough pointer space. Increase PTRSP.\n"); return(-1);}
164     i++; commlim++;
165   }
166   if (rpb+1-rpf<marg)
167   { printf("Running out of space.\n"); return(-1);}
168   ncommer=i;
169   for (subdp=subdp0+1;subdp<=substeps;subdp++)
170   { if (subdp<substeps)
171     { while (dpth[commst]<=subdp) commst++;
172       while (commlim<=exp && dpth[commlim]<=depth-substeps+subdp) commlim++;
173       commno=rpf-commst; rpf+= (commlim-commst);
174       if (rpb+1-rpf<marg)
175       { printf("Running out of space.\n"); return(-1);}
176       for (i=commst;i<commlim;i++) commno[i]=0;
177       comm=pcp-1; pcp+= (commlim-commst); ncomm=0;
178       if (pcp-pcptr>ptrsp)
179       { printf("Not enough pointer space. Increase PTRSP.\n"); return(-1);}
180     }
181 
182     for (i=1;i<=ncommer;i++) for (j=1;j<=subexp;j++)
183     { zero(expnt,eexpnt); zero(nexpnt,enexpnt);
184       wf=rpf; wc=wf-2;
185       enter(commer[i],-1);enter(subg[j],-1);enter(commer[i],1);enter(subg[j],1);
186       collect(wc,wf,1);
187       if (subdp==substeps)
188       { if (prnrel()== -1) return(-1);}
189       else
190       { ptr=expnt+commst-1; ptrlim=expnt+commlim;
191         while (++ptr<ptrlim) if ((l= *ptr)!=0)
192         { n=pinv[l]; k=ptr-expnt;
193           while (ptr<ptrlim) {*ptr*=n; *(ptr++)%=prime;}
194           if ((n=commno[k])!=0)
195           { ptr= comm[n]; l= *ptr; ptre=ptr+l;
196             while (++ptr<ptre)
197             { ptr2=expnt+ *ptr; *ptr2-= *(++ptr);
198               if (*ptr2<0) *ptr2+=prime;
199             }
200             ptr=expnt+k;
201           }
202           else
203           { ptr=expnt+k-1; ncomm++; commno[k]=ncomm; comm[ncomm]=rpf; l=0;
204             while (++ptr<ptrlim) if ((n= *ptr)!=0)
205             { l+=2; *(++rpf)=ptr-expnt; *(++rpf)=n;}
206             *comm[ncomm]=l; rpf++;
207             if (rpb+1-rpf<marg)
208             { printf("Running out of space in subact.\n");
209               if (wsp>marg) bgc(); else return(-1);
210             }
211             break;
212           }
213         }
214       }
215     }
216 
217     fscanf(ips,"%hd",&i);
218     if (i!=0)
219     { subexp=i; rpf=orpf;
220       if (subexp==facexp)
221       { k=0; for (i=1;i<=facexp;i++) if (wt[i]==1)
222         { k++; subg[k]=rpf; *(rpf++)=2; *(rpf++)=i; *(rpf++)=1;}
223         subexp=k;
224       }
225       else
226       for (i=1;i<=subexp;i++)
227       { subg[i]=rpf; fscanf(ips,"%hd",rpf); l= *rpf; rpf++;
228         for (j=1;j<=l;j++) {fscanf(ips,"%hd",rpf); rpf++;}
229       }
230       commer[1]=rpf;
231     }
232 
233     if (subdp<substeps)
234     {
235       rpf=commer[1]; ptr=comm[1];
236       for (i=1;i<=ncomm;i++)
237       { commer[i]=rpf; *(rpf++)=l= *(ptr++);
238         for (j=1;j<=l;j++) *(rpf++)= *(ptr++);
239       }
240       ncommer=ncomm; pcp=commer+ncommer+1;
241     }
242   }
243   printf("File %s.   ",infc);
244   printf("Reduced order at depth %d is: %d ^ %d\n",depth,prime,exp+nng);
245   fflush(stdout);
246   fclose(ips); subfile[1]++; goto restart;
247 }
248