# Try to beat these MRB constant records!

GROUPS:
 Marvin Ray Burns 1 Vote POSTED BY: Marvin Ray Burns .I think this important point got buried near the end.When it comes to mine and a few more educated people's passion to calculate many digits and the dislike possessed by a few more educated people; it is all a matter telling us that the human mind is multi faceted in giving passion, to person a, for one task and to person b for another task!The MRB constant is defined below. See http://mathworld.wolfram.com/MRBConstant.html Here are some record computations. If you know of any others let me know.. On or about Dec 31, 1998 I computed 1 digit of the (additive inverse of the) MRB constant with my TI-92's, by adding 1-sqrt(2)+3^(1/3)-4^(1/4) as far as I could. That first digit by the way is just 0. On Jan 11, 1999 I computed 3 digits of the MRB constant with the Inverse Symbolic Calculator. In Jan of 1999 I computed 4 correct digits of the MRB constant using Mathcad 3.1 on a 50 MHz 80486 IBM 486 personal computer operating on Windows 95. Shortly afterwards I computed 9 correct digits of the MRB constant using Mathcad 7 professional on the Pentium II mentioned below. On Jan 23, 1999 I computed 500 digits of the MRB constant with the online tool called Sigma. In September of 1999, I computed the first 5,000 digits of the MRB Constant on a 350 MHz Pentium II with 64 Mb of ram using the simple PARI commands \p 5000;sumalt(n=1,((-1)^n*(n^(1/n)-1))), after allocating enough memory. On June 10-11, 2003 over a period, of 10 hours, on a 450mh P3 with an available 512mb RAM: I computed 6,995 accurate digits of the MRB constant. Using a Sony Vaio P4 2.66 GHz laptop computer with 960 MB of available RAM, on 2:04 PM 3/25/2004, I finished computing 8000 digits of the MRB constant. On March 01, 2006 with a 3GH PD with 2GBRAM available, I computed the first 11,000 digits of the MRB Constant. On Nov 24, 2006 I computed 40, 000 digits of the MRB Constant in 33hours and 26min via my own program in written in Mathematica 5.2. The computation was run on a 32-bit Windows 3GH PD desktop computer using 3.25 GB of Ram. Finishing on July 29, 2007 at 11:57 PM EST, I computed 60,000 digits of MRB Constant. Computed in 50.51 hours on a 2.6 GH AMD Athlon with 64 bit Windows XP. Max memory used was 4.0 GB of RAM. Finishing on Aug 3 , 2007 at 12:40 AM EST, I computed 65,000 digits of MRB Constant. Computed in only 50.50 hours on a 2.66GH Core2Duo using 64 bit Windows XP. Max memory used was 5.0 GB of RAM. Finishing on Aug 12, 2007 at 8:00 PM EST, I computed 100,000 digits of MRB Constant. They were computed in 170 hours on a 2.66GH Core2Duo using 64 bit Windows XP. Max memory used was 11.3 GB of RAM. Median (typical) daily record of memory used was 8.5 GB of RAM. Finishing on Sep 23, 2007 at 11:00 AM EST, I computed 150,000 digits of MRB Constant. They were computed in 330 hours on a 2.66GH Core2Duo using 64 bit Windows XP. Max memory used was 22 GB of RAM. Median (typical) daily record of memory used was 17 GB of RAM. Finishing on March 16, 2008 at 3:00 PM EST, I computed 200,000 digits of MRB Constant using Mathematica 5.2. They were computed in 845 hours on a 2.66GH Core2Duo using 64 bit Windows XP. Max memory used was 47 GB of RAM. Median (typical) daily record of memory used was 28 GB of RAM. Washed away by Hurricane Ike -- on September 13, 2008 sometime between 2:00PM - 8:00PM EST an almost complete computation of 300,000 digits of the MRB Constant was destroyed. Computed for a long 4015. Hours (23.899 weeks or 1.4454*10^7 seconds) on a 2.66GH Core2Duo using 64 bit Windows XP. Max memory used was 91 GB of RAM. The Mathematica 6.0 code used follows: Block[{$MaxExtraPrecision = 300000 + 8, a, b = -1, c = -1 - d, d = (3 + Sqrt[8])^n, n = 131 Ceiling[300000/100], s = 0}, a[0] = 1; d = (d + 1/d)/2; For[m = 1, m < n, a[m] = (1 + m)^(1/(1 + m)); m++]; For[k = 0, k < n, c = b - c; b = b (k + n) (k - n)/((k + 1/2) (k + 1)); s = s + c*a[k]; k++]; N[1/2 - s/d, 300000]] On September 18, 2008 a computation of 225,000 digits of MRB Constant was started with a 2.66GH Core2Duo using 64 bit Windows XP. It was completed in 1072 hours. Memory usage is recorded in the attachment pt 225000.xls, near the bottom of this post. . 250,000 digits was attempted but failed to be completed to a serious internal error which restarted the machine. The error occurred sometime on December 24, 2008 between 9:00 AM and 9:00 PM. The computation began on November 16, 2008 at 10:03 PM EST. Like the 300,000 digit computation this one was almost complete when it failed. The Max memory used was 60.5 GB. On Jan 29, 2009, 1:26:19 pm (UTC-0500) EST, I finished computing 250,000 digits of the MRB constant. with a multiple step Mathematica command running on a dedicated 64bit XP using 4Gb DDR2 Ram on board and 36 GB virtual. The computation took only 333.102 hours. The digits are at http://marvinrayburns.com/250KMRB.txt . The computation is completely documented in the attached 250000.pd at bottom of this post. On Sun 28 Mar 2010 21:44:50 (UTC-0500) EST, I started a computation of 300000 digits of the MRB constant using an i7 with 8.0 GB of DDR3 Ram on board.; But it failed due to hardware problems. I computed 299,998 Digits of the MRB constant. The computation began Fri 13 Aug 2010 10:16:20 pm EDT and ended 2.23199*10^6 seconds later | Wednesday, September 8, 2010. I used Mathematica 6.0 for Microsoft Windows (64-bit) (June 19, 2007) That is an average of 7.44 seconds per digit.. I used my Dell Studio XPS 8100 i7 860 @ 2.80 GH 2.80 GH with 8GB physical DDR3 RAM. Windows 7 reserved an additional 48.929 GB virtual Ram. I computed exactly 300,000 digits to the right of the decimal point of the MRB constant from Sat 8 Oct 2011 23:50:40 to Sat 5 Nov 2011 19:53:42 (2.405*10^6 seconds later). This run was 0.5766 seconds per digit slower than the 299,998 digit computation even though it used 16GB physical DDR3 RAM on the same machine. The working precision and accuracy goal combination were maximized for exactly 300,000 digits, and the result was automatically saved as a file instead of just being displayed on the front end. Windows reserved a total of 63 GB of working memory of which at 52 GB were recorded being used. The 300,000 digits came from the Mathematica 7.0 command Quit; DateString[] digits = 300000; str = OpenWrite[]; SetOptions[str, PageWidth -> 1000]; time = SessionTime[]; Write[str, NSum[(-1)^n*(n^(1/n) - 1), {n, \[Infinity]}, WorkingPrecision -> digits + 3, AccuracyGoal -> digits, Method -> "AlternatingSigns"]]; timeused = SessionTime[] - time; here = Close[str] DateString[] 314159 digits of the constant took 3 tries do to hardware failure. Finishing on September 18, 2012 I computed 314159 digits, taking 59 GB of RAM. The digits are came from the Mathematica 8.0.4 code DateString[] NSum[(-1)^n*(n^(1/n) - 1), {n, \[Infinity]}, WorkingPrecision -> 314169, Method -> "AlternatingSigns"] // Timing DateString[] Where I have 10 digits to round off. (The command NSum[(-1)^n*(n^(1/n) - 1), {n, [Infinity]}, WorkingPrecision -> big number, Method -> "AlternatingSigns"] tends to give about 3 digits of error to the right.)The following records are due to the work of Richard Crandall found here. Sam Noble of Apple computed 1,000,000 digits of the MRB constant in 18 days 9 hours 11 minutes 34.253417 seconds Finishing on Dec 11, 2012 Ricard Crandall, an Apple scientist, computed 1,048,576 digits in a lighting fast 76.4 hours. That's on a 2.93 Ghz 8-core Nehalem I computed a little over 1,200,000 digits of the MRB constant in 11 days, 21 hours, 17 minutes, and 41 seconds,( finishing on on March 31 2013). I used a six core Intel(R) Core(TM) i7-3930K CPU @ 3.20 GHz 3.20 GHz. On May 17, 2013 I finished a 2,000,000 or more digit computation of the MRB constant, using only around 10GB of RAM. It took 37 days 5 hours 6 minutes 47.1870579 seconds. I used a six core Intel(R) Core(TM) i7-3930K CPU @ 3.20 GHz 3.20 GHz. Finally, I would like to announce a new unofficial world record computation of the MRB constant that was finished on Sun 21 Sep 2014 18:35:06. It took 1 month 27 days 2 hours 45 minutes 15 seconds. I computed 3,014,991 digits of the MRB constant with Mathematica 10.0. I Used my new version of Richard Crandall's code, below, optimized for my platform and large computations. I also used a six core Intel(R) Core(TM) i7-3930K CPU @ 3.20 GHz 3.20 GHz with 64 GB of RAM of which only 16 GB was used. Can you beat it (in more number of digits, less memory used, or less time taken)? This confirms that my previous "2,000,000 or more digit computation" was actually accurate to 2,009,993 digits. (They were used as MRBtest2M.) (**Fastest (at MRB's end) as of 25 Jul 2014*.*) DateString[] prec = 3000000; (**Number of required decimals.*.*)ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 12, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.005 pre]}, chunksize = cores*tsize; n = Floor[1.32 pr]; end = Ceiling[n/chunksize]; Print["Iterations required: ", n]; Print["end ", end]; Print[end*chunksize]; d = ChebyshevT[n, 3]; {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l; x = N[E^(Log[ll]/(ll)), iprec]; pc = iprec; While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(**N[Exp[Log[ll]/ll], pr]**)x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 4]]; ctab = ParallelTable[Table[c = b - c; ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); c, {l, chunksize}], Method -> "EvaluationsPerKernel" -> 2]; s += ctab.(xvals - 1); start += chunksize; Print["done iter ", k*chunksize, " ", SessionTime[] - T0];, {k, 0, end - 1}]; N[-s/d, pr]]; t2 = Timing[MRBtest2 = expM[prec];]; DateString[] Print[MRBtest2] MRBtest2 - MRBtest2M t2 From the computation was {1.961004112059*10^6, Null}.Here are a couple of graphs of my record computations in max digits/ year: Attachments: Answer 3 years ago 62 Replies  Daniel Lichtblau 1 Vote It is hard to be certain that c1 and c2 are correct to 77 digits even though they agree to that extent. I'm not saying that they are incorrect and presumably you have verified this. Just claiming that whatever methods NSum may be using to accelerate convergence, there is really no guarantee that they apply to this particular computation. So c1 aand c2 could agree to that many places because they are computed in a similar manner without all digits actually being correct. Answer 3 years ago  Marvin Ray Burns 1 Vote I could really use your opinions here! If your getting tired of my posts let me know.Power towers are defined below. See http://mathworld.wolfram.com/PowerTower.html . ;;; How about breaking my record for computing the infinite power tower of the MRB constant. Here's how I have to go about computing n digits of it: I compute a decimal expansion of the MRB constant, which I save as m. Then I compute l = Log[m]; N[-ProductLog[-l]/l, n]. My record is 3 million digits, in which I used my 3,014,991 digit computation of the MRB constant (mentioned in the first post of this thread) for m.It would be EXTREMELY helpful if anyone could find a way to compute the exact power tower without first computing the MRB constant! Such a revelation would be synonymous with finding a closed for solution to the MRB constant! I've been trying to find the exact power tower, lately, with no success. Even a high precision approximation to the power tower would generate a nearly equally precise approximation to the MRB constant. (Let ll be the infinite power tower of m, then ll^(1/ll)=m)Attached are 3 million digits of the power tower. Attachments: Answer 3 years ago  The following might help anyone serious about breaking my record.Richard Crandall wrote about couple of new formulae for computing the MRB constant faster on pp 28 and 29 here. By the way, I wrote Crandall about formula (44) -- it seems it should have a minus sign in front of it -- and he wrote back: Answer 3 years ago  The following email Crandall sent me before he died might be useful for anyone checking their results: Answer 3 years ago  Perhaps some of these speed records will be easier to beat digits seconds sec./digit 500 0.062 0.000124 1000 0.1230 0.000123 3000 0.8640 0.000288 5000 2.4801 0.00049602 10,000 9.8515 0.00098515 20,000 47.6900668 0.0023845 30,000 121.3201698 0.00404401 40,000 233.3583473 0.00583396 50,193 363.7078029 0.00724619 100,493 2647.6114347 0.0263462 301,492 32075.94123 . 0.106391 The ultimate speed record to beat, though, is Crandall's 1,048,576 digits in a blistering fast 76.4 hours. Sam Noble of Apple computed 1,000,000 digits of the MRB constant in 18 days 9 hours 11 minutes 34.253417 seconds Answer 3 years ago  I found the following program Crandall wrote to check his code. (it gives the MRB constant to the desired level of precision.) It's hard to follow! But here it is: t1 = 1/2 (2 EulerGamma Log[2] - Log[2]^2); t2 = 1/24 (\[Pi]^2 Log[2]^2 - 2 \[Pi]^2 Log[ 2] (EulerGamma + Log[2] - 12 Log[Glaisher] + Log[\[Pi]]) - 6 (Zeta^\[Prime]\[Prime])[2]); FLAG = 1; prec = 200;$MaxPrecision = Infinity; expM[pr_] := Module[{a, d, s, k, b, c}, n = Floor[1.3 pr]; Print[n]; d = N[(3 + Sqrt[8])^n, prec]; d = 1/2 (d + 1/d); {b, c, s} = {-1, -d, 0}; T0 = T00 = SessionTime[]; Do[c = b - c; s = s + c*N[E^(Log[k + 1]/(k + 1)) - 1 - FLAG * (Log[k + 1]/(k + 1) + 1/2! Log[k + 1]^2/(k + 1)^2), pr]; (* s=s+c*N[(k+1)^(1/(k+1))-1,pr]; *) b = (k + n) (k - n) b/((k + 1) (k + 1/2)); If[Mod[k, 1000] == 0, Print[{k, SessionTime[] - T0}]; T0 = SessionTime[]], {k, 0, n - 1}]; N[-s/d, pr]]; (*Begin test to given precision.*) (*t1=Timing[MRBtrue=NSum[(-1)^n*(n^(1/n)-1),{n,1,Infinity},\ WorkingPrecision\[Rule]prec,Method\[Rule]"AlternatingSigns"];];*) t3 = Timing[MRBtest = expM[prec];] SessionTime[] - T00 MRBtest + FLAG * (t1 + t2)
3 years ago
 The MRB constantIt is easy to compute the MRB constant using the LHS, but how about computing it from the RHS, ?How many correct digits can you compute using the using it? Here's how one can compute 173 accurate digits using it. First compute known correct digits. Then compute using the LHS using NSumTerms -> 2^22 and compare. c1 = NSum[(-1)^n*(n^(1/n) - 1), {n, 1, Infinity}, Method -> "AlternatingSigns", WorkingPrecision -> 200, PrecisionGoal -> 200]; c2 = NSum[(2 k)^(1/(2 k)) - (2 k - 1)^(1/(2 k - 1)), {k, 1, Infinity}, WorkingPrecision -> 600, PrecisionGoal -> 600, NSumTerms -> 2^22]; c1 - c2 gives -6.6444030376101061743388604*10^-173While, if you try it naively, NSum[(-1)^n*(n^(1/n) - 1), {n, 1, \[Infinity]}, Method -> "AlternatingSigns", WorkingPrecision -> 100, PrecisionGoal -> 100] - N[Sum[(2 k)^(1/(2 k)) - (2 k - 1)^(1/(2 k - 1)), {k, 1, \[Infinity]}], 100] you get 8.65025879128817600921262272409653679780852599347660397129681915033728\ 91*10^-28, with only 28 digits of accuracy! Attachments:
3 years ago
 C1, the approximation to and a special case of , is computed correctly by Henri Cohen, Fernando Rodriguez Villegas, and Don Zagier,s Convergence Acceleration of Alternating Series as found next. Below it is computed by that method to at least 1000 decimals. (* Newer loop with Newton interior. *) prec = 1000;(*Number of required decimals.*) expM[pr_] := Module[{a, d, s, k, b, c}, n = Floor[1.32 pr]; Print["Iterations required: ", n]; d = N[(3 + Sqrt[8])^n, pr + 10]; d = Round[1/2 (d + 1/d)]; {b, c, s} = {-1, -d, 0}; T0 = SessionTime[]; Do[c = b - c; x = N[E^(Log[k + 1]/(k + 1)), iprec = Ceiling[prec/128]]; pc = iprec; Do[nprec = Min[2 pc, pr]; x = SetPrecision[x, nprec];(*Adjust precision*) x = N[x - x/(k + 1) + 1/x^k, nprec]; pc *= 2; If[nprec >= pr, Break[]], {ct, 1, 19}]; s += c*(x - 1); b *= 2 (k + n) (k - n)/((k + 1) (2 k + 1)); If[Mod[k, 1000] == 0, Print["Iterations: ", k, " Cumulative time (sec): ", SessionTime[] - T0];], {k, 0, n - 1}]; N[-s/d, pr]]; MRBtest2 = expM[prec] During evaluation of In[17]:= Iterations required: 1320During evaluation of In[17]:= Iterations: 0 Cumulative time (sec): 0.*10^-8During evaluation of In[17]:= Iterations: 1000 Cumulative time (sec): 0.1872004Out[18]= 0. 1878596424620671202485179340542732300559030949001387861720046840894772 3156466021370329665443310749690384234585625801906123137009475922663043 8929348896184120837336626081613602738126379373435283212552763962171489 3217020762820621715167154084126804483635416719985197680252759893899391 4457983505561350964852107120784442309586812949768852694956420425558648 3670441042527952471060666092633974834103115781678641668915460034222258 8380025455396892947114212218910509832871227730802003644521539053639505 5332203470627551159812828039510219264914673176293516190659816018664245 8249506972033819929584209355151625143993576007645932912814517090824249 1588320416906640933443591480670556469280678700702811500938060693813938 5953360657987405562062348704329360737819564603104763950664893061360645 5280675151935082808373767192968663981030949496374962773830498463245634 7931157530028921252329181619562697369707486576547607117801719578736830 0965902260668753656305516567361288150201438756136686552210674305370591 0397357561914891 c1 = NSum[(-1)^n*(n^(1/n) - 1), {n, 1, \[Infinity]}, Method -> "AlternatingSigns", WorkingPrecision -> 100, PrecisionGoal -> 100] 0.18785964246206712024851793405427323005590309490013878617200468408947 723156466021370329665443310750
