ClearAll["Global`*"]
Needs["Statistics`NormalDistribution`"]
Needs["Statistics`LinearRegression`"]
Needs["ERAEvaluate`"]
Needs["Calculus`Pade`"]
Off[LinearSolve::luc]
We set up preliminary values for B, m, nmax, and the precision and read in the necessary files.
B=31;
order=30;
file1 = FileName[B,"02"];
file2 = FileName[B,"10"];
numprecision = 34;
infile1 = OpenRead[file1];
infile2 = OpenRead[file2];
Table[Read[infile1, String], {7}];
Table[Read[infile2, String], {7}];
coeffs1 = Table[SetPrecision[Read[infile1,Number],numprecision], {order+2}];
coeffs2 = Table[SetPrecision[Read[infile2,Number],numprecision], {order+2}];
nmax = Min[Length[coeffs1],Length[coeffs2]];
m=-2;
delta = N[1/(2+2 Abs[m]),20];
bcoeffs and ccoeffs are the coefficients of the parameters b and c in Tim's thesis, Eq. (?). 
bcoeffs = Table[SetPrecision[coeffs1[[i]]+coeffs2[[i]],numprecision],
    {i,1,nmax}];
bnumberfunc = SeriesData[x,0,bcoeffs,0,Length[bcoeffs],1]//Normal;
bnumberfunc /. x -> delta;
ccoeffs=Table[0,{l,1,nmax}];
k=1;
Do[
   ccoeffs[[k]]=
       Sum[N[coeffs1[[k-j+1]],
       numprecision]*N[coeffs2[[j]],numprecision],{j,1,k}];
   k++,
   {nmax}];
cnumberfunc = SeriesData[x,0,ccoeffs,0,nmax,1]//Normal;
cnumberfunc /. x -> delta;
epsilon[i_,n_,delta_] := EconomizedRationalApproximation[cnumberfunc,
    {x,{-i,i},n,n}] /. x -> delta ;
N[Table[epsilon[0,n,delta], {n, 1, 15}], 20]

i=0.00000000000000000000000;
istep = 0.0005000000000000000000;
k0 = 0.3;
kstep = .01;
diff[n_] := Abs[energy[[n]]-energy[[n-1]]]-(Abs[energy[[n-1]]-energy[[n-2]]]);
monotone[n_] := Sign[energy[[n]]-energy[[n-1]]];
fit[k_] := N[Fit[energy, {1, Exp[-k x], Exp[-k x^2]}, x], numprecision];
reg[k_] := N[Regress[energy,{1,Exp[-k x], Exp[-k x^2]}, x][[2]][[2]], numprecision];
Do[
    energy=Table[N[epsilon[i,n,delta],numprecision], {n, 12, 15}];
    k=k0;
    reg2=0;
    While[ k <= k0, 
      k+=kstep;
      If[diff[4]< 0 && diff[3] < 0 && monotone[4] > 0 && reg > .999,
        Print[N[i,4]," ",k-kstep,"  ",reg[k],
            "  ",energy];
        Print["Converged value is ",Limit[fit[k], x -> Infinity]]
      ];
      reg2=reg[k];
    ];
    i+=istep,
{1000}];
b=\[InvisibleSpace]13.694639960;
c=44.320860160;
energyp[b_,c_] := (1/2)(b + Sqrt[b^2-4c])//N;
energym[b_,c_] := (1/2)(b - Sqrt[b^2-4c])//N;
energyp[b,c]
energym[b,c]
EnergyToBindingEnergy[energyp[b,c],B,m]
EnergyToBindingEnergy[energym[b,c],B,m]

BTildeToBeta[32,-2]//N