Meissel–Mertens constant

Task

Calculate Meissel–Mertens constant up to a precision your language can handle.

Motivation

Analogous to Euler's constant, which is important in determining the sum of reciprocal natural numbers, Meissel-Mertens' constant is important in calculating the sum of reciprocal primes.

Example

We consider the finite sum of reciprocal natural numbers:

1 + 1/2 + 1/3 + 1/4 + 1/5 ... 1/n

this sum can be well approximated with:

log(n) + E

where E denotes Euler's constant: 0.57721...

log(n) denotes the natural natural logarithm of n.

Now consider the finite sum of reciprocal primes:

1/2 + 1/3 + 1/5 + 1/7 + 1/11 ... 1/p

this sum can be well approximated with:

log( log(p) ) + M

where M denotes Meissel-Mertens constant: 0.26149...

See

Details in the Wikipedia article: Meissel–Mertens constant

PARI/GP

Summation method

{
MM(t)=
  my(s=0);
  forprime(p = 2, t,
    s += log(1.-1./p)+1./p
  );
  Euler+s
};

Output:

Running 10^9 summations to get 9 valid digits:

? \p10
   realprecision = 19 significant digits (10 digits displayed)
? MM(1e9)
?
%1 = 0.2614972129
?
? ##
  ***   last result: cpu time 1min, 18,085 ms, real time 1min, 18,094 ms.
?

Wren

Library: Wren-math

Library: Wren-fmt

Wren's only native number type is a 64 bit float (15 or 16 digits accuracy) and it appears from the following results that, no matter how many primes we use, we're not going to be able to improve on 8 digit accuracy for M with this particular number type.

import "./math" for Int
import "./fmt" for Fmt

var euler = 0.57721566490153286
var primes = Int.primeSieve(2.pow(31))
var pc = primes.count
var sum = 0
var c = 0
System.print("Primes added         M")
System.print("------------  --------------")
for (p in primes) {
    var rp = 1/p
    sum = sum + (1-rp).log + rp
    c = c + 1
    if ((c % 1e7) == 0 || c == pc) Fmt.print("$,11d   $0.12f", c, sum + euler)
}

Output:

Primes added         M
------------  --------------
 10,000,000   0.261497213008
 20,000,000   0.261497213008
 30,000,000   0.261497213009
 40,000,000   0.261497213008
 50,000,000   0.261497213009
 60,000,000   0.261497213009
 70,000,000   0.261497213009
 80,000,000   0.261497213009
 90,000,000   0.261497213009
100,000,000   0.261497213009
105,097,565   0.261497213009