3 years ago
 For c1 I used a method for alternating series to get a very precise sum approximation: That is from the "AlternatingSigns" option.
3 years ago
 One more record, Richard Crandall was the first one to compute 319,000 digits of the MRB constant in one day*, as early as Nov 24, 2012. Here is a portion of his email:*He may have been counting CPU time alone (from the Timing[] comand which is an under estimate of the actual time taken).. My personal records are in actual time taken.
3 years ago
 Above where Cohen said, the coefficients of the numerators are furthered below. I know the pattern for the first and last coefficient of each line, but can you figure out a simple pattern for the rest of them? My hat's off to you if you can! 2 16, 8 98, 80, 32 576, 544, 384, 128 3362, 3312, 2912, 1792, 512 19600, 19528, 18688, 15104, 8192, 2048 114242, 114144, 112576, 103168, 76288, 36864, 8192 665856, 665728, 663040, 641536, 557056, 376832, 163840, 32768 3880898, 3880736, 3876416, 3832064, 3603968, 2945024, 1826816, 720896, 131072
3 years ago
 OK, you might think I'm real vain for saying the following! Nonetheless, here you go. The ultimate record to break: Richard Mathar and Richard Crandall both wrote 1 scholarly article that talks about the MRB constant to a good extent. Can you break that record and publish a full article or two about the MRB constant? (You might have to come up with some new ideas.) Mathar's and Crandall's papers can be found on Google Scholar: https://scholar.google.com/scholar?q=%22MRB+constant%22 Steven R Finch, Eric W.Weisstein,and i have written some about it in other genres.
3 years ago
 What Richard Crandall and maybe others did to come up with that method is really good and somewhat mysterious. I still don't really understand the inner workings, and I had shown him how to parallelize it. So the best I can say is that it's really hard to compete against magic. (I don't want to discourage others, I'm just explaining why I myself would be reluctant to tackle this. Someone less familiar might actually have a better chance of breaking new ground.)In a way this should be good news. Should it ever become "easy" to compute, the MRB number would lose what is perhaps its biggest point of interest. It just happens to be on that cusp of tantalizingly "close" to easily computable (perhaps as sums of zeta function and derivatives thereof), yet still hard enough that it takes a sophisticated scheme to get more than a few dozen digits.
3 years ago
 Daniel Lichtblau, Which one do you not understand? how he derived the eta' formula: from ?Or is it how and why he added the extra loop to Cohen's method:ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); ?I've been asking for a proof of the eta formula derivation. I think it could become the foundation to a great paper.If the extra loop changes Cohen's method for to a method for I think the code in that loop could help prove the derivation.
3 years ago
 I can't say I understand either. My guess is the Eta stuff comes from summing (-1)^k*(Log[k]/k)^n over k, as those are the terms that appear in the double sum you get from expanding k^(1/k)-1 in powers of Log[k]/k (use k^(1/k)=Exp[Log[k]/k] and the power series for Exp). Even if it does come from this the details remain elusive..
3 years ago
 Daniel Lichtblau and others, Richard Crandall did intend to explian his work on the MRB constant and his program to compute it. When I wrote him with a possible small improvement to his program he said, "It's worth observing when we write it up." See screenshot:
3 years ago
 Crandall is not using his eta formulas directly!!!!!!! He computes Sum[(-1)^k*(k^(1/k) - 1), {k, 1, Infinity}] directly!Going back to Crandall's code: (*Fastest (at RC's end) as of 30 Nov 2012.*)prec = 500000;(*Number of \ required decimals.*)ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 4, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.02 pre]}, chunksize = cores*tsize; n = Floor[1.32 pr]; end = Ceiling[n/chunksize]; Print["Iterations required: ", n]; Print["end ", end]; Print[end*chunksize]; d = N[(3 + Sqrt[8])^n, pr + 10]; d = Round[1/2 (d + 1/d)]; {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l; x = N[E^(Log[ll]/(ll)), iprec]; pc = iprec; While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(*N[Exp[Log[ll]/ll], pr]*)x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 1]]; ctab = Table[c = b - c; ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); c, {l, chunksize}]; s += ctab.(xvals - 1); start += chunksize; Print["done iter ", k*chunksize, " ", SessionTime[] - T0];, {k, 0, end - 1}]; N[-s/d, pr]]; t2 = Timing[MRBtest2 = expM[prec];]; MRBtest2 - MRBtest3 x = N[E^(Log[ll]/(ll)), iprec]; Gives k^(1/k) to only 1 decimal place; they are either 1.0, 1.1, 1.2, 1.3 or 1.4 (usually 1.1 or 1.0).. On the other hand, While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));], takes the short precision x and gives it the necessary precision and accuracy for k^(1/k) (k Is ll there.) It actually computes k^(1/k). Then he remarks, "(N[Exp[Log[ll]/ll], pr])."After finding a fast way to compute k^(1/k) to necessary precision he uses Cohen's algorithm 1 (See a screenshot in a previous post.) to accelerate convergence of Sum[(-1)^k*(k^(1/k) - 1), {k, 1, Infinity}]. That is his secret!!As I mentioned in a previous post the "MRBtest2 - MRBtest3" is for checking with a known-to-be accurate approximation to the MRB constant, MRBtest3I'm just excited that I figured it out! as you can tell.
3 years ago
 Nice work. Worth a bit of excitement, I' d say.
3 years ago
 Richard Crandall might of had some help in developing his method. He wrote one time: "Marvin I am working on a highly efficient method for your constant, and I've been in touch with other mathematics scholars.Please be patient...recSent from my iPhone."
3 years ago
 The result of k^(1/k) from the iteration in Crandall's code grows like the following example of 123^(1/123).. Timing[ll = 123; x = 1 + 1/10; Table[Block[{$MaxExtraPrecision = Infinity}, x = SetPrecision[x, 20 + a^6]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll)); N[ll^(1/ll) - x, a + 20]], {a, 1, 13}]] That gives the following timing and accuracies per iteration. {8.829657, {-0.04239426367665370783, -0.02519017014329854597401, \ -0.0097174389047505410455275, -0.000933400633967433388167272, \ -9.456045144810702667437515*10^-7, \ -9.8570302819631224677862339*10^-16, \ -1.11649387558952636447488240*10^-42, \ -1.622508976729293063326569833*10^-123, \ -4.9794273187141078980478269830*10^-366, \ -1.43931637321260733551982880202*10^-1093, \ -3.476056847637409626605825363818*10^-3276, \ -4.8964202322309104109089811118496*10^-9824, \ -1.36852937628135955635109182098985*10^-29467}} Here is a more extreme example from my big computer where we get about 64,000,000 digits of accuracy of 123^(1/123) in only 20 iterations! In[47]:= Timing[ll = 123; x = 1 + 1/10; Table[Block[{$MaxExtraPrecision = Infinity}, x = SetPrecision[x, 20 + a^6]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll)); N[ll^(1/ll) - x, a + 20]], {a, 1, 20}]] Out[47]= {210.180147, {-0.04239426367665370783, \ -0.02519017014329854597401, -0.0097174389047505410455275, \ -0.000933400633967433388167272, -9.456045144810702667437515*10^-7, \ -9.8570302819631224677862339*10^-16, \ -1.11649387558952636447488240*10^-42, \ -1.622508976729293063326569833*10^-123, \ -4.9794273187141078980478269830*10^-366, \ -1.43931637321260733551982880202*10^-1093, \ -3.476056847637409626605825363818*10^-3276, \ -4.8964202322309104109089811118496*10^-9824, \ -1.36852937628135955635109182098985*10^-29467, \ -2.987998984519890677800153644718129*10^-88398, \ -3.1099929855547029560727524646233719*10^-265190, \ -3.50668133180579791324881778853196092*10^-795566, \ -5.026977910888395698424738167604170974*10^-2386694, \ -1.4809452514436132712292881458571004047*10^-7160077, \ -3.78647491015226246115034702337662410645*10^-21480228, 0.*10^-64000020}} Furthermore this process can compute 6,379,358,787 digits of the Sqrt[2], from an initial guess of machine precision, in 18 iterations: In[90]:= Timing[ll = 2; x = 2.^(1/2); Table[Block[{$MaxExtraPrecision = Infinity}, x = SetPrecision[x, 20 + a^7.5]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll)); N[ll^(1/ll) - x, a + 20]], {a, 1, 18}]] Out[90]= {1913.882668, {0.*10^-21, -1.800462644283455490473*10^-148, \ -7.2956225729359623776786*10^-445, \ -4.85397000370369097837099*10^-1334, \ -1.429556345255675471801415*10^-4001, \ -3.6518576944392011682818847*10^-12004, \ -6.08767627470915271802076570*10^-36012, \ -2.820100999650988078782491560*10^-108035, \ -2.8035222068968058313783728204*10^-324105, \ -2.75436831997714080172674040062*10^-972315, \ -2.612017346393607427831858236332*10^-2916945, \ -2.2276049962448401256762035120120*10^-8750835, \ -1.38173437723345455709717486529180*10^-26252505, \ -3.297491628267321222864512567028780*10^-78757516, \ -4.4818892212256882776648121052875660*10^-236272548, \ -1.12536490316593766234907240580868063*10^-708817643, 0.*10^-1073741813, 0.*10^-1073741813}} The above pattern is 3 times the previous number of accurate digits per iteration.. 708817643 * 3^2= 6379358787 I believe it will only take 191+16 = 207 iterations to compute a google digits of the Sqrt[2]: In[55]:= Log[10.^100/708817643]/Log[3] Out[55]= 191.04 Answer 3 years ago  How about computing the MRB constant from Crandall's eta derivative formulas? They are mentioned in a previous post but here they are again: I computed and checked 500 digits, using the first eta derivative formula in 38.6 seconds. How well can you do? Can you improve my program? (It is a 51.4% improvement of one of Crandall's programs.) I would like a little competition in some of these records! (That formula takes just 225 summands, compared to 10^501 summands using -1^(1/1)+2^(1/2)-3^(1/3)+.... See http://arxiv.org/pdf/0912.3844v3.pdf for more summation requirements for other summation methods.) In[37]:= mm = 0.187859642462067120248517934054273230055903094900138786172004684089\ 4772315646602137032966544331074969038423458562580190612313700947592266\ 3043892934889618412083733662608161360273812637937343528321255276396217\ 1489321702076282062171516715408412680448363541671998519768025275989389\ 9391445798350556135096485210712078444230958681294976885269495642042555\ 8648367044104252795247106066609263397483410311578167864166891546003422\ 2258838002545539689294711421221891050983287122773080200364452153905363\ 9505533220347062755115981282803951021926491467317629351619065981601866\ 4245824950697203381992958420935515162514399357600764593291281451709082\ 4249158832041690664093344359148067055646928067870070281150093806069381\ 3938595336065798740556206234870432936073781956460310476395066489306136\ 0645528067515193508280837376719296866398103094949637496277383049846324\ 5634793115753002892125232918161956269736970748657654760711780171957873\ 6830096590226066875365630551656736128815020143875613668655221067430537\ 0591039735756191489093690777983203551193362404637253494105428363699717\ 0244185516548372793588220081344809610588020306478196195969537562878348\ 1233497638586301014072725292301472333336250918584024803704048881967676\ 7601198581116791693527968520441600270861372286889451015102919988536905\ 7286592870868754254925337943953475897035633134403826388879866561959807\ 3351473990256577813317226107612797585272274277730898577492230597096257\ 2562718836755752978879253616876739403543214513627725492293131262764357\ 3214462161877863771542054231282234462953965329033221714798202807598422\ 1065564890048536858707083268874877377635047689160983185536281667159108\ 41219342016438600025850842655643500695483283012054619321661.\ 273833491444; In[30]:= Timing[ etaMM[m_, pr_] := Module[{a, d, s, k, b, c}, a[j_] := Log[j + 1]^m/(j + 1)^m; n = Floor[1.32 pr]; d = Cos[n ArcCos[3]]; {b, c, s} = {-1, -d, 0}; Do[c = b - c; s = s + c a[k]; b = (k + n) (k - n) b/((k + 1) (k + 1/2)), {k, 0, n - 1}]; N[s/d, pr] (-1)^m]; eta[s_] := (1 - 2^(1 - s)) Zeta[s]; eta1 = Limit[D[eta[s], s], s -> 1]; MRBtrue = mm; prec = 500; MRBtest = eta1 - Sum[(-1)^m etaMM[m, prec]/m!, {m, 2, Floor[.45 prec]}]; MRBtest - MRBtrue] Out[30]= {36.831836, 0.*10^-502} Here is a short table of computation times with that program: Digits Seconds 500 36.831836 1000 717.308198 1500 2989.759165 2000 3752.354453 I just now retweaked the program. It is now Timing[etaMM[m_, pr_] := Module[{a, d, s, k, b, c}, a[j_] := N[(-PolyLog[1, -j]/(j + 1))^m, pr]; n = Floor[1.32 pr]; d = Cos[n ArcCos[3]]; {b, c, s} = {-1, -d, 0}; Do[c = b - c; s = s + c a[k]; b = N[(k + n) (k - n) b/((k + 1) (k + 1/2)), pr], {k, 0, n - 1}]; Return[N[s/d, pr] (-1)^m]]; eta[s_] := (1 - 2^(1 - s)) Zeta[s]; eta1 = Limit[D[eta[s], s], s -> 1]; MRBtrue = mm; prec = 1500; MRBtest = eta1 - Sum[(-1)^m etaMM[m, prec]/Gamma[m + 1], {m, 2, Floor[.45 prec]}, Method -> "Procedural"]; MRBtest - MRBtrue] Here are my best eta derivative records: Digits Seconds 500 9.874863 1000 62.587601 1500 219.41540 2000 1008.842867 2500 2659.208646 3000 5552.902395 3500 10233.821601 That is using V10.0.2.0 Kernel. Here is a sample Timing[etaMM[m_, pr_] := Module[{a, d, s, k, b, c}, a[j_] := N[(-PolyLog[1, -j]/(j + 1))^m, pr]; n = Floor[1.32 pr]; d = Cos[n ArcCos[3]]; {b, c, s} = {-1, -d, 0}; Do[c = b - c; s = s + c a[k]; b = N[(k + n) (k - n) b/((k + 1) (k + 1/2)), pr], {k, 0, n - 1}]; Return[N[s/d, pr] (-1)^m]]; eta[s_] := (1 - 2^(1 - s)) Zeta[s]; eta1 = Limit[D[eta[s], s], s -> 1]; MRBtrue = mm; prec = 500; MRBtest = eta1 - Sum[(-1)^m etaMM[m, prec]/Gamma[m + 1], {m, 2, Floor[.45 prec]}]; ] N::meprec: Internal precision limit$MaxExtraPrecision = 50. reached while evaluating -Cos[660 ArcCos[3]]. N::meprec: Internal precision limit $MaxExtraPrecision = 50. reached while evaluating -Cos[660 ArcCos[3]]. N::meprec: Internal precision limit$MaxExtraPrecision = 50. reached while evaluating -Cos[660 ArcCos[3]]. General::stop: Further output of N::meprec will be suppressed during this calculation. Out[1]= {9.874863, Null}
3 years ago
 OOPS! I did believe for a little while that I figured out how to rapidly compute AND CHECK a computation of the MRB constant! but it wasn't quite a good check! (The timing given is in processor time [for computing and checking] only. T0 can be used with another SessionTime[] call at the end to figure out all time expired during running of the program.) I used both of Crandall's methods for computing it and used for a check, the nontrivial identity ,where gamma is the Euler constant and M is the MRB constant.Below is my first version of the code with results. If nothing else, I thought, the code pits Crandall's 2 methods against each other to show if one is wrong they both are wrong. (So it is not a real proof.) But these are two totally different methods! (the first of which has been proven by Henry Cohen to be theoretically correct here). For a second check mm is a known approximation to the constant; over 3 million non checked (as of now) digits are found in the attached file 3M.nb. (You will have to change the Format/Style to Input to use the digits.)I at least did hope that I did more than to only show == , but I'm afraid I didn't! In[15]:= (*MRB constant computation with verification! The constant's \ decimal approximation is saved as MRBtest*)prec = 5000;(*Number of \ required decimals.*)ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 4, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.02 pre]}, chunksize = cores*tsize; n = Floor[1.32 pr]; end = Ceiling[n/chunksize]; d = N[(3 + Sqrt[8])^n, pr + 10]; d = Round[1/2 (d + 1/d)]; {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l; x = N[E^(Log[ll]/(ll)), iprec]; pc = iprec; While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(*N[Exp[Log[ll]/ll], pr]*)x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 1]]; ctab = Table[c = b - c; ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); c, {l, chunksize}]; s += ctab.(xvals - 1); start += chunksize;;, {k, 0, end - 1}]; etaMs = N[-s/d - (EulerGamma Log[2] - Log[2]^2/2), pr]]; t2 = Timing[MRBtest2 = expM[prec];]; Print["The MRB constant was computed and checked to ", prec, " digits \ in ", t1 = t2[[1]] + Timing[eta[s_] := (1 - 2^(1 - s)) Zeta[s]; eta1 = Limit[D[eta[s], s], s -> 1]; MRBtrue = mm; MRBtest = eta1 + etaMs; check = MRBtest - MRBtrue][[1]], " seconds"]; check During evaluation of In[15]:= The MRB constant was computed and checked to 5000 digits in 2.12161 seconds Out[18]= 0.*10^-5000 In[19]:= MRBtest - mm Out[19]= 0.*10^-5000 Attachments:
