Message Boards Message Boards

1
|
7816 Views
|
8 Replies
|
9 Total Likes
View groups...
Share
Share this post:

Record breaking direct summation of MRB constant

The MRB constant is approximated by

m = NSum[(-1)^k (k^(1/k) - 1), {k, 1, Infinity}, WorkingPrecision -> 200, 
  Method -> "AlternatingSigns"]

.

This method is pitifully inefficient in its direct summation of MRB constant ie. -1^(1/1)+2^(1/2)-3^(1/3)+... . Here is how much error there is after 10,000 terms:

 m - (NSum[(-1)^k (k^(1/k) - 1), {k, 1, 10^4}, 
   WorkingPrecision -> 200, Method -> "AlternatingSigns"])

(* -0.
0004607086148833843903798558880939839952579916174598636084830590133354\
3604452568023408196838217734393902621089862492299952504355243485749500\
9513789757123430499002402147185148514618052628511533818725*)

.

Due to the work of Daniel Lichtblau, me, and the late Richard Crandall, found in the middle of my most successful post,

http://community.wolfram.com/groups/-/m/t/366628?ppauth=zGnyw9gs

, in the second reply that starts with the phrase "Daniel Lichtblau and others," without the quotes; I think I found the optimum direct summation method of the m= -1^(1/1)+f1 +2^(1/2)+f2 -3^(1/3)+f3+ f(k)(x)+ ...- sum[f[k]]. That is

f[k_]=Log[(k - Log[k])/k]

. Then subtract Sum[(-1)^k Log[(k - Log[k])/k] . Here is an example showing the much smaller error in adding the first 10^4 terms:

     f[k_] := Log[(k - Log[k])/
        k]; m - (NSum[(-1)^k (k^(1/k) - 1 + f[k]), {k, 1, 10^4}, 
         WorkingPrecision -> 200, Method -> "AlternatingSigns"] - 
        NSum[(-1)^k *f[k], {k, 1, Infinity}, WorkingPrecision -> 200, 
         Method -> "AlternatingSigns"])

(* 6.
5176019395153886286218356168146776861946664225985899368189404704176944\
5768612992269488431506197821023123569647840336109019435539071457039552\
26445450414770330729207675435060098684928265988*10^-11*)

Why don't see if you can find a better algorithm!

POSTED BY: Marvin Ray Burns
8 Replies
POSTED BY: Henrik Schachner

Dear Henrik,

that's a beautiful piece of mathemagic.

Thanks for sharing,

Marco

POSTED BY: Marco Thiel

Dear Marco,

thank you for your nice and encouraging comments! Yes, it is almost magical - and it shows a sort of "beauty of simplicity"! I learned about Shanks transformation when reading - or trying to read - Benders book (which for me is somewhat challenging) and watching his lectures (which I find excellent!). Bender shows lots of most interesting methods I had never heard of before, particularly not in my math courses (about 30 years ago ...).

Best wishes -- Henrik

POSTED BY: Henrik Schachner

Henrik et al,

I have a big job ahead of me if I try to use Shanks to break my record of calculated digits of the MRB constant!! The terms needed to break my record of 3,014,991 digits seems impractical to list! Maybe I just need a new way of computing the Shanks' terms (someway of using known acceleration methods to accumulate them).

POSTED BY: Marvin Ray Burns
POSTED BY: Henrik Schachner

Thank you Henrik! I think the Shanks transformation will help me. It seems, however, Mathematica has an issue with 10's and 100's of iterations of Shanks.

POSTED BY: Marvin Ray Burns

Let g[n] be the error in using f, as you use 10^n terms, as shown in the original post; using the following piece of code:

m = NSum[(-1)^k (k^(1/k) - 1), {k, 1, Infinity}, 
   WorkingPrecision -> 1000, Method -> "AlternatingSigns"];

N::meprec: Internal precision limit $MaxExtraPrecision = 50.` reached while evaluating -E^(2 I Interval[{0,\[Pi]}])+E^(2 I Interval[{0,\[Pi]}]). >>

f[k_] := Log[(k - Log[k])/k];

ClearAll[g]; 
g[n_] := Module[{}, 
  m - (NSum[(-1)^k (k^(1/k) - 1 + f[k]), {k, 1, 10^n}, 
      WorkingPrecision -> Floor[3 n + 100], 
      Method -> "AlternatingSigns"] - 
     NSum[(-1)^k*f[k], {k, 1, Infinity}, 
      WorkingPrecision -> Floor[3 n + 100], 
      Method -> "AlternatingSigns"])]

. You find there is a pattern to g[n]

For[n = 2, n < 20, 
 Print[RootApproximant[
   MantissaExponent[g[10 n]/g[10 (n + 1)]][[1]]]]; n++]

8/27

27/64

64/125

125/216

216/343

343/512

512/729

729/1000

1000/1331

1331/1728

1728/2197

2197/2744

2744/3375

3375/4096

4096/4913

4913/5832

5832/6859

6859/8000

, and find g[10 n]/g[10 (n + 1)] ->n^3/(n+1)^3 *10 to a factor:

For[n = 2, n < 20, Print[
  RootApproximant[
    MantissaExponent[
      g[10 n]/g[10 (n + 1)]
      ][[1]]
    ] - n^3/(n + 1)^3
  ]; n++]

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0

0
POSTED BY: Marvin Ray Burns
POSTED BY: Marvin Ray Burns
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract