In[]:=
ClearAll[expM,MRB1,MRB2,err,goodDigits];pr=2000;steps=3;expM[pre_Integer,alpha_?NumericQ]:=Module[{d,bb,cc,s,n,pr2,block,start,stop,rng,u,xvals,ctab,m,dot},pr2=Floor[1.005pre];n=Floor[alphapr2];block=Ceiling[n/steps];Print["Iterations required: ",n];d=N[ChebyshevT[n,3],pr2+50];bb=SetPrecision[-1,pr2+50];cc=-d;s=SetPrecision[0,pr2+50];start=1;While[start<=n,stop=Min[start+block-1,n];Print["Starting block ",start," to ",stop];(*rng=Range[start,stop];*)u=ParallelTable[N[Log[rng]/rng,pr2],{rng,start,stop}];(*documented,stablereplacementforExpMinus1*)xvals=N[Exp[u]-1,pr2];ctab=Table[cc=bb-cc;m=start+j-2;bb*=2(m+n)(m-n)/((m+1)(2m+1));cc,{j,1,stop-start+1}];dot=ctab.xvals;If[Dimensions[dot]=!={},Print["ERROR: dot is not scalar."];Print["Head[dot] = ",Head[dot]];Print["Dimensions[dot] = ",Dimensions[dot]];Abort[];];s+=dot;Print[stop," iterations done."];start=stop+1;];N[-s/d,pre]];Print["Computing MRB1..."];Print["The first run took this many seconds:",AbsoluteTiming[MRB1=expM[pr,1.32]][[1]]];Print["Computing MRB2..."];Print["the second run took this many seconds:",AbsoluteTiming[MRB2=expM[pr,1.34]][[1]]];If[ValueQ[MRB1]&&ValueQ[MRB2],err=N[Abs[MRB2-MRB1],50];Print["Error estimate = ",ScientificForm[err,20]];goodDigits=If[NumericQ[err]&&err>0,Floor[-Log10[err]],"at least "<>ToString[pr]];Print["Estimated good digits = ",goodDigits],Print["MRB1 or MRB2 failed; no valid error estimate."]];Print[pr," digits are ",MRB2]
Computing MRB1...
Iterations required: 2651
Starting block 1 to 884
884 iterations done.
Starting block 885 to 1768
1768 iterations done.
Starting block 1769 to 2651
2651 iterations done.
The first run took this many seconds:0.494573
Computing MRB2...
Iterations required: 2692
Starting block 1 to 898
898 iterations done.
Starting block 899 to 1796
1796 iterations done.
Starting block 1797 to 2692
2692 iterations done.
the second run took this many seconds:0.476123
Error estimate = 0.×
-2000
10
Estimated good digits = at least 2000
2000" digits are "0.18785964246206...2209597744091137