3 years ago
 the identity in question is straightforward. Write n^(1/n) as Exp[Log[n]/n], take a series expansion at 0, and subtract the first term from all summands. That means subtracting off Log[n]/n in each summand. This gives your left hand side. We know it must be M - the sum of the terms we subtracted off. Now add all of them up, accounting for signs. Expand[Sum[(-1)^n*Log[n]/n, {n, 1, Infinity}]] (* Out[74]= EulerGamma Log[2] - Log[2]^2/2 *) So we recover the right hand side.I have not understood whether this identity helps with Crandall's iteration. One advantage it confers, a good one in general, is that it converts a conditionally convergent alternating series into one that is absolutely convergent. From a numerical computation point of view this is always good.
3 years ago
 Daniel Lichtblau and others, I just deciphered an Identity Crandall used for checking computations of the MRB constant just before he died. It is used in a previous post about checking, where I said it was hard to follow. The MRB constant is B here. B= In input form that is B= Sum[(-1)^(k + 1)*(-1 + (1 + k)^(1/(1 + k)) - Log[1 + k]/(1 + k) - Log[1 + k]^2/(2*(1 + k)^2)), {k, 0, Infinity}] + 1/24 (\[Pi]^2 Log[2]^2 - 2 \[Pi]^2 Log[ 2] (EulerGamma + Log[2] - 12 Log[Glaisher] + Log[\[Pi]]) - 6 (Zeta^\[Prime]\[Prime])[2]) + 1/2 (2 EulerGamma Log[2] - Log[2]^2) For 3000 digit numeric approximation, it is B=NSum[((-1)^( k + 1) (-1 + (1 + k)^(1/(1 + k)) - Log[1 + k]/(1 + k) - Log[1 + k]^2/(2 (1 + k)^2))), {k, 0, Infinity}, Method -> "AlternatingSigns", WorkingPrecision -> 3000] + 1/24 (\[Pi]^2 Log[2]^2 - 2 \[Pi]^2 Log[ 2] (EulerGamma + Log[2] - 12 Log[Glaisher] + Log[\[Pi]]) - 6 (Zeta^\[Prime]\[Prime])[2]) + 1/2 (2 EulerGamma Log[2] - Log[2]^2) It is anylitaclly straight forward too because Sum[(-1)^(k + 1)*Log[1 + k]^2/(2 (1 + k)^2), {k, 0, Infinity}] gives 1/24 (-\[Pi]^2 (Log[2]^2 + EulerGamma Log[4] - 24 Log[2] Log[Glaisher] + Log[4] Log[\[Pi]]) - 6 (Zeta^\[Prime]\[Prime])[2]) That is I wonder why he chose it?
