(* This program is for calculating the Q-value for any of the 175 modular functions *)
(* Schiff Research Project, 2002-2003, Joshua L. Wiczer *)
Remove["Global`*"]
(*Let md = The Label of the Modular Function you want to compute *)
md = 1A;
(*Let r=degree of numerator of the rational polynomial *)
r = 15;
(*Let s=degree of denominator of the rational polynomial *)
s=15;

n=300;
<<Coeff.txt
x=h[md];
lst={};
lst=Append[lst,x];
While[f[x]!=x, lst=Append[lst,f[x]]; x=f[x]];
lst=Sort[lst];
size=Length[lst];
For[et=1,et <=size, et=et+1,
    {x = lst[[et]],y=lst[[et-1]],
    If[et == 1, {
        c[k_]:=(c[(k/4)2+1]+((c[k/4])^2-c[k/4])/2+Sum[((c[j])*(c[2(k/4)-j])),{j,1,(k/4)-1}])/;(Mod[k,4]==0);
        c[k_]:=c[2((k-1)/4)+3]-(c[2])*(c[2((k-1)/4)])+((c[2((k-1)/4)])^2+c[2((k-1)/4)])/2+((c[((k-1)/4)+1])^2-c[((k-1)/4)+1])/2+Sum[((c[j])*c[2((k-1)/4)-j+2]),{j,1,((k-1)/4)}]+Sum[((c[j])*c[4((k-1)/4)-4j]),{j,1,((k-1)/4)-1}]+Sum[(((-1)^j)*(c[j])*c[4((k-1)/4)-j]),{j,1,2((k-1)/4)-1}]/;(Mod[k,4]==1);
        c[k_]:=(c[2((k-2)/4)+2]+Sum[(c[j]*c[2((k-2)/4)-j+1]),{j,1,(k-2)/4}])/;(Mod[k,4]==2);
        c[k_]:=(c[2((k-3)/4)+4]-(c[2])*(c[2((k-3)/4)+1])-(((c[2((k-3)/4)+1])^2)-(c[2((k-3)/4)+1]))/2+Sum[((c[j])*(c[2((k-3)/4)-j+3])),{j,1,((k-3)/4)+1}]+Sum[((c[j])*(c[4((k-3)/4)-4j+2])),{j,1,((k-3)/4)}]+Sum[(((-1)^j)*(c[j])*(c[4((k-3)/4)-j+2])),{j,1,2((k-3)/4)}])/;(Mod[k,4]==3);
        For[i=-1,i<=5,i=i+1,c[i]=d[x,i]];
        For[i=6,i<=n,i=i+1,l=c[i];c[i]=l; d[x,i]=l];
        Clear[c]}];
    If[et > 1, {
        c[k_]:=(c[(k/4)2+1]+((c[k/4])^2-e[k/4])/2+Sum[((c[j])*(c[2(k/4)-j])),{j,1,(k/4)-1}])/;(Mod[k,4]==0);
        c[k_]:=c[2((k-1)/4)+3]-(c[2])*(c[2((k-1)/4)])+((c[2((k-1)/4)])^2+e[2((k-1)/4)])/2+((c[((k-1)/4)+1])^2-e[((k-1)/4)+1])/2+Sum[((c[j])*c[2((k-1)/4)-j+2]),{j,1,((k-1)/4)}]+Sum[((e[j])*c[4((k-1)/4)-4j]),{j,1,((k-1)/4)-1}]+Sum[(((-1)^j)*(c[j])*c[4((k-1)/4)-j]),{j,1,2((k-1)/4)-1}]/;(Mod[k,4]==1);
        c[k_]:=(c[2((k-2)/4)+2]+Sum[(c[j]*c[2((k-2)/4)-j+1]),{j,1,(k-2)/4}])/;(Mod[k,4]==2);
        c[k_]:=(c[2((k-3)/4)+4]-(c[2])*(c[2((k-3)/4)+1])-(((c[2((k-3)/4)+1])^2)-(e[2((k-3)/4)+1]))/2+Sum[((c[j])*(c[2((k-3)/4)-j+3])),{j,1,((k-3)/4)+1}]+Sum[((e[j])*(c[4((k-3)/4)-4j+2])),{j,1,((k-3)/4)}]+Sum[(((-1)^j)*(c[j])*(c[4((k-3)/4)-j+2])),{j,1,2((k-3)/4)}])/;(Mod[k,4]==3);
        For[i=-1,i<=n,i=i+1,e[i]=d[y,i]];
        For[i=-1,i<=5,i=i+1,c[i]=d[x,i]];
        For[i=6,i<=n,i=i+1,l=c[i];c[i]=l; d[x,i]=c[i]]; Clear[c]; Clear[e];}];
    }]
t[q_]=(Sum[(d[x,i]*q^i),{i,-1,2s+5}])+(O[q])^(2s+6);
z[q_]=1/t[q];
dz1[q_]=q*D[z[q],q];
dz2[q_]=q*D[q*D[z[q],q],q];
dz3[q_]=q*D[q*D[q*D[z[q],q],q],q];
R[q_]=Sum[(Subscript[b,i])*(z[q])^i,{i,0,r}];
S[q_]=Sum[(Subscript[\[Theta],i])*(z[q])^i,{i,0,s}];
F[q_]=2*R[q]*((dz1[q])^2)*((dz1[q])^2)+(dz3[q])(dz1[q])*S[q]-(3/2)((dz2[q])^2)S[q];
F[q_]=Collect[Expand[F[q]],q];
Equations = {};
For[\[Tau]=2,\[Tau]<=r+s+1,\[Tau]=\[Tau]+1, {Equations=Append[Equations,Coefficient[F[q],q^(\[Tau])]==0]}];
solutions = Solve[Equations];
Simplify[((Sum[(Subscript[b,i])*(z^i),{i,0,r}])/(Sum[(Subscript[\[Theta],i])*(z^i),{i,0,s}]))/.solutions]