5
$\begingroup$

Computing $\ln(10)$ to 6 million digits of precision on my 2.5 GHz machine running Mathematica 12.1 takes about 23 seconds using the methods below. Wish to compute $\ln(x)$ with much higher precision. Is this the fastest I can compute $\ln(10)$ with 6 million digits on my machine?

  1. Use the built-in Log[x] function,
  2. Invert the exponential expression $e^y=x$ which results in a Newton iteration of the form: $$ y_{n+1}=y_n+2\frac{x-e^{y_n}}{x+e^{y_n}} $$
  3. Use the arithmetic-geometric expression: $$ \log(x)\approx \frac{\pi}{2M(1,4/s)}-m\log(2);\quad s=x 2^m>2^{p/2} $$ for $p$ bits of precision.

Unfortunately, these all take about the same amount of time. The code below is for 6 million digits and the best time is about 23 seconds. I would like to compute $\ln(x)$ to about 60 million digits of precision as fast as possible in regards to computing $\Phi(0,10^8)$ explained in the links to this problem: https://artofproblemsolving.com/community/c7h2320487_fixed_points_of_an_iterated_exponential

totalD = 6000000;
(*
 set up A-G mean parameters
*)
myPi = SetPrecision[Pi, totalD];
myLog2 = SetPrecision[Log[2], totalD];
pFun[x_] := Ceiling[x Log[10]/Log[2]];
mFun[p_] := Ceiling[1/Log[2] (p/2 Log[2] - Log[10])];
sFun[m_] := 10 2^m;
(*
  check built-in Log:
*)
AbsoluteTiming[
 actVal = SetPrecision[Log[10], totalD];
 ]
(*
 check arithmetic-geometric mean approach
*)
AbsoluteTiming[
 pVal = pFun[totalD];
 mVal = mFun[pVal];
 sVal = sFun[mVal];
 denom = SetPrecision[ArithmeticGeometricMean[1, 4/sVal], totalD];
 myVal = SetPrecision[myPi/(2 denom) - mVal myLog2, totalD];
 ]
(*
  check just the quotient expression of the inverted exp expression
*)
y0 = 23/10;
AbsoluteTiming[
 SetPrecision[(
   10 - Exp[y0])/(10 + Exp[y0]), totalD];
 ]

{23.911, Null}

{22.7023, Null}

{50.6573, Null}
$\endgroup$
4
  • $\begingroup$ What for? Isn't it art for art's sake? $\endgroup$ – user64494 Oct 30 '20 at 12:46
  • $\begingroup$ This: artofproblemsolving.com/community/… $\Phi(0,10^8)$ becomes problematic. $\endgroup$ – Dominic Oct 30 '20 at 12:52
  • 1
    $\begingroup$ Related: mathematica.stackexchange.com/questions/124418/… $\endgroup$ – Michael E2 Oct 30 '20 at 14:05
  • $\begingroup$ Thanks for that Michael. Seems from that thread Log[x] is already optimized and even though the arithmetic-geometric method is faster above, I don't have a good handle of how accurate the precision is that way. Only option then is just to either let it run for days or use a faster machine. $\endgroup$ – Dominic Oct 30 '20 at 14:35
3
$\begingroup$

Well It took a few minutes but I was able to run N[Log[10],6*10^6] on my laptop (Intel Core i7-8750H @ 2.2 GHz) in 330 seconds. Looking into the problem I found the OEIS sequence for the digits of Log(10) https://oeis.org/A002392 and in the references of this page the current world record of 1,200,000,000,100 digits: https://ehfd.github.io/world-record/natural-logarithm-of-10-log10/. This page includes a link to download the digits (more then 400 GB of data) and details on how they were computed. I have not found a good solution (in the rather limited time I spend researching this) to compute the digits faster in Mathematica BUT:

Reading about this I came across the tool y-cruncher (can be downloaded from multiple sources) which is purpose build to compute digits of mathematical constants. Using this command line tool to compute digits of Log[10] is very easy:

.\y-cruncher.exe custom log:10 -dec:100m -o "c:/log10/"

(on Windows) computes 100 million decimal digits of Log[10] and saves them to "c:/log10/". This took only 28 seconds on my machine (more then 10 times faster then Mathematica where I only computed 60 million digits). I uploaded the decimal digits here (43 MB compressed).

$\endgroup$
4
  • $\begingroup$ Thanks for that Nova. Wouldn't be able to use OEIS or y-cruncher as I need all the data to do at least two Newton Iterations. Took 10 minutes for 60 million digits on my 2.5 GHz machine and I need at least 60 similar 60-million digit precision operations to confirm $\Phi(0,10^8)$ so about 10 hours, I suspect I would run out of memory (have 16 GBytes) however. $\endgroup$ – Dominic Oct 31 '20 at 10:17
  • $\begingroup$ Well OEIS has only 80 digits online but I do not really see the problem regarding runtime and memory: 10 minutes times 60 is ten hours (one night) and those digits should not take up 16 GB: 100 million digits are around 100MB in a plain text file which would lead me to assume that 60 would take up 6GB. $\endgroup$ – N0va Oct 31 '20 at 15:09
  • $\begingroup$ Ok, will try to set it up tonight. Actually have 10 GB available and 18 virtual. $\endgroup$ – Dominic Oct 31 '20 at 15:51
  • $\begingroup$ Went over the precision requirements and decided I only needed 32 million digits of working precision. Computed and verified $\Phi(0,10^8)$ in about 2 hours. $\Phi(0,10^9)$ would require 360 million. Suspect that's not reasonably attainable on my machine so will leave it for someone else to solve. $\endgroup$ – Dominic Oct 31 '20 at 19:07
0
$\begingroup$

UPDATE: As per comments below, this approach will not work.

cn = Compile[{{x, _Real}}, Log[x]];
AbsoluteTiming[res = SetPrecision[cn[10], 6000000];]   (* {0.003, Null} *}
$\endgroup$
2
  • 6
    $\begingroup$ Afraid that won't work. Compile uses machine precision. Try comparing the results of that to just a N[Log[10],100]. $\endgroup$ – Dominic Oct 30 '20 at 13:20
  • $\begingroup$ Good point! I was fooled by Precision[res] giving the required precision. Passed machine precison we just have a bunch of zeros. $\endgroup$ – Jean-Pierre Oct 30 '20 at 13:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.