3 years ago
 The new sum is this. Sum[(-1)^(k + 1)*(-1 + (1 + k)^(1/(1 + k)) - Log[1 + k]/(1 + k) - Log[1 + k]^2/(2*(1 + k)^2)), {k, 0, Infinity}] That appears to be the same as for MRB except now we subtract two terms from the series expansion at the origin of k^(1/k). For each k these terms are Log[k]/k + 1/2*(Log[k]/k)^2. Accounting for the signs (-1)^k and summing, as I did earlier for just that first term, we get something recognizable. Sum[(-1)^(k)*(Log[k]/(k) + Log[k]^2/(2*k^2)), {k, 1, Infinity}] (* Out[21]= 1/24 (24 EulerGamma Log[2] - 2 EulerGamma \[Pi]^2 Log[2] - 12 Log[2]^2 - \[Pi]^2 Log[2]^2 + 24 \[Pi]^2 Log[2] Log[Glaisher] - 2 \[Pi]^2 Log[2] Log[\[Pi]] - 6 (Zeta^\[Prime]\[Prime])[2]) *) So what does this buy us? For one thing, we get even better convergence from brute force summation, because now our largest terms are O((logk/k)^3) and alternating (which means if we sum in pairs it's actually O~(1/k^4) with O~ denoting the "soft-oh" wherein one drops polylogarithmic factors).How helpful is this? Certainly it cannot hurt. But even with 1/k^4 size terms, it takes a long time to get even 40 digits, let alone thousands. So there is more going on in that Crandall approach.
3 years ago
 If you copy and paste the code for Zeta''[2] you get a non evaluable input!I get the following results comparing two methods: In[655]:= Timing[MRBtrue = NSum[(-1)^n*(n^(1/n) - 1), {n, 1, Infinity}, WorkingPrecision -> prec, Method -> "AlternatingSigns"]][[1]] Out[655]= 8.049652 In[656]:= Timing[B = NSum[((-1)^(k + 1) (-1 + (1 + k)^(1/(1 + k)) - Log[1 + k]/(1 + k) - Log[1 + k]^2/(2 (1 + k)^2))), {k, 0, Infinity}, Method -> "AlternatingSigns", WorkingPrecision -> 3000] + 1/24 (\[Pi]^2 Log[2]^2 - 2 \[Pi]^2 Log[ 2] (EulerGamma + Log[2] - 12 Log[Glaisher] + Log[\[Pi]]) - 6 (Zeta''[2])) + 1/2 (2 EulerGamma Log[2] - Log[2]^2)][[1]] Out[656]= 12.277279
3 years ago
 Daniel Lichtblau,Notice Log[k]/(k) and Log[k]^2/(2k^2) are of the form Log[k]^n/(nk^n). We can take that summing of Log[k]/(k) + Log[k]^2/(2*k^2)+... to its end and write B= or compute NSum[(-1)^(k)*(k^(1/k) - 1 + Sum[Log[k]^n/(n*k^n), {n, 1, Infinity}]), {k, 1, Infinity}, WorkingPrecision -> 1000, Method -> "AlternatingSigns"] - NSum[(-1)^(1 + k) Log[(k - Log[k])/k], {k, 1, Infinity}, WorkingPrecision -> 1000, Method -> "AlternatingSigns"] Simplified,.B=,What does that do to our "better convergence from brute force summation?"P.S. I only used 2 terms in my sample data, so I might not have the correct pattern for Log[k]/(k) + Log[k]^2/(2*k^2)+... .The reason that I suspect that I've got the wrong pattern is, in the following code the summands Log[k]/(k) and Log[k]^2/(2k^2) give 0 where the Log[k]^n/(nk^n).for n>2 does not. (etaMM[m, n] Gives a numeric approximation for the mth derivative of m.) etaMM[m_, pr_] := Module[{a, d, s, k, b, c}, a[j_] := Log[j + 1]^m/(j + 1)^m; n = Floor[1.32 pr]; d = Cos[n ArcCos[3]]; {b, c, s} = {-1, -d, 0}; Do[c = b - c; s = s + c a[k]; b = (k + n) (k - n) b/((k + 1) (k + 1/2)), {k, 0, n - 1}]; Return[N[s/d, pr] (-1)^m]]; Table[-Sum[(-1)^m*etaMM[m, 30]/m!, {m, 1, a}] - Sum[(-1)^(k)*(Sum[Log[k]^n/(n*k^n), {n, 1, a}]), {k, 1, Infinity}], {a, 1, 15}]
3 years ago
 What about records of computing the integral analog of the MRB constant,: ?Richard Mathar did a lot of work on it here.I've gotten Matheamtica to compute 125 digits. However, they are not proven to be correct yet! They are 0.68765236892769436980931240936544016493963738490362254179507101010743\ 366253478493706862729824049846818873192933433546612328629 . First we compute the real part as far as Mathematica will allow. a1 = NIntegrate[Cos[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 100] 0.07077603931152880353952802183028200136575469620336299759658471973672\ 987938741600037745028756981434374 a2 = NIntegrate[Cos[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 120] a2 - a1 0.07077603931152880353952802183028200136575469620336302758317278266053\ 31986618615110244568060496758380620699811570793175408 2.998658806292380331927444551064700651847986149432*10^-53 a3 = NIntegrate[Cos[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 150] a3 - a2 0.07077603931152880353952802183028200136575469620336302758317278816361\ 8457264385970709799491401005081151056924116255307801983594127144525095\ 5653544005192 5.5030852586025244596853426853513292430889869429591759902612*10^-63 a4 = NIntegrate[Cos[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 200] a4 - a3 0.07077603931152880353952802183028200136575469620336302758317278816361\ 8457264382036580831881266177238210031756216795402920795214039271485948\ 634659563768084109747493815003439875479076850383786911941519465 -3.9341289676101348278429410251678994599048811883800878730391469306948\ 367511*10^-78 a5 = NIntegrate[Cos[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 250] a5 - a4 0.07077603931152880353952802183028200136575469620336302758317278816361\ 8457264382036580831881266177238209440733969109717926999044694539086929\ 3857095687266500964737783523859835124762555195276023702167529617039725\ 7261177753806842756198742365511173334813888 -5.9102224768568499379616934473239901924894999504143401327371546261745\ 6363002821330856184541724766503*10^-103 Next we compute the imaginary part to the same precision. b1 = NIntegrate[I Sin[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 100] - I/Pi 0.*10^-117 - 0.6840003894379321291827444599926611267109914826550016181302726087470\ 544306934833279937664708191960468 I b2 = NIntegrate[I Sin[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 120] - I/Pi b2 - b1 0.*10^-137 - 0.6840003894379321291827444599926611267109914826549994343226304054256\ 46767722886537984405858512438464223325361496951820797 I 0.*10^-117 + 2.1838076422033214076629705967900093606123067575826*10^-51 I b3 = NIntegrate[I Sin[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 150] - I/Pi b3 - b2 0.*10^-167 - 0.6840003894379321291827444599926611267109914826549994343226303771381\ 5305812568206208637713014270949108628424796532117557865488349026349505\ 4352728287677 I 0.*10^-137 + 2.8287493709597204475898028728369728973137041113531630645218*10^-62 I b4 = NIntegrate[I Sin[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 200] - I/Pi b4 - b3 0.*10^-218 - 0.6840003894379321291827444599926611267109914826549994343226303771381\ 5305812497663815095983421272147867735031056071869477552727290571462108\ 208123698276619850397331432861469605963724235550107655309644965 I 0.*10^-167 + 7.0542393541729592998801240893393740460248080312761058454887397227149\ 1304910*10^-76 I b5 = NIntegrate[I Sin[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 250] - I/Pi b5 - b4 0.*10^-268 - 0.6840003894379321291827444599926611267109914826549994343226303771381\ 5305812497663815095983421272147867223796451609148860995867828496814126\ 9810848570299802270095261060286697622600207986034863822997401942304753\ 4951409792726050072747412751162199808963072 I 0.*10^-218 + 5.1123460446272061655685946207464798122703884124663962338780532683279\ 9843703703436946621273009904771*10^-102 I b6 = NIntegrate[I Sin[x Pi] x^(1/x), {x, 1, Infinity}, WorkingPrecision -> 300] - I/Pi b6 - b5 0.*10^-318 - 0.6840003894379321291827444599926611267109914826549994343226303771381\ 5305812497663815095983421272147867223796451609148860995867804988314557\ 9408739051911924508290758754789975176921766748245229306743723292030351\ 1357229649514450909272015113199881208930542548540913212596310791355732\ 04151474091653439098975 I 0.*10^-268 + 2.3508499569040210951838787776180450230549672244567844123778963451625\ 36786502744023594180143211599163475397637962318600032530*10^-127 I Notice that WorkingPrecision->100 gave 51 consistant (correct) digits, WorkingPrecision->120 gave 62 correct digits, WorkingPrecision->150 gave 76 correct digits, WorkingPrecision->200 gave 102 correct digits, so it is not too much of a stretch to believe WorkingPrecision->250 gave 125 correct digits. In[78]:= c = N[Abs[a5 + b5], 125] Out[78]= 0.\ 6876523689276943698093124093654401649396373849036225417950710101074336\ 6253478493706862729824049846818873192933433546612328629
3 years ago
 Going back to the Crandall equations for the MRB constant: In particular, this time we will work on equation 44, which should have a negative sign in front.The naive approach is about as good as any other, so far:: Clear[a, c]; eta[s_] := (1 - 2^(1 - s)) Zeta[s]; a[i_] := Derivative[i][eta][0]; c[j_] := Sum[Binomial[j, d]*(-1)^d*d^(j - d), {d, 1, j}]; -N[Sum[c[m]/m!*a[m], {m, 1, 40}], 6] giving, after some warning, 0.187859 Nonetheless, I found a shortcut for computing the m'th derivatives of the eta function( with an argument of 0). For an example below you can replace 70 with any other n>3. Timing[Derivative[70][DirichletEta][0]]] can be computed in a lot less time by a = EulerGamma^2/2 - Pi^2/24 - (1/2)*(Log[2] + Log[Pi])^2 + StieltjesGamma[1]; b = EulerGamma^3 + (3/2)*EulerGamma^2*Log[2*Pi] - (1/8)*Pi^2* Log[2*Pi] - (1/2)*Log[2*Pi]^3 + 3*EulerGamma*StieltjesGamma[1] + 3*Log[2*Pi]*StieltjesGamma[1] + (3*StieltjesGamma[2])/2 - Zeta[3]; Timing[((-1)^n*Log[2]^n + (-1)^(n + 1)*n*Log[2]^(n - 1)* Log[2*Pi] + (-1)^(n + 1)*n*(n - 1)*Log[2]^(n - 2)*a) + (-1)^n*2*Binomial[n, 3]*Log[2]^(n - 3)* b - (Derivative[n][Zeta][0] + Sum[(-1)^x*2* Binomial[n, x]*(Derivative[n - x][Zeta][0]*Log[2]^x), {x, 1, n - 4}]) /. n -> 70] . That shortcut leads to the following raw program to compute the MRB constant. (Unfortunately, it is severely limited in its precision due to Mathematica's refusal to give high precision results for zeta derivatives of order >3 with argument of 0. For example try N[Derivative[4][Zeta][0], 40] ).Be careful reading the program. (There is a dot product in it. And it could use a lot of cleanup!) a = EulerGamma^2/2 - Pi^2/24 - (1/2)*(Log[2] + Log[Pi])^2 + StieltjesGamma[1]; b = EulerGamma^3 + (3/2)*EulerGamma^2*Log[2*Pi] - (1/8)*Pi^2* Log[2*Pi] - (1/2)*Log[2*Pi]^3 + 3*EulerGamma*StieltjesGamma[1] + 3*Log[2*Pi]*StieltjesGamma[1] + (3*StieltjesGamma[2])/2 - Zeta[3]; -N[With[{prec = 40}, (c = CoefficientList[Series[Exp[(-x)*Exp[x]], {x, 0, prec}], x])[[2 ;; 4]] . Table[Derivative[n][DirichletEta][0], {n, 1, 3}] + c[[5 ;; prec]] . Table[((-1)^n*Log[2]^n + (-1)^(n + 1)*n*Log[2]^(n - 1)* Log[2*Pi] + (-1)^(n + 1)*n*(n - 1)*Log[2]^(n - 2)*a) + (-1)^ n*2*Binomial[n, 3]*Log[2]^(n - 3)*b - (N[Derivative[n][Zeta][0], prec] + Sum[(-1)^x*2* Binomial[n, x]*(N[Derivative[n - x][Zeta][0], prec]*Log[2]^x), {x, 1, n - 4}]), {n, 4, prec - 1}]], 40] Does anyone care to work on this?The zeta derivatives of order n with argument of 0 do have exact forms not involving the zeta function.(I just haven't figured out the pattern to them.) For examples, Derivative[6][Zeta][0] is the same as Limit[D[2 (2 Pi)^(s - 1) Sin[Pi*s/2] Gamma[1 - s] Zeta[1 - s], {s, 6}], s -> 0] which gives the following that doesn't use zeta. (1/2688)*(6720*EulerGamma^6 + 5040*EulerGamma^4*Pi^2 + 1596*EulerGamma^2*Pi^4 - 275*Pi^6 + 32256*EulerGamma^5*Log[2*Pi] + 13440*EulerGamma^3*Pi^2*Log[2*Pi] + 60480*EulerGamma^4*Log[2*Pi]^2 + 10080*EulerGamma^2*Pi^2*Log[2*Pi]^2 - 1596*Pi^4*Log[2*Pi]^2 + 53760*EulerGamma^3*Log[2*Pi]^3 + 20160*EulerGamma^2*Log[2*Pi]^4 - 1680*Pi^2*Log[2*Pi]^4 - 1344*Log[2*Pi]^6 - 53760*EulerGamma^3*PolyGamma[2, 1] - 80640*EulerGamma^2*Log[2*Pi]*PolyGamma[2, 1] + 6720*Pi^2*Log[2*Pi]*PolyGamma[2, 1] + 26880*Log[2*Pi]^3*PolyGamma[2, 1] - 13440*PolyGamma[2, 1]^2 + 8064*Log[2*Pi]*PolyGamma[4, 1] + 40320*EulerGamma^4*StieltjesGamma[1] + 20160*EulerGamma^2*Pi^2*StieltjesGamma[1] + 3192*Pi^4*StieltjesGamma[1] + 161280*EulerGamma^3*Log[2*Pi]*StieltjesGamma[1] + 40320*EulerGamma*Pi^2*Log[2*Pi]*StieltjesGamma[1] + 241920*EulerGamma^2*Log[2*Pi]^2*StieltjesGamma[1] + 20160*Pi^2*Log[2*Pi]^2*StieltjesGamma[1] + 161280*EulerGamma*Log[2*Pi]^3*StieltjesGamma[1] + 40320*Log[2*Pi]^4*StieltjesGamma[1] - 161280*EulerGamma*PolyGamma[2, 1]*StieltjesGamma[1] - 161280*Log[2*Pi]*PolyGamma[2, 1]* StieltjesGamma[1] + 80640*EulerGamma^3*StieltjesGamma[2] + 20160*EulerGamma*Pi^2*StieltjesGamma[2] + 241920*EulerGamma^2*Log[2*Pi]*StieltjesGamma[2] + 20160*Pi^2*Log[2*Pi]*StieltjesGamma[2] + 241920*EulerGamma*Log[2*Pi]^2*StieltjesGamma[2] + 80640*Log[2*Pi]^3*StieltjesGamma[2] - 80640*PolyGamma[2, 1]*StieltjesGamma[2] + 80640*EulerGamma^2*StieltjesGamma[3] + 6720*Pi^2*StieltjesGamma[3] + 161280*EulerGamma*Log[2*Pi]*StieltjesGamma[3] + 80640*Log[2*Pi]^2*StieltjesGamma[3] + 40320*EulerGamma*StieltjesGamma[4] + 40320*Log[2*Pi]*StieltjesGamma[4] + 8064*StieltjesGamma[5]) In traditional form that result isEven Derivative[8][Zeta][0] is the same as Limit[D[2 (2 Pi)^(s - 1) Sin[Pi*s/2] Gamma[1 - s] Zeta[1 - s], {s, 8}], s -> 0] gives -(Pi^8/4608) + (1/6)* ... (an exact form not involving zeta). Derivative[n][Zeta][0] is approximately -n!Try In[656]:= ss = Table[N[Limit[ D[2 (2 Pi)^(s - 1) Sin[Pi*s/2] Gamma[1 - s] Zeta[1 - s], {s, n}], s -> 0] + n!, 20], {n, 1, 9}] Out[656]= {0.081061466795327258220, -0.0063564559085848512101, \ -0.0047111668622544477611, 0.0028968119862920410128, \ -0.00023290755845472453599, -0.00093682513005092950428, \ 0.00084982376500166915171, -0.00023243173551155958286, \ -0.00033058966361229644526}, so we can switch them in our ugly program to compute the MRB constant and get a rapidly converging rough approximation as follows. a = EulerGamma^2/2 - Pi^2/24 - (1/2)*(Log[2] + Log[Pi])^2 + StieltjesGamma[1]; b = EulerGamma^3 + (3/2)*EulerGamma^2*Log[2*Pi] - (1/8)*Pi^2* Log[2*Pi] - (1/2)*Log[2*Pi]^3 + 3*EulerGamma*StieltjesGamma[1] + 3*Log[2*Pi]*StieltjesGamma[1] + (3*StieltjesGamma[2])/2 - Zeta[3]; -N[With[{prec = 50}, (c = CoefficientList[Series[Exp[(-x)*Exp[x]], {x, 0, prec}], x])[[2 ;; 4]].Table[ Derivative[n][DirichletEta][0], {n, 1, 3}] + c[[5 ;; prec]].Table[((-1)^n*Log[2]^n + (-1)^(n + 1)*n* Log[2]^(n - 1)*Log[2*Pi] + (-1)^(n + 1)*n*(n - 1)* Log[2]^(n - 2)*a) + (-1)^n*2*Binomial[n, 3]*Log[2]^(n - 3)* b - (-n! + Sum[(-1)^x*2*Binomial[n, x]*(-(n - x)!*Log[2]^x), {x, 1, n - 4}]), {n, 4, prec - 1}]], 40] , which unfortunately is off from the true value of the MRB constant by about 0.00112931. Does anyone have suggestions on how to adjust this new program to give a better rapidly converging approximation to the MRB constant? or how to proceed from here?
3 years ago
 Using the same programs to compute the MRB constant, MMA V10.1 is a little slower than 10.0.0.Here is a sample of such programs: (*Fastest (at MRB's end) as of 24 Dec 2014.*) prec = 10000;(*Number of required decimals.*)ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 12, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.005 pre]}, chunksize = cores*tsize; n = Floor[1.32 pr]; end = Ceiling[n/chunksize]; Print["Iterations required: ", n]; Print["end ", end]; Print[end*chunksize]; d = Cos[n ArcCos[3]]; {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l; h = Log[ll]/ll; x = N[Exp[h], iprec]; pc = iprec; While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(*N[Exp[Log[ll]/ll], pr]*)x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 4]]; ctab = ParallelTable[Table[c = b - c; ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); c, {l, chunksize}], Method -> "EvaluationsPerKernel" -> 2]; s += ctab.(xvals - 1); start += chunksize; Print["done iter ", k*chunksize, " ", SessionTime[] - T0];, {k, 0, end - 1}]; N[-s/d, pr]]; t2 = Timing[MRBtest2 = expM[prec];]; MRBtest2 I'm in the process of tweaking that program to perform better in V 10.1.I tweaked the program a little. I computed 301493 digits in 31266 seconds, or 8 hours 41 minutes 6 seconds. (See attachment.) It's pretty hard to tweak without loosing accuracy!!!Some other 10.1 records: digits seconds 10,000 9.7 20,000 44 30,000 119 40,000 228 50,193 360 100,493 2573 301,493 31266 Attachments:
3 years ago
 Going back to integral analog of the MRB constant':Using formula 5 on page 3 of http://arxiv.org/pdf/0912.3844v3.pdf.We can compute a great deal of digits of the integral analog of the MRB constant' (I once called it the MKB constant, named after Marsha Kell-Burns my, now ex, wife.) In the paper Mathar simply calls it M1.Until further notice in this post when we compute the imaginary part of M1, we will be concerned with the imaginary part's absolute value only,This time we will compute the Imaginary part first to at least 500 digits: a[1] = 0; For[n = 1, n < 11, a[n] = N[2/Pi - 1/Pi*NIntegrate[ Cos[Pi*x]*x^(1/x)*(1 - Log[x])/x^2, {x, 1, Infinity}, WorkingPrecision -> 100*n], 50 n]; Print[a[n] - a[n - 1]], n++]; Print[a[11]] \ giving 0.6840003894379321291827444599926611267109914826549994343226303771381530581249766381509598342127214787 0.*10^-101 0.*10^-151 0.*10^-201 0.*10^-251 0.*10^-301 0.*10^-351 0.*10^-401 0.*10^-451 0.*10^-501 0.6840003894379321291827444599926611267109914826549994343226303771381530581249766381509598342127214786722379645160914886099586780498831455794087390519118879988351918366211827085883779918191195794251385436100844782462528597869421390620796113023053439642582325892202911183326091512210367124716901047132601108752764946385830438156754378694878046808312868541961166205744280461776232345922905313658259576212809654022016030244583148587352474339130505540080799774619683572540292971258866450201101870835703060314349396491402064932644813564545345219868887520120 . Likewise the real part: b[1] = 0; For[n = 1, n < 11, b[n] = N[-1/Pi* NIntegrate[Sin[Pi*x]*x^(1/x)*(1 - Log[x])/x^2, {x, 1, Infinity}, WorkingPrecision -> 100*n], 50 n]; Print[b[n] - b[n - 1]], n++]; Print[b[11]] giving 0.07077603931152880353952802183028200136575469620336302758317278816361845726438203658083188126617723821 0.*10^-102 0.*10^-152 0.*10^-202 0.*10^-252 0.*10^-302 0.*10^-352 0.*10^-402 0.*10^-452 0.*10^-502 0.07077603931152880353952802183028200136575469620336302758317278816361845726438203658083188126617723820944073396910971792699904464538475364292258443860652193330471222906120205483985764336623434898438270710499897053952312269178485299032185072743545220051257328105422174249313177670295863771714489658779291185716175115405623656039914848817528200250723061535734571065031458992196831648681239079549382556509741967588147362548743205919028695774572411439927516593391029992733107982746794845130889328251307263102570083031527430861023428334369104098217022622689 . Then the magnitude: N[Sqrt[a[11]^2 + b[11]^2], 500] giving 0.68765236892769436980931240936544016493963738490362254179507101010743\ 3662534784937068627298240498468188731929334335466123286287665409457565\ 9577211580255650416284625143925097120589697986500952590195706813170472\ 5387265069668971286335322245474865156721299946377659227025219748069576\ 0895993932096027520027641920489863095279507385793449828250341732295653\ 3809181101532087948181335825805498812728097520936901677028741356923292\ 2644964771090329726483682930417491673753430878118054062296678424687465\ 624513174205 . That checks with the 200 digits computed by the quadosc command in mpmath by FelisPhasma at https://github.com/FelisPhasma/MKB-Constant .The function is defined here: http://mpmath.googlecode.com/svn/trunk/doc/build/calculus/integration.html#oscillatory-quadrature-quadoscP.S.I just now finished 750 digits, (about the max with formula 5 from the paper, as far as Mathematica is concerned).Here is the work: a[1] = 0; For[n = 1, n < 16, a[n] = N[2/Pi - 1/Pi*NIntegrate[ Cos[Pi*x]*x^(1/x)*(1 - Log[x])/x^2, {x, 1, Infinity}, WorkingPrecision -> 100*n], 50 n]; Print[a[n] - a[n - 1]], n++]; Print[a[16]]; b[1] = 0; For[n = 1, n < 16, b[n] = N[-1/Pi* NIntegrate[Sin[Pi*x]*x^(1/x)*(1 - Log[x])/x^2, {x, 1, Infinity}, WorkingPrecision -> 100*n], 50 n]; Print[b[n] - b[n - 1]], n++]; Print[b[16]]; Print[N[Sqrt[a[16]^2 + b[16]^2], 750]] 0.6840003894379321291827444599926611267109914826549994343226303771381530581249766381509598342127214787 0.*10^-101 0.*10^-151 0.*10^-201 0.*10^-251 0.*10^-301 0.*10^-351 0.*10^-401 0.*10^-451 0.*10^-501 0.*10^-551 0.*10^-601 3.*10^-650 -4.*10^-700 -2.6*10^-749 0.68400038943793212918274445999266112671099148265499943432263037713815\ 3058124976638150959834212721478672237964516091488609958678049883145579\ 4087390519118879988351918366211827085883779918191195794251385436100844\ 7824625285978694213906207961130230534396425823258922029111833260915122\ 1036712471690104713260110875276494638583043815675437869487804680831286\ 8541961166205744280461776232345922905313658259576212809654022016030244\ 5831485873524743391305055400807997746196835725402929712588664502011018\ 7083570306031434939649140206493264481356454534521986888752011950353818\ 1776359577265099302389566135475579468144849763261779452665955246258699\ 8679271659049208654746533234375478909962633090080006358213908728990850\ 5026759549928935029206442637425786005036048098598304092996753145589012\ 64547453361707037686708654522699 0.07077603931152880353952802183028200136575469620336302758317278816361845726438203658083188126617723821 0.*10^-102 0.*10^-152 0.*10^-202 0.*10^-252 0.*10^-302 0.*10^-352 0.*10^-402 0.*10^-452 0.*10^-502 2.*10^-551 -1.*10^-600 1.8*10^-650 1.27*10^-699 4.34*10^-749 0.07077603931152880353952802183028200136575469620336302758317278816361\ 8457264382036580831881266177238209440733969109717926999044645384753642\ 9225844386065219333047122290612020548398576433662343489843827071049989\ 7053952312269178485299032185072743545220051257328105422174249313177670\ 2958637717144896587792911857161751154056236560399148488175282002507230\ 6153573457106503145899219683164868123907954938255650974196758814736254\ 8743205919028695774572411439927516593391029992733107982746794845130889\ 3282513072631025700830315274308610234283343691040982170226226904594029\ 7055093272952022662549075225941956559080574835998923469310063614655255\ 0629713179601483134045038416878054929072981851045829413286377842843667\ 5378730394247519728064887287780998671021887797977772522419765594172569\ 277490031071938177749184834961300 0.687652368927694369809312409365440164939637384903622541795071010107433662534784937068627298240498468188731929334335466123286287665409457565957721158025565041628462514392509712058969798650095259019570681317047253872650696689712863353222454748651567212999463776592270252197480695760895993932096027520027641920489863095279507385793449828250341732295653380918110153208794818133582580549881272809752093690167702874135692329226449647710903297264836829304174916737534308781180540622966784246874656245131742049004832216427665542900559350289936114782223424261285828326467186036500189315374147638489679365569122714398706519530651330568884655048857998738535162606116788633540389660052822237449082894798620397228331715198160243676576563833057235963591510865254600 Using formula 7 from http://arxiv.org/pdf/0912.3844v3.pdf, . (Treating it as we did formula 5), First, the imaginary part to at least 1000 digits:: a[1] = 0; For[n = 1, n < 21, a[n] = N[2/Pi + 1/Pi^2 NIntegrate[ Sin[x Pi] x^(1/x) (1 - 3 x + 2 (x - 1) Log[x] + Log[x]^2)/ x^4, {x, 1, Infinity}, WorkingPrecision -> 100 n], 50 n]; Print[a[n] - a[n - 1]], n++]; Print[a[21]] 0.6840003894379321291827444599926611267109914826549994343226303771381530581249766381509598342127214787 0.*10^-101 0.*10^-151 0.*10^-201 0.*10^-251 0.*10^-301 0.*10^-351 0.*10^-401 0.*10^-451 0.*10^-501 0.*10^-551 0.*10^-601 0.*10^-651 0.*10^-701 0.*10^-751 0.*10^-801 0.*10^-851 0.*10^-901 -2.*10^-950 5.*10^-1000 0.684000389437932129182744459992661126710991482654999434322630377138153058124976638150959834212721478672237964516091488609958678049883145579408739051911887998835191836621182708588377991819119579425138543610084478246252859786942139062079611302305343964258232589220291118332609151221036712471690104713260110875276494638583043815675437869487804680831286854196116620574428046177623234592290531365825957621280965402201603024458314858735247433913050554008079977461968357254029297125886645020110187083570306031434939649140206493264481356454534521986888752011950353818177635957726509930238956613547557946814484976326177945266595524625869986792716590492086547465332343754789099626330900800063582139087289908505026759549928935029206442637425786005036048098598304092996753145589012645474533617070376867086545228223060940434935219252885333298390272342234952870883304116640409421452765284609364941205344122569781634782508368641126766528707019957340895061936246645065753101916781254557006989818409283317145837167345971516970849116096077030635788389165381066055992688 Then the real part to at least 1000 digits: b[1] = 0; For[n = 1, n < 21, b[n] = N[1/Pi^2 - 1/Pi^2 NIntegrate[ Cos[Pi x] x^(1/x) (1 - 3 x + 2 (x - 1) Log[x] + Log[x]^2)/ x^4, {x, 1, Infinity}, WorkingPrecision -> 100 n], 50 n]; Print[b[n] - b[n - 1]], n++]; Print[b[21]] 0.07077603931152880353952802183028200136575469620336302758317278816361845726438203658083188126617723821 0.*10^-102 0.*10^-152 0.*10^-202 0.*10^-252 0.*10^-302 0.*10^-352 0.*10^-402 0.*10^-452 0.*10^-502 0.*10^-552 0.*10^-602 0.*10^-652 0.*10^-702 0.*10^-752 0.*10^-802 0.*10^-852 -3.*10^-901 8.*10^-951 -4.6*10^-1000 0.0707760393115288035395280218302820013657546962033630275831727881636184572643820365808318812661772382094407339691097179269990446453847536429225844386065219333047122290612020548398576433662343489843827071049989705395231226917848529903218507274354522005125732810542217424931317767029586377171448965877929118571617511540562365603991484881752820025072306153573457106503145899219683164868123907954938255650974196758814736254874320591902869577457241143992751659339102999273310798274679484513088932825130726310257008303152743086102342833436910409821702262269045940297055093272952022662549075225941956559080574835998923469310063614655255062971317960148313404503841687805492907298185104582941328637784284366753787303942475197280648872877809986710218877979777725224197655941725692774900310719381777491848349627938468198411955193898347075098152638657614980900350262780319142430252921925131515239611841070722530473939496294305264627977744876814858325335947117076721493110160508928494597906728688873533031986215124467678736429981544321187124269147141804397293341613 Then the magnitude: In[97]:= N[Sqrt[a[21]^2 + b[21]^2], 1000] Out[97]= 0.\ 6876523689276943698093124093654401649396373849036225417950710101074336\ 6253478493706862729824049846818873192933433546612328628766540945756595\ 7721158025565041628462514392509712058969798650095259019570681317047253\ 8726506966897128633532224547486515672129994637765922702521974806957608\ 9599393209602752002764192048986309527950738579344982825034173229565338\ 0918110153208794818133582580549881272809752093690167702874135692329226\ 4496477109032972648368293041749167375343087811805406229667842468746562\ 4513174204900483221642766554290055935028993611478222342426128582832646\ 7186036500189315374147638489679365569122714398706519530651330568884655\ 0488579987385351626061167886335403896600528222374490828947986203972283\ 3171519816024367657656383305723596359151086525460036387486837632622334\ 2987257095524637683005910353149353985736118868884201748241906260834981\ 7303422370398413326428269921074045506558966667483453656748906071577744\ 4147548424388220133662816274116986724576330176058912438027319979840883\ 05950589130911719199 PPS. I just now finished a 1500 digit computation of the integral analog of the MRB constant, but I don't have any way of checking it other than to see that it confirms smaller computations. Which thing it does. In[99]:= aa = N[2/Pi + 1/Pi^2 NIntegrate[ Sin[x Pi] x^(1/x) (1 - 3 x + 2 (x - 1) Log[x] + Log[x]^2)/ x^4, {x, 1, Infinity}, WorkingPrecision -> 3000], 1500] Out[99]= 0.\ 6840003894379321291827444599926611267109914826549994343226303771381530\ 5812497663815095983421272147867223796451609148860995867804988314557940\ 8739051911887998835191836621182708588377991819119579425138543610084478\ 2462528597869421390620796113023053439642582325892202911183326091512210\ 3671247169010471326011087527649463858304381567543786948780468083128685\ 4196116620574428046177623234592290531365825957621280965402201603024458\ 3148587352474339130505540080799774619683572540292971258866450201101870\ 8357030603143493964914020649326448135645453452198688875201195035381817\ 7635957726509930238956613547557946814484976326177945266595524625869986\ 7927165904920865474653323437547890996263309008000635821390872899085050\ 2675954992893502920644263742578600503604809859830409299675314558901264\ 5474533617070376867086545228223060940434935219252885333298390272342234\ 9528708833041166404094214527652846093649412053441225697816347825083686\ 4112676652870701995734089506193624664506575310191678125455700698981840\ 9283317145837167345971516970849116096077030635788389165381066055992708\ 4284702473154303800276803908560080204997803241058414188902018357202062\ 9532415382916822796942734253441520784640814155687968986766443021927163\ 6249354786973717955004441549085673392105556692081075647388204227896978\ 1483978754685921758294318270385312597177598977912650715548994562461701\ 1553879109152932039370312241134127950112036269188660519350584627066913\ 4925878278209048717316088629321353274101519307401594635990058104175474\ 300641475776727955287474213040 In[98]:= bb = N[1/Pi^2 - 1/Pi^2 NIntegrate[ Cos[Pi x] x^(1/x) (1 - 3 x + 2 (x - 1) Log[x] + Log[x]^2)/ x^4, {x, 1, Infinity}, WorkingPrecision -> 3000], 1500] Out[98]= 0.\ 0707760393115288035395280218302820013657546962033630275831727881636184\ 5726438203658083188126617723820944073396910971792699904464538475364292\ 2584438606521933304712229061202054839857643366234348984382707104998970\ 5395231226917848529903218507274354522005125732810542217424931317767029\ 5863771714489658779291185716175115405623656039914848817528200250723061\ 5357345710650314589921968316486812390795493825565097419675881473625487\ 4320591902869577457241143992751659339102999273310798274679484513088932\ 8251307263102570083031527430861023428334369104098217022622690459402970\ 5509327295202266254907522594195655908057483599892346931006361465525506\ 2971317960148313404503841687805492907298185104582941328637784284366753\ 7873039424751972806488728778099867102188779797777252241976559417256927\ 7490031071938177749184834962793846819841195519389834707509815263865761\ 4980900350262780319142430252921925131515239611841070722530473939496294\ 3052646279777448768148583253359471170767214931101605089284945979067286\ 8887353303198621512446767873642998154432118712426914714180439729334146\ 8345902382977472975053271988386946291215512340931334841526712825988330\ 6521193975174379922254198045615178994412133135553490942451521573377205\ 4086429300485891441696490339106907723915822537813700713422515725943626\ 7756749980892097547020923938358076198570370106085596863039832425037481\ 4946826330552459256977035009973219582010379262683780372730214991685800\ 3676611833579648850161974289307066295385292264148146789532534018500663\ 1153014589399140567464592864024 In[109]:= c1500 = Sqrt[aa^2 + bb^2] Out[109]= \ 0.68765236892769436980931240936544016493963738490362254179507101010743\ 3662534784937068627298240498468188731929334335466123286287665409457565\ 9577211580255650416284625143925097120589697986500952590195706813170472\ 5387265069668971286335322245474865156721299946377659227025219748069576\ 0895993932096027520027641920489863095279507385793449828250341732295653\ 3809181101532087948181335825805498812728097520936901677028741356923292\ 2644964771090329726483682930417491673753430878118054062296678424687465\ 6245131742049004832216427665542900559350289936114782223424261285828326\ 4671860365001893153741476384896793655691227143987065195306513305688846\ 5504885799873853516260611678863354038966005282223744908289479862039722\ 8331715198160243676576563833057235963591510865254600363874868376326223\ 3429872570955246376830059103531493539857361188688842017482419062608349\ 8173034223703984133264282699210740455065589666674834536567489060715777\ 4441475484243882201336628162741169867245763301760589124380273199798408\ 8305950589130911719198776146941477264898934365742508503405073273852990\ 3546587114217499635584514475429656959327732862489935076490012861232249\ 2446704232200904844779690044774489466704342791971033325818579375177198\ 9865742583276770011926585495711579480114327818546199372349313180236079\ 1389248808154759564302727311223193005229640892474022665093207969297797\ 9723087954832182561714039165214592519432072341006090867558444590500046\ 6707963346545638317950978935794173691635274461184852166407791838662429\ 40408834876470623546535579027725
3 years ago
 Mathar gives a simple scheme to find better formulas at http://arxiv.org/pdf/0912.3844v3.pdf . I could use some help in programming it: (I keep getting erroneous results!) Does anyone get the right results here?
3 years ago
 Below, where the upper limit of the following integrals shows Infinity, it is meant to be the (Ultraviolet limit of the sequence) as mentioned by Mathar here:Until further notice in this post when we compute the imaginary part of M1, we will be concerned with the imaginary part's absolute value only,I derived a new formula for computing the integral analog of the MRB constant': f[x_]:=x^(1/x);-((2 I)/\[Pi]^3) + 1/\[Pi]^2 - ( 2 I)/\[Pi] + (I/Pi)^3* Integrate[(-1)^x*D[f[x], {x, 3}], {x, 1, Infinity}] In traditional form that is M1= Using it I computed 2000 digits in only 10.8 minutes: In[131]:= Timing[f[x_] = x^(1/x); a = N[1/\[Pi]^2 + (1/Pi)^3* NIntegrate[Sin[Pi*x]*D[f[x], {x, 3}], {x, 1, Infinity}, WorkingPrecision -> 4000], 2000]; b = N[2/\[Pi]^3 + 2/\[Pi] + (1/Pi)^3* NIntegrate[Cos[Pi x]*D[f[x], {x, 3}], {x, 1, Infinity}, WorkingPrecision -> 4000], 2000]; Print[N[Sqrt[a^2 + b^2], 2000]]] During evaluation of In[131]:= 0.68765236892769436980931240936544016493963738490362254179507101010743366253478493706862729824049846818873192933433546612328628766540945756595772115802556504162846251439250971205896979865009525901957068131704725387265069668971286335322245474865156721299946377659227025219748069576089599393209602752002764192048986309527950738579344982825034173229565338091811015320879481813358258054988127280975209369016770287413569232922644964771090329726483682930417491673753430878118054062296678424687465624513174204900483221642766554290055935028993611478222342426128582832646718603650018931537414763848967936556912271439870651953065133056888465504885799873853516260611678863354038966005282223744908289479862039722833171519816024367657656383305723596359151086525460036387486837632622334298725709552463768300591035314935398573611886888420174824190626083498173034223703984133264282699210740455065589666674834536567489060715777444147548424388220133662816274116986724576330176058912438027319979840883059505891309117191987761469414772648989343657425085034050732738529903546587114217499635584514475429656959327732862489935076490012861232249244670423220090484477969004477448946670434279197103332581857937517719898657425832767700119265854957115794801143278185461993723493131802360791389248808154759564302727311223193005229640892474022665093207969297797972308795483218256171403916521459251943207234100609086755844459050004667079633465456383179509789357941736916352744611848521664077918386624294040883487647062354653558109265769644276994369741555722263494599492834558291937955573706480722982389806312472239746286527176248883116124285469947303667188075506826507811479428582807366599407544908560990699866167233307144245764835741501174979679166078765231145175411199825822532170091858833628202128777966026600647843068442894310401343003939117236867245656732686719139206716028255819141802331701942027248337771633882445225049334329008827371320849006472846226868011129149192754883153995560921671208059671732704499253517327447921147157 Out[131]= {653.145, Null} I am presently computing 10,000 digits using that formula. Come back here for results!That formula didn't work out; I will try one of the following formulas.Here are 2 more, more advanced formulas; remember f(x) is x^(1/x):I did finish a 5,000 digit computation using M1=in 48.11 minutes.Here are the 5000 digits:of the magnitude: 0.68765236892769436980931240936544016493963738490362254179507101010743366253478493706862729824049846818873192933433546612328628766540945756595772115802556504162846251439250971205896979865009525901957068131704725387265069668971286335322245474865156721299946377659227025219748069576089599393209602752002764192048986309527950738579344982825034173229565338091811015320879481813358258054988127280975209369016770287413569232922644964771090329726483682930417491673753430878118054062296678424687465624513174204900483221642766554290055935028993611478222342426128582832646718603650018931537414763848967936556912271439870651953065133056888465504885799873853516260611678863354038966005282223744908289479862039722833171519816024367657656383305723596359151086525460036387486837632622334298725709552463768300591035314935398573611886888420174824190626083498173034223703984133264282699210740455065589666674834536567489060715777444147548424388220133662816274116986724576330176058912438027319979840883059505891309117191987761469414772648989343657425085034050732738529903546587114217499635584514475429656959327732862489935076490012861232249244670423220090484477969004477448946670434279197103332581857937517719898657425832767700119265854957115794801143278185461993723493131802360791389248808154759564302727311223193005229640892474022665093207969297797972308795483218256171403916521459251943207234100609086755844459050004667079633465456383179509789357941736916352744611848521664077918386624294040883487647062354653558109265769644276994369741555722263494599492834558291937955573706480722982389806312472239746286527176248883116124285469947303667188075506826507811479428582807366599407544908560990699866167233307144245764835741501174979679166078765231145175411199825822532170091858833628202128777966026600647843068442894310401343003939117236867245656732686719139206716028255819141802331701942027248337771633882445225049334329008827371320849006472846226868011129149192754883153995560921671208059671732704499253517327447529208297180672654123457301218758892278525894167935930983363218877512533994251978272092700003994136520699813263053327399132641690231179063314931546906927612775633995348209911166678724589467821767106592498663827057034363632241807121831546175498178011687284590439293322231263406301066863589072717290630291441982684113819198880100231182613587798104863611185433976009254862585527222843445901958943153561148829083242874018226480554274231391324767376148485531787767908124831873688579979114662856184612164534836370699371440464263768724668291617743681719766849740663590277737977490693183461320266666793472116774276618408124767965369796362732668987556797338128876129264558867657737417548617146808592137056879602982206609613881069490166381528825180204703315896719667069923077454352649723496033985893188309150391579573916059639453655188856334980355047281560296288150836680499821806918067869468571687709518088408966653716009356556714281694904914038988996962213833530636987279769672200413448893419914190954063100962251649102614676944333201213024711868954772741991675045198246947499574872027800654821823797116399297131866662866832215332914761325880983081211272181775518951539503852063119472301382766303820851467743266039356123495461914463960644386394228342211998370152351720235034997434035743513051754761571835043769475528640144621307760159481496713401409374957729200400650100318226988524015127382509490642900236553851499823658269458873976032051355393161653806016080446394196719312454167915154602448638624354575153334932298393406734174580316934939632892851077461038399470015366439910136971186909599331204517462262508377673477745789645309425145559198802530351403897927622891172233239135167420567162398873965477371498335087310395422796362380227536212159184529243644094285328763286873653399867593200891823468738537356817916009007206857590792983184556882143118383332812491747733056313117179696094921120670802012310012864110800437831852620698327457619035904268498030693438632685623213366864129523404256345542376567721287706234359125016588483777876970236084456277023948551334490591022594253744077631232660869593809453087749830900393202787736482133628148979992109544954840067942735030391105496026321872468122542495017023785810605820545392820104069279893067324597299043883381251767370331206913429284614563732308018369972360638019778425246546329838131639355043236388708044857300408692365733932897876809202025693305332974091411983635619038514442263783801745983300121464879550146672827072002317686396598587702487509572349422593441184802476344187280014450860069307120621758277552124841158659386176036703247124389223327008210072318671884895179305778728051888524412158486781863155034447221379906386062559915129172725833420555901857729690605950941678587057025641848365090809750870051863842805803189784976076099574956436664131457150096711473033060684065060747340764998195621425524824611657787212347497307297184843276100338110267863618974154272345482369968216663233417338501929114697679974461999040589290327155974468087040862022522065912789 I'm getting closer to 10K digits of M1: Using , where f(x)=x^(1/x).I got approx. 10K digits of the imaginary part, but the real part was a little garbled. Finally using , where f(x) = x^(1/x) ,.I got about 10000 digits of M1 in about 12 hours. (It showed that the 5000 digit computation was only correct to 4979 digits, though.) Here is a rough program to get it: f[x_]:=x^(1/x); Print[DateString[]]; Print[T0 = SessionTime[]]; prec = 10000; Timing[Print[ a = N[Re[-(136584/Pi^10) - (34784*I)/Pi^9 + 5670/Pi^8 + (786*I)/Pi^7 - 90/Pi^6 - (4*I)/Pi^5 - 3/Pi^4 - (2*I)/Pi^3 + 1/Pi^2 - (2*I)/Pi] - (1/Pi)^10* NIntegrate[Cos[Pi*x]*D[f[x], {x, 10}], {x, 1, Infinity}, WorkingPrecision -> prec, PrecisionGoal -> prec], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[b = -Im[-(136584/Pi^10) - (34784*I)/Pi^9 + 5670/Pi^8 + (786*I)/Pi^7 - 90/Pi^6 - (4*I)/Pi^5 - 3/Pi^4 - (2*I)/Pi^3 + 1/Pi^2 - (2*I)/Pi] + (1/Pi)^10* NIntegrate[Sin[Pi*x]*D[f[x], {x, 10}], {x, 1, Infinity}, WorkingPrecision -> prec, PrecisionGoal -> prec], prec]]]; Print[ SessionTime[] - T0, " seconds"]; Print[N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; See attached 10000MKB.pdf and 10KMKB.nb for work and digits.On May 5, I computed another 10,000 digits in 9.55 hours see attached faster10KMKB.On May 6, I computed another 10,000 digits in a blistering fast 5.1 hours see attached fastest10KMKB.nb.On May 9, I improved that timing to 4.8 hours (17355 seconds). Here is the code I used: d = 15; f[x_] = x^(1/x); ClearAll[a, b, h]; h[n_] := Sum[ StirlingS1[n, k]*Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; h[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}]; Print[ DateString[]]; Print[T0 = SessionTime[]]; prec = 10000; Print[N[a = -Re[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = Im[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; Attachments:
3 years ago
 FelisPhasma has been helpful in providing me with a little competition in computing the integral analog of the MRB constant. See https://github.com/FelisPhasma/MKB-Constant.I've never done this before. But I so much would like to see others breaking these records that I'm going to give away a program that is practically guaranteed to break my record of 10,000 digits, for the integral analog of the MRB constant in a day or so. The program could use some "clean up" if you care to go that far. (The imaginary part is given as a positive, real constant: it actually starts with a negative sign and of course ends with I.)Here it is: f[x_] = x^(1/x); Print[DateString[]]; Print[ T0 = SessionTime[]]; prec = 11000; Timing[ Print[a = N[Re[(633666648 I)/\[Pi]^13 - 33137280/\[Pi]^12 - ((824760 I)/\[Pi]^11) - 136584/\[Pi]^10 - (34784 I)/\[Pi]^9 + 5670/\[Pi]^8 + (786 I)/\[Pi]^7 - 90/\[Pi]^6 - (4 I)/\[Pi]^5 - 3/\[Pi]^4 - (2 I)/\[Pi]^3 + 1/\[Pi]^2 - (2 I)/\[Pi]] + (1/Pi)^12* NIntegrate[Cos[Pi x]*D[f[x], {x, 12}], {x, 1, Infinity}, WorkingPrecision -> prec, PrecisionGoal -> prec], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = -Im[(633666648 I)/\[Pi]^13 - 33137280/\[Pi]^12 - ((824760 I)/\[Pi]^11) - 136584/\[Pi]^10 - (34784 I)/\[Pi]^9 + 5670/\[Pi]^8 + (786 I)/\[Pi]^7 + 90/\[Pi]^6 - (4 I)/\[Pi]^5 - 3/\[Pi]^4 - (2 I)/\[Pi]^3 + 1/\[Pi]^2 - (2 I)/\[Pi]] - (1/Pi)^13* NIntegrate[Cos[Pi x]*D[f[x], {x, 13}], {x, 1, Infinity}, WorkingPrecision -> prec, PrecisionGoal -> prec], prec]]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; Will anyone let me know you are running this program to break my record?Edit: On Sat 2 May 2015 19:03:45 I started a 15,000 digit, new record computation of the real and imaginary parts and magnitude of the integral analog of the MRB constant, (where the imaginary part is given as a positive, real constant), using the following code. f[x_] = x^(1/x); ClearAll[a]; h[n_] := Sum[ StirlingS1[n, k]*Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; h[0] = 1; g = -2 I/Pi + Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, 18}]; Print[DateString[]]; Print[T0 = SessionTime[]]; prec = 15000; Print[N[a = Re[g] + (1/Pi)^19* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, 19}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = -Im[g] + (1/Pi)^19* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, 19}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; The formula behind this computation is Edit: The program took 33.75 hours,The full run is attached in 15KMKB3.nb.Edit May 9, 2015: I better than halved my time! I computed 15000 digits in 14.83 hours. See fastestMKB15K.nb/. The faster formula isIf you still want me to write out a code for more digits, for you to break that record, let me know. Attachments:
3 years ago
 Still talking about the integral analog of the MRB constant:Here are my speed records -- can you beat any of them?Here is a graph of those speed records with a trendline:The 20,000 digit run is attached as MKB20K.nb, and MKB20K.pdf, Here is the algorithm used:Here is the code: d = 30; f[x_] = x^(1/x); ClearAll[a, b, h]; a[n_] := Sum[ StirlingS1[n, k]*Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; a[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) a[n]/Pi^(n + 1), {n, 1, d}]; Print[ DateString[]]; Print[T0 = SessionTime[]]; prec = 20000; Print[N[a = -Re[g] - (1/Pi)^(d)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = Im[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; I just now completed a 25,000 digit computation. It took 63.7 hours and confirmed the 20,000 digits. I updated MKB20K.nb and MKB20K.pdf. Here is the algorithm and the code I used: d = 35; f[x_] = x^(1/x); ClearAll[a, b, h]; h[n_] := Sum[ StirlingS1[n, k]*Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; h[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}]; Print[ DateString[]]; Print[T0 = SessionTime[]]; prec = 25000; Print[N[a = -Re[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = Im[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; Here is new a graph of those speed records with a trendline: Edit:On Tue 26 May 2015 06:21:00, I started a 30,000 digit computation using the following code. Does any one else want to try to break the record? $MaxExtraPrecision = 100; d = 43; f[x_] = x^(1/x); ClearAll[a, b, h]; h[n_] := Sum[ StirlingS1[n, k]*Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; h[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}]; Print[ DateString[]]; Print[T0 = SessionTime[]]; prec = 30000; Print[N[a = -Re[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = Im[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]]; Edit: My first full 30,0000 run finished on Sun 31 May 2015 00:45:09.Time span: {"4.767 days", "114.4 hours", "6864 minutes", "411849 seconds"} See attached MKB30K2.nb worksheet.Here is an updated speed record plot, with trendline. (I think the 30,000 digit run can be done faster.) Attachments: Answer 3 years ago I think I came up with a rough program that computes any "prec" digits of the integral analog of the MRB constant. It chooses, d, the best (or close to the best) order of derivative to use in Mathar's algorithm mentioned in a previous post (formula (12) at http://arxiv.org/pdf/0912.3844v3.pdf ), Then uses the appropriate code that integrates the integral analog of the constant. It shows the real and imaginary parts as postive real constants and the value the integral, and gives some timings. It could use a lot of cleanup! I hope someone can help me test it with varying values of prec. Please reply your intentions to use it and results.If no one else can clean it up I will after I tested it more. prec = 2000; d = Ceiling[0.264086 + 0.00143657 prec]; If[ Mod[d, 4] == 0, f[x_] = x^(1/x); ClearAll[a, b, h]; a[n_] := Sum[StirlingS1[n, k]* Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; a[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) a[n]/Pi^(n + 1), {n, 1, d}]; Print[DateString[]]; Print[T0 = SessionTime[]]; Print[N[a = -Re[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[b = Im[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]];, If[Mod[d, 4] == 1, f[x_] = x^(1/x); ClearAll[a, b, h]; h[n_] := Sum[StirlingS1[n, k]* Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; h[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}]; Print[DateString[]]; Print[T0 = SessionTime[]]; Print[N[ a = -Re[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[ b = Im[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]];, If[Mod[d, 4] == 2, f[x_] = x^(1/x); ClearAll[a, b, h]; a[n_] := Sum[StirlingS1[n, k]* Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; a[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) a[n]/Pi^(n + 1), {n, 1, d}]; Print[DateString[]]; Print[T0 = SessionTime[]]; Print[N[ a = -Re[g] - (1/Pi)^(d)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[ b = Im[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]];, If[Mod[d, 4] == 3, f[x_] = x^(1/x); ClearAll[a, b, h]; h[n_] := Sum[StirlingS1[n, k]* Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; h[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}]; Print[DateString[]]; Print[T0 = SessionTime[]]; Print[ N[a = -Re[g] + (1/Pi)^(d + 1)* NIntegrate[ Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[ N[b = Im[g] - (1/Pi)^(d + 1)* NIntegrate[ Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity}, WorkingPrecision -> prec*(105/100), PrecisionGoal -> prec*(105/100)], prec]]; Print[SessionTime[] - T0, " seconds"]; Print[N[Sqrt[a^2 + b^2], prec]]; Print[DateString[]];]]]] Here are some of my best timings to compare with the program's results: digits seconds 1000 38.8650545 2000 437.4906125 3000 889.473875 4000 1586.000714 5000 2802.591704 6000 4569.41586 7000 6891.057587 8000 9659.318566 9000 13491.43967 10000 17355 11000 12000 13000 14000 15000 53385.02323 16000 17000 18000 19000 20000 123876.4331 21000 22000 23000 24000 25000 229130.3088 26000 27000 28000 29000 30000 411848.6322 Edit: On Fri 5 Jun 2015 20:41:45 I started a 35,000 digit computation with the above "automated" program. Edit: The 35,000 digit computation should be done by 10:24:38 am EDT | Sunday, June 14, 2015. In the above "automated" program I forgot to adjust the MaxExtraPrecision, but that shouldn't affect the accuracy in that program. It already computed the real part of the integral to 35,000 digits and the first 30,000 of those are the same as the real part of my previously mentioned 30,000 digit calculation. I will keep you posted. Edit: The 35,000 digit computation finished on Sun 14 Jun 2015 06:52:29. It is attached as 35KMKB.nb. The first 30,000 digits of those are the same as the ones of my previously mentioned 30,000 digit calculation. (That shows the computation didn't take any "wild" turns because of the lack of MaxExtraPrecision.) Further it is a good check of the 30,000 digit run, as all of the bigger computations are of the smaller, because they all are calculated with distinct formulas using different orders of the derivative of x^(1/x). # 727844 seconds Attachments: Answer 2 years ago  Here is an extensive record of records of computing the integral analog of the MRB constant: Here is a graph of those records. (The progression of computed digits is so extreme, it is almost unbelievable!) Answer 2 years ago For 2000 digits Mathematica 10. 2.0 shows some remarkable improvement over 10.1.2 with the above "automated program" for computing the digits of the integral analog of the MRB constant. I will post some speed records that are strictly what the program produces in V 10.2.0, below, no picking and choosing of the methods by a human being. Some results will naturally be slower than my previously mentioned speed records, because I tried so very many methods and recorded only the fastest results. digits seconds 2000 256.3853590 3000 794.4361122 4000 1633.5822870 5000 2858.9390025 10000 17678.7446323 20000 121431.1895170 40000 I got error msg to be continued # I have to change the program for 40,000 digits! I'll post the new program when I get 40,000 to work. Answer 2 years ago As of Wed 29 Jul 2015 11:40:10, one of my computers was happily and busily churning away at 40,000 digits of the integral analog of the MRB constant, using the following formula. Edit: Mathematica crashed at 11:07 PM 8/4/2016 (I used MKB as a symbol for the integral analog because it is called the MKB constant. You can find the name MKB constant at http://www.ebyte.it/library/educards/constants/MathConstants.html .) If you can weed through my code, at the bottom of this reply, you might want to check the formula for the placement of pluses, minuses and imaginary units!!! A little hint when checking if the formula matches the code, d is 80 so Mod[d,4] =0. f[x_] = x^(1/x) : a[n_] := Sum[StirlingS1[n, k]* Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}]; a[0] = 1; g = 2 I/Pi - Sum[-I^(n + 1) a[n]/Pi^(n + 1), {n, 1, 80}] MKB = -g + (I/Pi)^81* Integrate[f[x]*D[f[x], {x, 81}], {x, 1, Infinity}] In traditional form that looks like Here is the code,cleaned up a little: This is the version from Aug 6, 2015 452 pm; for the first time the imaginary part is signed and shown to be multiplied by the imaginary unit! Block[{$MaxExtraPrecision = 200}, prec = 4000; f[x_] = x^(1/x);
ClearAll[a, b, h];
Print[DateString[]];
Print[T0 = SessionTime[]];

If[prec > 35000, d = Ceiling[0.002 prec],
d = Ceiling[0.264086 + 0.00143657 prec]];

h[n_] :=
Sum[StirlingS1[n, k]*
Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}];

h[0] = 1;
g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}];

sinplus1 :=
NIntegrate[
Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity},
WorkingPrecision -> prec*(105/100),
PrecisionGoal -> prec*(105/100)];

cosplus1 :=
NIntegrate[
Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity},
WorkingPrecision -> prec*(105/100),
PrecisionGoal -> prec*(105/100)];

middle := Print[SessionTime[] - T0, " seconds"];

end := Module[{}, Print[SessionTime[] - T0, " seconds"];
Print[c = Abs[a + b]]; Print[DateString[]]];

If[Mod[d, 4] == 0,
Print[N[a = -Re[g] - (1/Pi)^(d + 1)*sinplus1, prec]];
middle;
Print[N[b = -I (Im[g] - (1/Pi)^(d + 1)*cosplus1), prec]];
end];

If[Mod[d, 4] == 1,
Print[N[a = -Re[g] - (1/Pi)^(d + 1)*cosplus1, prec]];
middle;
Print[N[b = -I (Im[g] + (1/Pi)^(d + 1)*sinplus1), prec]]; end];

If[Mod[d, 4] == 2,
Print[N[a = -Re[g] + (1/Pi)^(d + 1)*sinplus1, prec]];
middle;
Print[N[b = -I (Im[g] + (1/Pi)^(d + 1)*cosplus1), prec]];
end];

If[Mod[d, 4] == 3,
Print[N[a = -Re[g] + (1/Pi)^(d + 1)*cosplus1, prec]];
middle;
Print[N[b = -I (Im[g] - (1/Pi)^(d + 1)*sinplus1), prec]];
end];]

Come back to see if I decided whether to try the 40K run again.

EDIT: It looks like I've only got one more test for the program (if it passes) before I retry the 40,000 digit calculation!

EDIT:On Thu 6 Aug 2015 17:23:18, I restarted the 40K run with Windows 10.

EDIT: My first thought was the program took up too much RAM, apparently over 115 GB! ( I have 64GB installed and a 51 GB paging file; nevertheless, Windows 10 closed the Mathematica kernel to keep from the computer from loosing data. Can someone else try the 40K run on their computer? It should take 2 weeks on a fast one. Please let me know if you try it and let me know the results, so I will know I don't have a problem with my computer., If two weeks is too great of a commitment, can you try taking note on the RAM used for two progressively larger runs, like 20K and 30 K? I will do the same, and we can compare notes. Thank You!

EDIT: I've been monitoring memory usage for smaller runs and found the program only uses minimal memory! This makes the action of Windows 10 (closing Mathematica kernel to avoid data loss) all the more a mystery! Could the 40K run really use up all of that RAM?

2 years ago
 I think I've reached an impasse in calculating more digits of the integral analog of the MRB constant. So I believe my best bet in setting a new record is to compute 4,000,000 digits of the MRB constant proper!Here is the data of my previous records on my big computer, the six core Intel(R) Core(TM) i7-3930K CPU @ 3.20 GHz 3.20 GHz with 64 GB of RAM. (Digits are rounded down, slightly.)1.2 million in 120,360 seconds2 million 1,870,579 second3 million 4,15,7400 seconds Remarkably they are very linear; so using the Fit command as shown next, Fit[data, {1, x}, x] I get y, in seconds, = -2.5894022213114775*10^6 + 2.2446041393442626 x.Substituting 4,000,000 for x gives 6,389,014;seconds or almost 2 and a half months. (I would round that to 3 months to be safe!)Before I start, I will listen to your suggestions, or encouragement; if you have any, (I do get lonely when I don't get any feedback.)EDIT:At 10:00 PM Friday Sept 04, 2015 I started the 4,000,000 run of the MRB constant proper. Does anyone want to work on it with me?Here is the code I am using: (*Fastest (at MRB's end) as of 24 dEC 2014.*) Block[{$MaxExtraPrecision = 50}, prec = 4000000;(*Number of required decimals.*)ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 16, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.005 pre]}, chunksize = cores*tsize; n = Floor[1.32 pr]; end = Ceiling[n/chunksize]; Print["Iterations required: ", n]; Print["end ", end]; Print[end*chunksize]; d = Cos[n ArcCos[3]]; {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l; h = Log[ll]/ll; x = N[Exp[h], iprec]; pc = iprec; While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(*N[Exp[Log[ll]/ll], pr]*)x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 4]]; ctab = ParallelTable[Table[c = b - c; ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); c, {l, chunksize}], Method -> "EvaluationsPerKernel" -> 2]; s += ctab.(xvals - 1); start += chunksize; Print["done iter ", k*chunksize, " ", SessionTime[] - T0];, {k, 0, end - 1}]; N[-s/d, pr]]; t2 = Timing[MRBtest2 = expM[prec];]; Print["Difference from 3014991 known digits is ", MRBtest2 - m3M]; MRBtest2] 3014991 digits are saved as m3M.I will post several updates as edits here as I go.EDIT: As of Sat, Sept 5, 2015, the program is indicating that the run will take 3.3 months.EDIT: Having to run the program for over 3 months, there is a chance that I would suffer a power outage longer than the 15 minutes or so that my battery backup will run my computer! Of course I'll keep you posted if that happens. (I've tried calculations that took longer than that before!)EDIT I didn't know that upgrading to Windows 10 caused my Windows update option to read, automatically install updates, so my computer restarted. I started the 4,000,000 digit run again at 8:34 PM, Thursday Sept 10, 2015. This time I made sure the deffer updates option was chosen!!EDIT As of Sept 19,2015, the program says it will take 94 more days to complete the 4 million digit run.EDIT As of Sept 23,2015, 9:20 PM, the program says it will take 90.9563 more days to complete the 4 million digit run.EDIT As of Sept 27,2015, 9:20 PM, the program says it will take 86.9283 more days to complete the 4 million digit run Answer 2 years ago  Windows 10 restarted my machine again.The Windows updates were installed again, arrggghhhh!I don't know if I am going to retry it again real soon. (It really ties up my best computer!)This might give other people a better chance to break my record. Answer 2 years ago  I am waiting for an improvement in Mathematica to try to break any new records, especially because the goal of 4,000,000 digits of the MRB constant requires, in the neighborhood of, 3 months worth of computing time with the best version so far! Also, I thought I would take advantage of this opportunity to introduce this discussion to any new members, that they too may learn or be encouraged by its voluminous replies, as manifold experienced members have -- and a few have added to it-- contributing to its success. Thank you all!In this post we have listed record calculations of the MRB constant (MRB) and its integral analog (MKB); found new formulas for computing each (many, many formulas for computing MKB with each more complicated and powerful that the previous like ); looked at records using only special formulas for MRB; and looked at a few related constants, like the infinite power tower of MRB. I haven't decided whether to start a new post when I try my next record, or not; if i decide to start a new one I will post a link to it here. Please feel free to continue to use this post for your contributions or questions! Answer 2 years ago  I've been working on Crandall's code for many digits of the MRB constant. On my big computer Windows 10 is just a little slower than Windows 7 was. This new code runs on Windows 10 at the same speed the my last mention in this blog program did in Windows 7.Here is what it now looks like in traditional form: Here it is in input form: Quiet[prec = 30000; ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 16, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.005*pre]}, chunksize = cores*tsize; n = Floor[1.32*pr]; end = Ceiling[n/chunksize]; Print["Iterations required: ", n]; Print["end ", end]; Print[end*chunksize]; d = (1/2)*(3 + 2*Sqrt[2])^n*(1 + (3 + 2*Sqrt[2])^ (-2*n)); {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[ Table[ll = start + j*tsize + l; x = N[ll^(1/ll), iprec]; pc = iprec; While[pc < pr, pc = Min[3*pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x - (2*x*y)/(y + ll*(2*ll + y))]; x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 4]]; ctab = ParallelTable[Table[c = b - c; ll = start + l - 2; b *= (2*(ll^2 - n^2))/ (1 + 3*ll + 2*ll^2); c, {l, chunksize}], Method -> "EvaluationsPerKernel" -> 2]; s += ctab . (xvals - 1); start += chunksize; Print["done iter ", k*chunksize, " ", SessionTime[] - T0]; , {k, 0, end - 1}]; N[-s/d, pr]]; t2 = Timing[MRBtest2 = expM[prec]; ]; MRBtest2] Also, here is my new fastest program for calculating digits of the MRB constant via the eta formula, This program is nearly twice as fast as the one in I gave "best eta timings" for, above.. prec = 1500; to = SessionTime[]; etaMM[m_, pr_] := Module[{a, d, s, k, b, c}, a[j_] := N[(Log[j + 1]/(j + 1))^m, pr]; n = Floor[132 /100 pr]; d = N[(1/2)*(3 + 2*Sqrt[2])^n*(1 + (3 + 2*Sqrt[2])^(-2*n)), prec]; {b, c, s} = {-1, -d, 0}; Do[c = b - c; s = s + c a[k]; b = (k + n) (k - n) b/((k + 1) (k + 1/2)), {k, 0, n - 1}]; Return[N[s/d, pr] (-1)^m]]; eta1 = N[EulerGamma Log[2] - Log[2]^2/2, prec]; MRBtest = eta1 - ParallelSum[(Cos[Pi m]) etaMM[m, prec]/ N[Gamma[m + 1], prec], {m, 2, Floor[.40 prec]}]; Print[MRBtest]; SessionTime[] - to gives (* During evaluation of In[188]:= 0.18785964246206712024851793405427323005590309490013878617200468408947723156466021370329665443310749690384234585625801906123137009475922663043892934889618412083733662608161360273812637937343528321255276396217148932170207628206217151671540841268044836354167199851976802527598938993914457983505561350964852107120784442309586812949768852694956420425558648367044104252795247106066609263397483410311578167864166891546003422225883800254553968929471142122189105098328712277308020036445215390536395055332203470627551159812828039510219264914673176293516190659816018664245824950697203381992958420935515162514399357600764593291281451709082424915883204169066409334435914806705564692806787007028115009380606938139385953360657987405562062348704329360737819564603104763950664893061360645528067515193508280837376719296866398103094949637496277383049846324563479311575300289212523291816195626973697074865765476071178017195787368300965902260668753656305516567361288150201438756136686552210674305370591039735756191489093690777983203551193362404637253494105428363699717024418551654837279358822008134480961058802030647819619596953756287834812334976385863010140727252923014723333362509185840248037040488819676767601198581116791693527968520441600270861372286889451015102919988536905728659287086875425492533794395347589703563313440382638887986656195980733514739902565778133172261076127975852722742777308985774922305970962572562718836755752978879253616876739403543214513627725492293131262764357321446216187786377154205423128223 Out[192]= 105.7159526*) Answer 2 years ago  On a recent post, I derived a bunch of approximations of the MRB constant.It was MRB constant approximations using TranscendentalRecognize -- http://community.wolfram.com/groups/-/m/t/819837?p_p_auth=m7XRds7k.Some of the last ones consist of arbitrary approximations to any n digits, while still using small coefficients, compared with the thousand of digits of the target. I wonder if anyone could show me how to extract the coefficient like terms which come before the summations, in a list like this:So someone could find a closed form consisting of those coefficient like terms Answer 1 year ago ## Repeated on top of first message When it comes to mine and a few more educated people's passion to calculate many digits and the dislike of a few more educated people; it is all a matter telling us that the human mind is multi faceted in giving passion, to person a, for one task and to person b for another task! Answer 1 year ago  My plan in the previous was to find a numeric list to dot product(.) a set of sums that equal the MRB constant, list or at least gave an improvable result. In the postI have made several numeric list and sums that are inner products, but I haven't found any good pattern to the numeric lists!Below is a sample of what I want to investigating with just a sample value for a, b , c..So from m = NSum[(-1)^n*(n^(1/n) - 1), {n, 1, Infinity},Method -> "AlternatingSigns", WorkingPrecision -> 1000] N[m - (-129858773922357615372945307143544254 - 332618118135196201861563173048187520*3^(1/5) + 360989141074787715535168098417020609*2^(2/5)*3^(1/5))/ 261711912538111957032871243762602971, 200] Out[116]= 0.*10^-99 you could have In[112]:= a = {-129858773922357615372945307143544254, 332618118135196201861563173048187520, 360989141074787715535168098417020609}; In[113]:= b = {1, -3^(1/5), 2^(2/5) 3^(1/5)}; In[114]:= c = 261711912538111957032871243762602971; In[115]:= N[m - (a.b/c), 200] Out[115]= 0.*10^-99 If only there was a formula for the elements of a,b, and c; up to the precision given by TranscendentalRecognize in a polynomial like from,you could derive a,b,c without TranscendentalRecognize getting arbitrary digits of m of the simple formula of a.b/c . Answer 1 year ago  I just reasoned a great reason to compute the MRB constant!As we saw above at https://web.archive.org/web/20130430193005/http://www.perfscipress.com/papers/UniversalTOC25.pdf , Richard Crandall wrote about Unified algorithms for polylogarithm, L-series, and zeta variants. Since this paper was cited 8 times in its first year (how much longer he lived), I believe this paper is suitable for an introduction to most of what is said in the subject mentioned above. In chapter 7 i.e. "Key fundamental constants." I searched for the meaning of Key fundamental constants, on google and found it referred to different constants, all depending on the subject that uses the constant. So a whole section of that chapter was called The MRB constant. So I wonder does considering the MRB constant to be a "Key fundamental constants." make to calculating of its digits any more interesting? Even if you don't want to use all of the computer power to break a record here, If anyone can find a faster algorithm that I could use, that would bee terrific! Answer 1 year ago  I think this article might be starting a second life for,since the efficient analysis of equation 44 is missing:http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.695.5959&rep=rep1&type=pdf <- New occasionally working reference, so subject not dead! This reference, on my site, is up a little more: http://www.marvinrayburns.com/UniversalTOC25.pdf .We touched formula 44 but didn't dig out any of the good stuff that Crandall believed to be in it!Here B is the MRB constant. Plus Crandall acknowledged that there should be a negative sign before first summation. Answer 1 year ago V. 11 is about 1.25 times faster with my newest program for calculating MKB, (the integral analog of the MRB sum). V 10. 4 calculated 20,000 digits in 121431.1895170 seconds and V 11 did it in 96979.6545388 seconds. I've got a little more testing to do, (about 1 day's worth), then I'll try 40,000 digits again, which should take about 12 days. I will post all my updates here, so you might want to save this message as a favorite so you won't loose it. # Update 1 The 40 K automatically started against my wishes on Thu 11 Aug 2016 15:42:08, (due to my pasting two codes at once). I'll keep you informed, how it goes. # Update 2 Windows 10 is pushing an update. Wednesday is the latest it will let me restart. I will restart now with all the updates I can get, Then deffer further ones and hopefully get 12 restart free days to do my 40K. # Update 3 I ran all the updates I could find, differed further ones and restarted 40K on Sun 14 Aug 2016 10:32:40. # Update 4 Widows 10 stopped the calculation! AGAIN! Can anyone else try it and see if you get anywhere? Here is my latest code: (*Other program:For large calculations.Tested for 1000-35000 digits-- \ see post at \ http://community.wolfram.com/groups/-/m/t/366628?p_p_auth=KA7y1gD4 \ and search for "analog" to find pertinent replies.Designed to include \ 40000 digits.A157852 is saved as c,the real part as a and the \ imaginary part as b.*)Block[{$MaxExtraPrecision = 200},
prec = 40000(*Replace 40000 with number of desired digits.40000 \
digits should take two weeks on a 3.5 GH Pentium processor.*);
f[x_] = x^(1/x);
ClearAll[a, b, h];
Print[DateString[]];
Print[T0 = SessionTime[]];
If[prec > 35000, d = Ceiling[0.002 prec],
d = Ceiling[0.264086 + 0.00143657 prec]];
h[n_] :=
Sum[StirlingS1[n, k]*
Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}];
h[0] = 1;
g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}];
sinplus1 :=
NIntegrate[
Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity},
WorkingPrecision -> prec*(105/100),
PrecisionGoal -> prec*(105/100)];
cosplus1 :=
NIntegrate[
Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity},
WorkingPrecision -> prec*(105/100),
PrecisionGoal -> prec*(105/100)];
middle := Print[SessionTime[] - T0, " seconds"];
end := Module[{}, Print[SessionTime[] - T0, " seconds"];
Print[c = Abs[a + b]]; Print[DateString[]]];
If[Mod[d, 4] == 0,
Print[N[a = -Re[g] - (1/Pi)^(d + 1)*sinplus1, prec]];
middle;
Print[N[b = -I (Im[g] - (1/Pi)^(d + 1)*cosplus1), prec]];
end];
If[Mod[d, 4] == 1,
Print[N[a = -Re[g] - (1/Pi)^(d + 1)*cosplus1, prec]];
middle;
Print[N[b = -I (Im[g] + (1/Pi)^(d + 1)*sinplus1), prec]]; end];
If[Mod[d, 4] == 2,
Print[N[a = -Re[g] + (1/Pi)^(d + 1)*sinplus1, prec]];
middle;
Print[N[b = -I (Im[g] + (1/Pi)^(d + 1)*cosplus1), prec]];
end];
If[Mod[d, 4] == 3,
Print[N[a = -Re[g] + (1/Pi)^(d + 1)*cosplus1, prec]];
middle;
Print[N[b = -I (Im[g] - (1/Pi)^(d + 1)*sinplus1), prec]];
end];] (*Marvin Ray Burns,Aug 06 2015*)
1 year ago
 V 11 has a significant improvement in my new most recently mentioned fastest program for calculating digits of the MRB constant via the eta formula, Here are some timings: Digits seconds 1500 42.6386632 2000 127.3101969 3000 530.4442911 4000 1860.1966540 5000 3875.6978162 6000 8596.9347275 10,000 53667.6315476 From an previous message that starts with "How about computing the MRB constant from Crandall's eta derivative formulas?" here are my first two sets of records to compare with the just mentioned ones. You can see that I increased time efficiency by 10 to 29 to even 72 fold for select computations! In the tests used in that "previous message," 4000 or more digit computations produced a seemingly indefinitely long hang-on. Digits Seconds 500 36.831836 1000 717.308198 1500 2989.759165 2000 3752.354453 Digits Seconds 500 9.874863 1000 62.587601 1500 219.41540 2000 1008.842867 2500 2659.208646 3000 5552.902395 3500 10233.821601 Comparing first of the just mentioned 2000 digit computations with the "significant improvement" one we get the following. 3752/127 ~=29. And from the slowest to the fastest 1500 digit run we get2989/42 ~=72,
1 year ago
 To try to get windows 10 from closing Mathematica during the computation I tried the instructions found at https://www.autoitscript.com/forum/topic/177749-stopping-windows-10-from-auto-closing-programs-to-free-up-ram/ . I will record progress in this spot as I did before.UPDATE I followed the memory usage on my computer and it did use around 64 GB of RAM. And then Windows closed down the Mathematica kernel. I assume that If I can ever afford to maximize my RAM to its 128GB limit the computation will be successful!Anyone have better luck?
1 year ago
 I have two versions of MRB constant code that I presently use, a short one (non-parallel) and a long one. Both are given below.Short (*Newer loop with Newton interior.*)prec = 5000;(*Number of required \ decimals.*)expM[pr_] := Module[{a, d, s, k, b, c}, n = Floor[1.32 pr]; Print["Iterations required: ", n]; d = N[(3 + Sqrt[8])^n, pr + 10]; d = Round[1/2 (d + 1/d)]; {b, c, s} = {-1, -d, 0}; T0 = SessionTime[]; Do[c = b - c; x = N[E^(Log[k + 1]/(k + 1)), iprec = Ceiling[prec/128]]; pc = iprec; Do[nprec = Min[2 pc, pr]; x = SetPrecision[x, nprec];(*Adjust precision*) x = N[x - x/(k + 1) + 1/x^k, nprec]; pc *= 2; If[nprec >= pr, Break[]], {ct, 1, 19}]; s += c*(x - 1); b *= 2 (k + n) (k - n)/((k + 1) (2 k + 1)); If[Mod[k, 1000] == 0, Print["Iterations: ", k, " Cumulative time (sec): ", SessionTime[] - T0];], {k, 0, n - 1}]; N[-s/d, pr]]; MRBtest2 = expM[prec] Long (**Fastest (at MRB's end) as of 25 Jul 2014*.*)DateString[] prec = 7000; (**Number of required decimals.*.*)ClearSystemCache[]; T0 = SessionTime[]; expM[pre_] := Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 12, tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.005 pre]}, chunksize = cores*tsize; n = Floor[1.32 pr]; end = Ceiling[n/chunksize]; Print["Iterations required: ", n]; Print["end ", end]; Print[end*chunksize]; d = ChebyshevT[n, 3]; {b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0}; iprec = Ceiling[pr/27]; Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l; x = N[E^(Log[ll]/(ll)), iprec]; pc = iprec; While[pc < pr, pc = Min[3 pc, pr]; x = SetPrecision[x, pc]; y = x^ll - ll; x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(**N[Exp[Log[ll]/ ll],pr]**)x, {l, 0, tsize - 1}], {j, 0, cores - 1}, Method -> "EvaluationsPerKernel" -> 4]]; ctab = ParallelTable[Table[c = b - c; ll = start + l - 2; b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1)); c, {l, chunksize}], Method -> "EvaluationsPerKernel" -> 2]; s += ctab.(xvals - 1); start += chunksize; Print["done iter ", k*chunksize, " ", SessionTime[] - T0];, {k, 0, end - 1}]; N[-s/d, pr]]; t2 = Timing[MRBtest2 = expM[prec];]; DateString[] Print[MRBtest2] MRBtest2 - MRBtest3
1 year ago
 Try to beat these MRB constant records!
8 months ago

As mentioned a few times above Richard Crandall made the following observations about the MRB constant:

His last comment is that we might find some algorithm advantages by focusing on the fixed argument of 0 for the eta derivatives. He gives Chapter 3 as a possible source to find the algorithm advantages. (Here is another link to his paper: http://marvinrayburns.com/UniversalTOC25.pdf .)

Here is what he was looking for.

c[0] = 0;
c[1] = Log[2]; Table[
c[n] = -2*Derivative[n - 1][Zeta][0], {n, 2, 22}];
l = CoefficientList[Normal[Series[Exp[(-x)*E^x], {x, 0, 21}]],
x]; N[-Sum[
l[[n + 1]]*((-1)^
n*(Sum[(-1)^x*(c[x + 1]*Binomial[n, x]*c[1]^(n - x)), {x, 1,
n - 1}] + c[1]^n +
c[0]) + 1/2 c[n + 1]), {n, 1, 17}]]

Mathematica still has trouble with higher zeta derivatives of 0, however. I'm going to look for a procedure for computing them better.

## 

I computed the zeta derivatives of 0 by

zeta0[p_] := (-1)^p Sum[StieltjesGamma[p + k]/k!, {k, 0, Infinity}] -
p!

and got the same results Mathematica got for the higher zeta derivatives of 0.

I also checked the CoefficientList of

CoefficientList[Normal[Series[Exp[(-x)*E^x], {x, 0, 21}]]

by the following, and haven't found where Mathematica is having trouble!

c[0] = 0; c[1] = Log[2]; Table[c[n] = -2*Derivative[n - 1][Zeta][0],
{n, 2,
18}]; N[-Sum[(Sum[(-1)^d*d^(n - d)*Binomial[n, d], {d, 1, n}]*
((-1)^
n*(Sum[(-1)^x*(c[x + 1]*Binomial[n, x]*c[1]^(n - x)),
{x, 1, n - 1}] + c[1]^n + c[0]) + (1/2)*
c[n + 1]))/n!, {n, 1, 17}]]

I even checked the binomials by their formula. Also the StieltjesGamma and Derivative[n ][Zeta][0]] results appear smooth for discrete operations:

Table[N[Derivative[n ][Zeta][0]], {n, 10, 25}]

and

N[Table[StieltjesGamma[x], {x, 1, 30}]]

.

## 

I found the problem as to why Mathematica doesn't increase its accuracy in the above formula after 17 iterations: Mathematica will give a numeric solution to many zeta derivatives of 0, only to machine precision!

For machine precision:

Table[N[Derivative[n ][Zeta][0]], {n, 1, 30}]

gives the machine sized numeric to many derivatives:

{-0.918939, -2.00636, -6.00471, -23.9971, -120., -720.001, \
-5040., -40320., -362880., -3.6288*10^6, -3.99168*10^7, \
-4.79002*10^8, -6.22702*10^9, -8.71783*10^10, -1.30767*10^12, \
-2.09228*10^13, -3.55687*10^14, -6.40237*10^15, -1.21645*10^17, \
-2.4329*10^18, -5.10909*10^19, -1.124*10^21, -2.5852*10^22, \
-6.20448*10^23, -1.55112*10^25, -4.03291*10^26, -1.08889*10^28, \
-3.04888*10^29, -8.84176*10^30, -2.65253*10^32}.

Then

Table[N[Derivative[n][Zeta][0], 10], {n, 1, 30}]

gives the high precision numeric for only the first 3 derivatives:

{-0.9189385332046727417803297364056176398610., \
-2.006356455908584851210100026729960438210.,
-6.0047111668622544478`10., Derivative[4][Zeta][0],
Derivative[5][Zeta][0], Derivative[6][Zeta][0],
Derivative[7][Zeta][0], Derivative[8][Zeta][0],
Derivative[9][Zeta][0], Derivative[10][Zeta][0],
Derivative[11][Zeta][0], Derivative[12][Zeta][0],
Derivative[13][Zeta][0], Derivative[14][Zeta][0],
Derivative[15][Zeta][0], Derivative[16][Zeta][0],
Derivative[17][Zeta][0], Derivative[18][Zeta][0],
Derivative[19][Zeta][0], Derivative[20][Zeta][0],
Derivative[21][Zeta][0], Derivative[22][Zeta][0],
Derivative[23][Zeta][0], Derivative[24][Zeta][0],
Derivative[25][Zeta][0], Derivative[26][Zeta][0],
Derivative[27][Zeta][0], Derivative[28][Zeta][0],
Derivative[29][Zeta][0], Derivative[30][Zeta][0]}.

With some confidence I can say

8 months ago
 I fixed it!!!I forced Mathematica to give up the needed digits of the zeta derivatives.My new program does not really compete with the previously mentioned record-breakers. But it is a fairly efficient use of Crandall's formula (44). First it calculates many digits on MRB the old fashion way to use as a check. So here is my new program: m = NSum[(-1)^n (n^(1/n) - 1), {n, 1, Infinity}, Method -> "AlternatingSigns", WorkingPrecision -> 120]; Block[{$MaxExtraPrecision = 1000}, s = 100(*number of terms*); c[1] = Log[2]; c[n_] := N[1 - 2*Derivative[n - 1][Zeta][0], Floor[3/2 s]] - 1; mtest = (-Sum[(Sum[(-1)^d*d^(n - d)*Binomial[n, d], {d, 1, n}]*((-1)^ n*(Sum[(-1)^x*(c[x + 1]*Binomial[n, x]*c[1]^(n - x)), {x, 1, n - 1}] + c[1]^n) + (1/2)*c[n + 1]))/n!, {n, 1, s}])]; Print[s, " terms give ", mtest, ", and the error is ", m - mtest] 100 terms give 0.18785964486341229249700337222, and the error is -2.40134517224848543816*10^-9Changing s to 200 shows all correct digits:200 terms give 0.187859642465, and the error is -3.*10^-12 Answer 8 months ago  I noticed that accuracy is not directly related to the number of terms used in the programs of Crandall's formula (44).I made a slight improvement in my program that computes Crandall's eta derivative of 0 formula (44):First compute a known value of the MRB constant: m = NSum[(-1)^n (n^(1/n) - 1), {n, 1, Infinity}, Method -> "AlternatingSigns", WorkingPrecision -> 120]; followed by Block[{$MaxExtraPrecision = 1000}, s = 129(*number of terms*); sc[t_] := StieltjesGamma[ t];(*-Pi/(2 n+2) NIntegrate[Log[1/2+I x]^(n+1)/Cosh[Pi \ x]^2,{x,-Infinity,Infinity},WorkingPrecision\[Rule]50];*) zeta0[p_] := N[(-1)^p Sum[sc[p + k]/k!, {k, 0, Infinity}] - p!, 3/2 s]; c[1] = Log[2]; c[n_] := N[1 - 2*zeta0[n - 1], Floor[3/2 s]] - 1; mtest = (-Sum[(Sum[(-1)^d*d^(n - d)*Binomial[n, d], {d, 1, n}]*((-1)^ n*(Sum[(-1)^x*(c[x + 1]*Binomial[n, x]*c[1]^(n - x)), {x, 1, n - 1}] + c[1]^n) + (1/2)*c[n + 1]))/n!, {n, 1, s}])]; Print[s, " terms give ", mtest, ", and the error is ", m - mtest] producing129 terms give 0.187859642466672738539355017, and the error is -4.605618290837083*10^-12.
7 months ago

# Reply starts with UPDATE on 8/13/2017 below.

Search 8/13/2017 below for where this reply was placed.

6 months ago

# UPDATE on 8/13/2017

Search 8/13/2017 in the above for this reply's place in Que. I placed it here to make it easy to find on this late of a date.

Here is another great improvement in calculating digits of MRB though this first Crandall eta formula

Here is an example calculating 1500 digits in less than 23 seconds. First calculate known accurate digits:

mTrue = NSum[(-1)^n (n^(1/n) - 1), {n, 1, Infinity},
Method -> "AlternatingSigns", WorkingPrecision -> 2000];

Then use this program:

prec = 1500;
to = SessionTime[];
etaMM[m_, pr_] :=
Module[{a, s, k, b, c}, a[j_] := (Log[j + 1]/(j + 1))^m;
{b, c, s} = {-1, -d, 0};
Do[c = b - c;
s = s + c a[k];
b = (k + n) (k - n) b/((k + 1) (k + 1/2)), {k, 0, n - 1}];
Return[N[s/d, pr] (-1)^m]];
eta1 = N[EulerGamma Log[2] - Log[2]^2/2, prec]; n =
Floor[132/100 prec]; d = N[ChebyshevT[n, 3], prec];
MRBtest =
eta1 - ParallelSum[(Cos[Pi m]) etaMM[m, prec]/
N[Gamma[m + 1], prec], {m, 2, Floor[.365 prec]}];
Print[mTrue - MRBtest]; SessionTime[] - to

The larger prec (precission) you want to choose, the smaller number you can replace .365 with in
ParallelSum[(Cos[Pi m]) etaMM[m, prec]/N[Gamma[m + 1], prec], {m, 2, Floor[.365 prec]}. For prec=5000, change .365 to .31 (approximately). This advantage for smaller upper limits being needed in the sum is due to the Gamma (!) in the denominator. On my big computer 5000 digits now gives a timing of around 1800 seconds (1/2 an hour). See message above the previous occurrence of 8/13/2017 for previous records and records before those.

6 months ago

## Update Feb 12, 2018

In celebration of 72000 views, this reply moved to what is presently the bottom of the Que. It contains many updates and will gain even more.

4 months ago

## Concentrating on integral analog of the MRB constant:

Search "integral analog" in the above messages for understanding of the integral anaolg of the MRB constant. And search "For 2000 digits Mathematica 10. 2.0" for my history of calculating 40,000 digits of it.

. The basic program I wrote to calculate the Integral analog of the MRB constant is

(*Other program:For large calculations.Tested for 1000-35000 digits-- \
see post at \
http://community.wolfram.com/groups/-/m/t/366628?p_p_auth=KA7y1gD4 \
and search for "analog" to find pertinent replies.Designed to include \
40000 digits.A157852 is saved as c,the real part as a and the \
imaginary part as b.*)Block[{\$MaxExtraPrecision = 200},
prec = 40000(*Replace 40000 with number of desired digits.40000 \
digits should take two weeks on a 3.5 GH Pentium processor.*);
f[x_] = x^(1/x);
ClearAll[a, b, h];
Print[DateString[]];
Print[T0 = SessionTime[]];
If[prec > 35000, d = Ceiling[0.264086 + 0.0017 prec],
d = Ceiling[0.264086 + 0.00143657 prec]];
h[n_] :=
Sum[StirlingS1[n, k]*
Sum[(-j)^(k - j)*Binomial[k, j], {j, 0, k}], {k, 1, n}];
h[0] = 1;
g = 2 I/Pi - Sum[-I^(n + 1) h[n]/Pi^(n + 1), {n, 1, d}];
sinplus1 :=
NIntegrate[
Simplify[Sin[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity},
WorkingPrecision -> prec*(105/100),
PrecisionGoal -> prec*(105/100)];
cosplus1 :=
NIntegrate[
Simplify[Cos[Pi*x]*D[f[x], {x, d + 1}]], {x, 1, Infinity},
WorkingPrecision -> prec*(105/100),
PrecisionGoal -> prec*(105/100)];
middle := Print[SessionTime[] - T0, " seconds"];
end := Module[{}, Print[SessionTime[] - T0, " seconds"];
Print[c = Abs[a + b]]; Print[DateString[]]];
If[Mod[d, 4] == 0,
Print[N[a = -Re[g] - (1/Pi)^(d + 1)*sinplus1, prec]];
middle;
Print[N[b = -I (Im[g] - (1/Pi)^(d + 1)*cosplus1), prec]];
end];
If[Mod[d, 4] == 1,
Print[N[a = -Re[g] - (1/Pi)^(d + 1)*cosplus1, prec]];
middle;
Print[N[b = -I (Im[g] + (1/Pi)^(d + 1)*sinplus1), prec]]; end];
If[Mod[d, 4] == 2,
Print[N[a = -Re[g] + (1/Pi)^(d + 1)*sinplus1, prec]];
middle;
Print[N[b = -I (Im[g] + (1/Pi)^(d + 1)*cosplus1), prec]];
end];
If[Mod[d, 4] == 3,
Print[N[a = -Re[g] + (1/Pi)^(d + 1)*cosplus1, prec]];
middle;
Print[N[b = -I (Im[g] - (1/Pi)^(d + 1)*sinplus1), prec]];
end];] (*Marvin Ray Burns,

I think I found out why the integral analog of the MRB constant is so hard to calculate to prec=40000 digits! I've been using too high of an order of the derivative of x^(1/x). I've been running out of memory because of using the 80th derivative from d = Ceiling[0.002 prec], because the 58th derivative from Ceiling[0.264086 + 0.00143657 prec] was apparently to small leaving an error statement. I just now asked myself, why make such a big jump? When my big computer gets back from its tuneup I think I will try Ceiling[0.00146 prec] = 59th derivative.

## EDIT

I tried Ceiling[0.00146 prec] and Ceiling[0.00145 prec] in Mathematica 11.0 and lost the kernel both times after 6 - 12 hours!

I'm now trying Ceiling[0.0017 prec] with v 10.4. It's been over 12 hours and I've not lost the kernel yet. Wish me luck!

## EDIT

I got the following error message and a real part that does not agree with previous computations.

"NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in x near {x} = {<<42008>>}. NIntegrate obtained -<<42012>> and <<42014>> for the integral and error estimates."

I'm now trying Ceiling[0.0018 prec] with v 10.4....

Ceiling[0.0018 prec] with v 10.4 gave the same error.

I'm working on a new program that uses less memory; stay tuned!

2 months ago

I just calculated 100,500 digits of the MRB constant in 2,359 seconds and a processor time (from the Timing[] command) of 1,196 seconds in V11.2. Here is a copy of an email from the late Richard Crandall showing his best for 100,000 digits using V9:

Re: 100,000 digits Richard Crandall crandall@reed.edu 11/30/12 at 11:00 PM To Marvin Burns Message body Something is really awry... On my humble Mac Book (2.2 GHz, Intel Core i7, 4 GB RAM) the attached program---my best so far takes 3600 seconds. That's right, 100K digitsw in one hour.

What timing do you get for the attached? (If you again get like 8000 seconds, tell me what kind of RAM etc. you have...)

-r

(* Fastest (at RC's end) as of 30 Nov 2012. *)

prec = 100000; (* Number of required decimals. *)
ClearSystemCache[];
T0 = SessionTime[];
expM[pre_] :=
Module[{a, d, s, k, bb, c, n, end, iprec, xvals, x, pc, cores = 4,
tsize = 2^7, chunksize, start = 1, ll, ctab, pr = Floor[1.02 pre]},
chunksize = cores*tsize;
n = Floor[1.32 pr];
end = Ceiling[n/chunksize];
Print["Iterations required: ", n];
Print["end ", end];
Print[end*chunksize];
d = N[(3 + Sqrt[8])^n, pr + 10];
d = Round[1/2 (d + 1/d)];
{b, c, s} = {SetPrecision[-1, 1.1*n], -d, 0};
iprec = Ceiling[pr/27];
Do[xvals = Flatten[ParallelTable[Table[ll = start + j*tsize + l;
x = N[E^(Log[ll]/(ll)), iprec];
pc = iprec;
While[pc < pr, pc = Min[3 pc, pr];
x = SetPrecision[x, pc];
y = x^ll - ll;
x = x (1 - 2 y/((ll + 1) y + 2 ll ll));];(*N[Exp[Log[ll]/ll],pr]*)
x, {l, 0, tsize - 1}], {j, 0, cores - 1},
Method -> "EvaluationsPerKernel" -> 1]];
ctab = Table[c = b - c;
ll = start + l - 2;
b *= 2 (ll + n) (ll - n)/((ll + 1) (2 ll + 1));
c, {l, chunksize}];
s += ctab.(xvals - 1);
start += chunksize;
Print["done iter ", k*chunksize, " ", SessionTime[] - T0];, {k, 0,
end - 1}];
N[-s/d, pr]];

t2 = Timing[MRBtest2 = expM[prec];];
MRBtest2 - MRBtest3

(*===============*)

On Nov 30, 2012, at 5:29 PM, Marvin Burns wrote:

On my 6 core AMD: 8511 seconds Error: 1.83506*10^-101989

----- Forwarded Message ---- From: Marvin Burns marvin@marvinrayburns.com To: Richard Crandall crandall@reed.edu Sent: Fri, November 30, 2012 6:09:09 PM Subject: Re: 100,000 digits

It took 8777 seconds. I have my laptop set not to go to sleep.

Error:1.83506*10^-101989.

I will try it on my 6 core AMD. From: Richard Crandall crandall@reed.edu To: Marvin Burns marvin@marvinrayburns.com Sent: Fri, November 30, 2012 2:49:26 AM Subject: Re: 100,000 digits

You report 8500 seconds on your 100K-digit run. Can you please try this version, to see how long it takes? (The output format is different, but it should have the same 100K-digit accuracy...

-r

PS. This has parallism but of a different style that might work well on your machine...

<MRB100Ktest[6 core AMD].nb>

(Above Is where I told Mr Crandall that my AMD was a 6 core, but if I remember right it was only a quad core.)

I would like to beat Crandall's million+ digit run. But I won't be able to right away. I am computing a fresh million+ digits using V 11.2, but it will take about 5 days vs. Crandall's 3.18 days. At least that's better than my first run which took 7.9 days and was started on ‎Tuesday, ‎January ‎22, ‎2013, ‏‎6:53:42 PM. It is attached as "1M." I used this same six core Intel(R) Core(TM) i7-3930K CPU @ 3.20 GHz 3.20 GHz and V9.0.)

## EDIT

1004993 digits finished in 420,123. sec or 4.86254 days! Processor time was 188,939 sec or 2.18679 days. See attached 1M2017

Attachments:
1 month ago

Search "I think I've reached an impasse" in the above messages for my history of calculating 4 million digits of the MRB constant.

Well I started the 4 million digit search again. This time I installed the Windows 10 Fall Creators Update to try to avoid the automatic restarts! Per my habit above I will post updates in this message. The code I'm using is the same mentioned twice above that starts with

## UPDATE

The Windows 10 Fall Creators Update restarted my machine on Thursday, Nov 16, 2017 at 5:33 pm!

## Dec 19, 2017

I applied some new update policies via LogMeIn, for windows to download but not apply updates. That should work!

The new run looks like

1 : 15 : 58 pm EST | Sunday, December 17, 2017
...  done iter 96256 165140.1164436.

Then

(165140.1164436) sec/(96256 iterations) (5306880 iterations)/(3600sec/hour)/(24 hour/day)

gives 105.378 day.

Entering "December 17, 2017+106 days" into Wolfram Alpha tells us the computation should be done around Monday, April 2, 2018.

## Jan 02, 2018

My computer went to sleep sometime between 10:15:51 pm EST | Thursday, December 28, 2017 and 9:22:00 am EST | Sunday, December 31, 2017 . It seems that my system needs someone to visit every other day. Otherwise it falls into a sleep that requires a restart. Does that make any sense? I have the settings for it to never sleep.

## Jan 07, 2018

I always remote into my big computer with my laptop and the big computer was registering "offline." I thought it had gone to sleep because I couldn't remote on to it. It wasn't asleep!

I hooked up a working monitor to it and found it only needed to be logged on to. So, that is another lesson for me!

With that learned I'm trying the 4 million digit again. I'll post some stats soon.

## Jan 09, 2018

The output starts out with

"Sun 7 Jan 2018 15:14:50"

Iterations required: 5306398

end 3455

5306880

done iter 0     2374.6810043.

The latest result is done iter 92160 153642.6180409.

Entering

(153642.6180409) sec/(92160 iterations) (5306880iterations)/(3600 sec/hour)/(24 hour/day)

gives

102.399 day or 3.4 months for the computation to run.

Entering "Sun 7 Jan 2018 15:14:50 EST+103 days" into Wolfram Alpha tells us the computation should be done around 3:14:50 pm EST | Friday, April 20, 2018. That is 101 days from now.

## Jan 17, 2018, 10:00 pm

Here is the latest output: done iter 531456 888255.5435780.

N[531456/5306880*100] says the computation is 1/10th of the way done!!!

## Jan 20, 2018, 5:39 am

Here is the latest output: done iter 654336 1.0955238835298*10^6.

Then

(1.0955238835298*10^6 sec) /(654336 iterations) (5306880 iterations)/(3600 sec/hour)/(24 hour/day)

gives 102.836 day.

Entering "Sun 7 Jan 2018 15:14:50 EST+102.836 days" into Wolfram Alpha tells us the computation should be done around "11:18:40 pm EST | Friday, April 20, 2018". That is 90.14 days from now.

## Feb 04, 2018,8:04 am

Here is the latest output: done iter done iter done iter 1420800 2.3920453426489*10^6.

Then

(2.3920453426489*10^6 sec) /(1420800 iterations) (5306880iterations)/(3600 sec/hour)/(24 hour/day)

gives 103.41 day.

Entering "7 Jan 2018 15:14:50 EST+103.41 days" into Wolfram Alpha tells us the computation should be done around "1:05:14 am EST | Saturday, April 21, 2018". That is 75.71 days from now.

## Feb 11, 2018, 8:13 am

Here is the latest output: done iter 1778688 2.9979034719762*10^6.

N[1778688/5306880*100] says the computation is more than 1/3 done!!!

## Feb 16, 2018,8:20 am

Here is the latest output: done iter 2032128 3.4292035877594*10^6.

Then

(3.4292035877594*10^6 sec) /(2032128 iterations) (5306880 iterations)/(3600 sec/hour)/(24 hour/day)

gives 103.65 day

Entering "7 Jan 2018 15:14:50 EST+103.65 days" into Wolfram Alpha tells us the computation should be done around "6:50:50 am EST | Saturday, April 21, 2018". That is 63.94 days from now.

## Feb 20, 2018, 5:06 am

Here is the latest output: done iter 2228736 3.7632326187330*10^6.

N[2228736/5306880*100] says the computation is 42% done.

## Feb 24, 2018, 11:10 am

Here is the latest output: done iter 2443776 4.1295732917069*10^6.

Then

(4.1295732917069*10^6 sec) /(2443776 iterations) (5306880 iterations)/(3600 sec/hour)/(24 hour/day) gives 103.793 day.

Entering "7 Jan 2018 15:14:50 EST+103.793 days" into Wolfram Alpha tells us the computation should be done around "10:16:45 am EST | Saturday, April 21, 2018".

That is 55.97 days from now.