Rosetta Code

Meissel–Mertens constant: Difference between revisions

Page Discussion

Meissel–Mertens constant (view source)

Revision as of 19:16, 19 February 2024

48,169 bytes added ,  3 months ago
Added Easylang
(Removed duplicate entry and fixed syntax highlighting.)
(Added Easylang)
 
(48 intermediate revisions by 20 users not shown)
Line 21:
where ''E'' denotes Euler's constant: 0.57721...
''log(n)'' denotes the natural natural logarithm of ''n''.
 
 
Line 37:
;See:
:*   Details in the Wikipedia article:  [https://en.wikipedia.org/wiki/Meissel%E2%80%93Mertens_constant Meissel–Mertens constant]
<br /><br />
 
=={{header|11l}}==
{{trans|Python}}
 
<syntaxhighlight lang="11l">
F primes_up_to_limit(Int limit)
[Int] r
I limit >= 2
r.append(2)
 
V isprime = [1B] * ((limit - 1) I/ 2)
V sieveend = Int(sqrt(limit))
L(i) 0 .< isprime.len
I isprime[i]
Int p = i * 2 + 3
r.append(p)
I i <= sieveend
L(j) ((p * p - 3) >> 1 .< isprime.len).step(p)
isprime[j] = 0B
R r
 
V euler = 0.57721566490153286
V m = 0.0
L(x) primes_up_to_limit(10'000'000)
m += log(1 - (1 / x)) + (1 / x)
 
print(‘MM = #.16’.format(euler + m))
</syntaxhighlight>
 
{{out}}
<pre>
MM = 0.2614972157776471
</pre>
 
=={{header|ALGOL 68}}==
Sieving odd primes only, not using extended precision.<br>
Note that to run this with Algol 68G, a large heap size must be specified, e.g.: by adding <code>-heap 256M</code> to the command line.
<syntaxhighlight lang="algol68">
BEGIN # compute an approximation to the Meissel-Mertens constant #
# construct a sieve of odd primes #
[ 0 : 10 000 000 ]BOOL primes;
BEGIN
FOR i TO UPB primes DO primes[ i ] := TRUE OD;
INT ip := 1;
FOR i WHILE i + ( ip +:= 2 ) <= UPB primes DO
IF primes[ i ] THEN
FOR s FROM i + ip BY ip TO UPB primes DO primes[ s ] := FALSE OD
FI
OD
END;
# sum the reciprocals of the primes #
INT p count := 1;
INT last p := 0;
LONG REAL sum := long ln( 0.5 ) + 0.5;
INT p := 1;
INT p10 := 10;
# Euler's constant from the wikipedia, truncated for LONG REAL #
LONG REAL eulers constant = 0.5772156649015328606 # 0651209008240243104215933593992 #;
FOR i TO UPB primes DO
p +:= 2;
IF primes[ i ] THEN
LONG REAL rp = 1 / LENG p;
sum +:= long ln( 1 - rp ) + rp;
p count +:= 1;
last p := p;
IF p count = p10 THEN
print( ( "after ", whole( p count, -8 ), " primes, the approximation is: "
, fixed( sum + eulers constant, -14, 12 )
, ", last prime considered: ", whole( last p, 0 )
, newline
)
);
p10 := IF p10 < 1 000 000 THEN p10 * 10 ELSE p10 + 1 000 000 FI
FI
FI
OD;
print( ( "after ", whole( p count, -8 ), " primes, the approximation is: "
, fixed( sum + eulers constant, -14, 12 )
, ", last prime considered: ", whole( last p, 0 )
, newline
)
)
END
</syntaxhighlight>
{{out}}
<pre>
after 10 primes, the approximation is: 0.265160104017, last prime considered: 29
after 100 primes, the approximation is: 0.261624626173, last prime considered: 541
after 1000 primes, the approximation is: 0.261503563617, last prime considered: 7919
after 10000 primes, the approximation is: 0.261497594498, last prime considered: 104729
after 100000 primes, the approximation is: 0.261497238454, last prime considered: 1299709
after 1000000 primes, the approximation is: 0.261497214692, last prime considered: 15485863
after 1270607 primes, the approximation is: 0.261497214255, last prime considered: 19999999
</pre>
 
=={{header|Arturo}}==
<syntaxhighlight lang="arturo">meisselMertens: function [depth][
Euler: 0.57721566490153286
m: (1//2) + ln 1-1//2
loop range.step:2 3 depth 'x ->
if prime? x ->
m: m + (1//x) + ln 1-1//x
 
return m + Euler
]
 
print meisselMertens 10000000</syntaxhighlight>
 
{{out}}
 
<pre>0.2614972157776471</pre>
 
=={{header|BASIC}}==
==={{header|BASIC256}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="freebasic">Euler = 0.5772156649
m = 0
for x = 2 to 1e6 # more prime numbers do not add more precision
if isPrime(x) then m += log(1-(1/x)) + (1/x)
next x
print "MM = "; Euler + m
print Euler
end
 
function isPrime(v)
if v < 2 then return False
if v mod 2 = 0 then return v = 2
if v mod 3 = 0 then return v = 3
d = 5
while d * d <= v
if v mod d = 0 then return False else d += 2
end while
return True
end function</syntaxhighlight>
{{out}}
<pre>MM = 0.26149724673</pre>
 
==={{header|Run BASIC}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="lb">function isPrime(n)
if n < 2 then isPrime = 0 : goto [exit]
if n = 2 then isPrime = 1 : goto [exit]
if n mod 2 = 0 then isPrime = 0 : goto [exit]
isPrime = 1
for i = 3 to int(n^.5) step 2
if n mod i = 0 then isPrime = 0 : goto [exit]
next i
[exit]
end function
 
e = 0.5772156
 
for x = 2 to 100000 ' more prime numbers do not add more precision
if isPrime(x) then m = m + log(1-(1/x)) + (1/x)
next x
print "MM = "; e + m</syntaxhighlight>
{{out}}
<pre>MM = 0.261274</pre>
 
==={{header|PureBasic}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="PureBasic">Procedure isPrime(v.i)
If v <= 1 : ProcedureReturn #False
ElseIf v < 4 : ProcedureReturn #True
ElseIf v % 2 = 0 : ProcedureReturn #False
ElseIf v < 9 : ProcedureReturn #True
ElseIf v % 3 = 0 : ProcedureReturn #False
Else
Protected r = Round(Sqr(v), #PB_Round_Down)
Protected f = 5
While f <= r
If v % f = 0 Or v % (f + 2) = 0
ProcedureReturn #False
EndIf
f + 6
Wend
EndIf
ProcedureReturn #True
EndProcedure
 
OpenConsole()
Euler.d = 0.5772156649 ;0153286
 
For x.i = 2 To 1e8
If isPrime(x)
m.d = m + Log(1-(1/x)) + (1/x)
EndIf
Next x
PrintN("MM = " + StrD(Euler + m))
PrintN(#CRLF$ + "--- terminado, pulsa RETURN---"): Input()
CloseConsole()
</syntaxhighlight>
{{out}}
<pre>MM = 0.2614972129</pre>
 
==={{header|True BASIC}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="qbasic">FUNCTION isPrime (n)
IF n = 2 THEN
LET isPrime = 1
ELSEIF n <= 1 OR REMAINDER(n, 2) = 0 THEN
LET isPrime = 0
ELSE
LET isPrime = 1
FOR i = 3 TO INT(SQR(n)) STEP 2
IF REMAINDER(n, i) = 0 THEN
LET isPrime = 0
EXIT FUNCTION
END IF
NEXT i
END IF
END FUNCTION
 
LET e = .5772156649
 
FOR x = 2 to 1e6 !more prime numbers do not add more precision
IF isPrime(x) = 1 THEN
LET m = m + LOG(1-(1/x)) + (1/x)
END IF
NEXT x
PRINT "MM ="; e + m
END</syntaxhighlight>
{{out}}
<pre>MM = 0.26149725</pre>
 
==={{header|Yabasic}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="freebasic">e = 0.5772156
 
for x = 2 to 1e6 // more prime numbers do not add more precision
if isPrime(x) m = m + log(1-(1/x)) + (1/x)
next x
print "MM = ", e + m
end
 
sub isPrime(v)
if v < 2 return False
if mod(v, 2) = 0 return v = 2
if mod(v, 3) = 0 return v = 3
d = 5
while d * d <= v
if mod(v, d) = 0 then return False else d = d + 2 : fi
wend
return True
end sub</syntaxhighlight>
{{out}}
<pre>MM = 0.261497</pre>
 
=={{header|C++}}==
<syntaxhighlight lang="c++">
#include <cmath>
#include <cstdint>
#include <iomanip>
#include <iostream>
#include <vector>
 
std::vector<double> list_prime_reciprocals(const int32_t& limit) {
const int32_t half_limit = ( limit % 2 == 0 ) ? limit / 2 : 1 + limit / 2;
std::vector<bool> composite(half_limit);
for ( int32_t i = 1, p = 3; i < half_limit; p += 2, ++i ) {
if ( ! composite[i] ) {
for ( int32_t a = i + p; a < half_limit; a = a + p ) {
composite[a] = true;
}
}
}
 
std::vector<double> result(composite.size());
result[0] = 0.5;
for ( int32_t i = 1, p = 3; i < half_limit; p += 2, ++i ) {
if ( ! composite[i] ) {
result.emplace_back(1.0 / p);
}
}
return result;
}
 
int main() {
std::vector<double> prime_reciprocals = list_prime_reciprocals(100'000'000);
const double euler = 0.577'215'664'901'532'861;
double sum = 0.0;
for ( double reciprocal : prime_reciprocals ) {
sum += reciprocal + log(1.0 - reciprocal);
}
 
const double meissel_mertens = euler + sum;
std::cout << "The Meissel-Mertens constant = " << std::setprecision(8) << meissel_mertens << std::endl;
}
</syntaxhighlight>
{{ out }}
<pre>
The Meissel-Mertens constant = 0.26149721
</pre>
 
=={{header|Delphi}}==
{{works with|Delphi|6.0}}
{{libheader|SysUtils,StdCtrls}}
This program doesn't test the limits of numbers in Delphi. The limit is the processing time. The time to process 10^7 is 3 minutes. The time to process 10^8 is nearly two hours. As a result, process beyond 10^8 is pretty much impossible and a different strategy would be required
 
<syntaxhighlight lang="Delphi">
function IsPrime(N: int64): boolean;
{Fast, optimised prime test}
var I,Stop: int64;
begin
if (N = 2) or (N=3) then Result:=true
else if (n <= 1) or ((n mod 2) = 0) or ((n mod 3) = 0) then Result:= false
else
begin
I:=5;
Stop:=Trunc(sqrt(N+0.0));
Result:=False;
while I<=Stop do
begin
if ((N mod I) = 0) or ((N mod (I + 2)) = 0) then exit;
Inc(I,6);
end;
Result:=True;
end;
end;
 
 
function GetNextPrime(var Start: integer): integer;
{Get the next prime number after Start}
{Start is passed by "reference," so the
{original variable is incremented}
begin
repeat Inc(Start)
until IsPrime(Start);
Result:=Start;
end;
 
 
function MeisselMertens(Depth: integer; Prog: TProgress): extended;
{Calculate MM value up a certain Depth}
var I,P: integer;
const Euler = 0.57721566490153286;
begin
Result:=0;
P:=1;
for I:=1 to Depth do
begin
P:=GetNextPrime(P);
Result:=Result+Ln(1-(1/P)) + (1/P);
if Assigned(Prog) and ((I mod 10000)=0) then
HandleProgress(MulDiv(100,I,Depth));
end;
Result:=Result+Euler;
end;
 
 
procedure ShowMeisselMertens(Memo: TMemo; Prog: TProgress);
var I,IT,Digits: integer;
var M,Last: extended;
begin
Memo.Lines.Add('Primes Digits M');
Memo.Lines.Add('-----------------------------------------');
Last:=0;
IT:=1;
{Calculate MM to specified Power of 10}
for I:=1 to 7 do
begin
IT:=IT*10;
M:=MeisselMertens(IT,Prog);
{Calculated Digits of accuracy}
Digits:=Trunc(abs(Log(abs(M-Last))));
Memo.Lines.Add(Format('10^%2d %7d %25.18f',[I,Digits,M]));
Last:=M
end;
end;
 
 
</syntaxhighlight>
{{out}}
<pre>
Primes Digits M
-----------------------------------------
10^ 1 0 0.265160104017311294
10^ 2 2 0.261624626172936147
10^ 3 3 0.261503563616850571
10^ 4 5 0.261497594498432488
10^ 5 6 0.261497238454297260
10^ 6 7 0.261497214692305277
10^ 7 8 0.261497212987251183
10^ 8 9 0.261497212858582276
 
</pre>
 
=={{header|Dart}}==
{{trans|FreeBASIC}}
<syntaxhighlight lang="dart">import 'dart:math';
 
bool isPrime(var n) {
if (n <= 1) return false;
if (n == 2) return true;
for (var i = 2; i <= sqrt(n); i++) if (n % i == 0) return false;
return true;
}
 
void main() {
const double euler = 0.57721566490153286;
 
double m = 0.0;
for (var x = 2; x <= 1e8; x++)
if (isPrime(x)) m += log(1 - (1 / x)) + (1 / x);
 
print('MM = ${euler + m}');
}</syntaxhighlight>
{{out}}
<pre>MM = 0.2614972131057144</pre>
 
=={{header|EasyLang}}==
{{trans|BASIC256}}
<syntaxhighlight>
fastfunc isprim num .
if num mod 2 = 0
if num = 2
return 1
.
return 0
.
i = 3
while i <= sqrt num
if num mod i = 0
return 0
.
i += 2
.
return 1
.
func log x .
return log10 x / log10 2.7182818284590452354
.
euler = 0.5772156649
for x = 2 to 1e6
if isprim x = 1
m += log (1 - (1 / x)) + (1 / x)
.
.
numfmt 11 0
print "mm = " & euler + m
</syntaxhighlight>
 
{{out}}
<pre>
mm = 0.26149724673
</pre>
 
=={{header|FreeBASIC}}==
On my CPU i5 it takes about 3.5 minutes
<syntaxhighlight lang="freebasic">'#include "isprime.bas"
 
Const As Double Euler = 0.57721566490153286
 
Dim As Double m = 0
For x As Ulongint = 2 To 1e8
If isPrime(x) Then m += Log(1-(1/x)) + (1/x)
Next x
Print "MM ="; Euler + m
Sleep</syntaxhighlight>
{{out}}
<pre>MM = 0.2614972131104154</pre>
 
=={{header|Go}}==
{{trans|Wren}}
{{libheader|Go-rcu}}
<syntaxhighlight lang="go">package main
 
import (
"fmt"
"math"
"rcu"
)
 
func contains(a []int, f int) bool {
for _, e := range a {
if e == f {
return true
}
}
return false
}
 
func main() {
const euler = 0.57721566490153286
primes := rcu.Primes(1 << 31)
pc := len(primes)
sum := 0.0
fmt.Println("Primes added M")
fmt.Println("------------ --------------")
for i, p := range primes {
rp := 1.0 / float64(p)
sum += math.Log(1.0-rp) + rp
c := i + 1
if (c%1e7) == 0 || c == pc {
fmt.Printf("%11s %0.12f\n", rcu.Commatize(c), sum+euler)
}
}
}</syntaxhighlight>
 
{{out}}
<pre>
Primes added M
------------ --------------
10,000,000 0.261497212987
20,000,000 0.261497212912
30,000,000 0.261497212889
40,000,000 0.261497212878
50,000,000 0.261497212871
60,000,000 0.261497212867
70,000,000 0.261497212864
80,000,000 0.261497212862
90,000,000 0.261497212861
100,000,000 0.261497212859
105,097,565 0.261497212858
</pre>
 
=={{header|J}}==
'''Implementation'''
<syntaxhighlight lang="j">Euler=: 0.57721566490153286
MM=: {{ Euler + +/ (+ ^.@-.)@% p: i. y }}</syntaxhighlight>
'''Task example'''
<syntaxhighlight lang="j"> 0j13 ": MM 1e8
0.2614972128591</syntaxhighlight>
 
That said, a more literal implementation, like
 
<syntaxhighlight lang=J>mm=: (% 10x(^#)10#.inv]) 26149721284764278375542683860869585905156664826119920619206421392x</syntaxhighlight>
 
would be more precise:
 
<syntaxhighlight lang=J> 0j65":mm
0.26149721284764278375542683860869585905156664826119920619206421392</syntaxhighlight>
 
=={{header|Java}}==
<syntaxhighlight lang="java">
import java.util.ArrayList;
import java.util.BitSet;
import java.util.List;
 
public final class MeisselMertensConstant {
public static void main(String[] aArgs) {
List<Double> primeReciprocals = listPrimeReciprocals(1_000_000_000);
final double euler = 0.577_215_664_901_532_861;
double sum = 0.0;
for ( double reciprocal : primeReciprocals ) {
sum += reciprocal + Math.log(1.0 - reciprocal);
}
final double meisselMertens = euler + sum;
System.out.println(String.format("%s%.9f", "The Meissel-Mertens constant = ", meisselMertens));
}
private static List<Double> listPrimeReciprocals(int aLimit) {
BitSet sieve = new BitSet(aLimit + 1);
sieve.set(2, aLimit + 1);
for ( int i = 2; i <= Math.sqrt(aLimit); i = sieve.nextSetBit(i + 1) ) {
for ( int j = i * i; j <= aLimit; j += i ) {
sieve.clear(j);
}
}
List<Double> result = new ArrayList<Double>(sieve.cardinality());
for ( int i = 2; i >= 0; i = sieve.nextSetBit(i + 1) ) {
result.add(1.0 / i);
}
return result;
}
 
}
</syntaxhighlight>
{{ out }}
<pre>
The Meissel-Mertens constant = 0.261497213
</pre>
 
=={{header|Julia}}==
Off by one in the 11th digit after 10^8 primes.
<syntaxhighlight lang="julia">
using Base.MathConstants # sets constant γ = 0.5772156649015...
using Primes
 
""" Approximate the Meissel-Mertons constant. """
function meissel_mertens(iterations = 100_000_000)
return mapreduce(p ->(d = 1/p; log(1 - d) + d), +, primes(prime(iterations))) + γ
end
 
@show meissel_mertens(100_000_000) # meissel_mertens(100000000) = 0.2614972128591237
</syntaxhighlight>
 
=={{header|Maxima}}==
<syntaxhighlight lang="maxima">
meissel_mertens:%gamma+lsum(log(1-(1/i))+(1/i),i,primes(2,10000000)),numer;
</syntaxhighlight>
{{out}}
<pre>
0.2614972157776471
</pre>
 
=={{header|Nim}}==
{{trans|Wren}}
<syntaxhighlight lang="Nim">import std/[math, strformat, strutils]
 
proc initPrimes(N: static int): seq[int] =
## Initialize the list of primes.
 
const M = 2 * N - 1
var composite = newSeq[bool](N)
composite[0] = true # 1 is not prime.
 
# Conversions from index to value and value to index.
template index(n: int): int = (n - 1) shr 1
template value(idx: int): int = idx shl 1 + 1
 
# Fill the sieve.
var n = 3
while n * n <= M:
if not composite[n.index]:
for k in countup(n * n, M, 2 * n):
composite[k.index] = true
inc n, 2
 
# Build list of primes.
result = @[2]
for idx in 0..composite.high:
if not composite[idx]:
result.add idx.value
 
 
const N = 2^30
let primes = initPrimes(N)
 
echo "Primes added M"
echo "──────────── ──────────────"
const γ = 0.57721566490153286 # Euler–Mascheroni constant.
let primeCount = primes.len
var sum = 0.0
var count = 0
for p in primes:
let rp = 1 / p
sum += ln(1 - rp) + rp
inc count
if count mod 10_000_000 == 0 or count == primeCount:
echo &"{insertSep($count):>11} {sum+γ:.12}"
</syntaxhighlight>
 
{{out}}
<pre>Primes added M
──────────── ──────────────
10_000_000 0.261497212987
20_000_000 0.261497212912
30_000_000 0.261497212889
40_000_000 0.261497212878
50_000_000 0.261497212871
60_000_000 0.261497212867
70_000_000 0.261497212864
80_000_000 0.261497212862
90_000_000 0.261497212861
100_000_000 0.261497212859
105_097_565 0.261497212858
</pre>
 
=={{header|PARI/GP}}==
Line 50 ⟶ 713:
};</syntaxhighlight>
{{out}}
10 valid digits, prime summation up to 10^9.
Running 10^9 summations to get 9 valid digits:
<pre>? \p10p 12
realprecision = 19 significant digits (1012 digits displayed)
? MM(1e9)
?
%1 = 0.2614972129261497212874
?
? ##
*** last result: cpu time 1min, 18,085 ms, real time 1min, 18,094 ms.
?</pre>
===Analytic method===
The Analytic method requires high precision calculation of Riemann zeta function.
<syntaxhighlight lang="parigp">{
Meissel_Mertens(d)=
default(realprecision, d);
my(prec = default(realprecision), z = 0, y = 0, q);
forprime(p = 2 , 7,
z += log(1.-1./p)+1./p
);
for(k = 2, prec,
q = 1;
forprime(p = 2, 7,
q *= 1.-p^-k
);
y += moebius(k)*log(zeta(k)*q)/k
);
Euler+z+y
};</syntaxhighlight>
{{out}}
1000 valid digits. For some reason the last digit is in some rare cases (rounded) wrong. Whatever, the last digit is per definition the one with the lowest weight.
<pre>? Meissel_Mertens(1001)
%1 = 0.26149721284764278375542683860869585905156664826119920619206421392492451089736820971414263143424665105161772887648602199778339032427004442454348740197238640666194955709392581712774774211985258807266272064144464232590023543105177232173925663229980314763831623758149059290382284758265972363422015971458785446941586825460538918007031787714156680620570605257601785334398970354507934530971953511716888598019955346947142883673537117910619342522616975101911159537244599605203558051780574237201332999961769676911386909654186249097435916294862238555389898241954857937738258646582212506260380084370067541379219020626760709633535981989783010762417792511961619355361391684002933280522289185167238258837930443067100391254985761418536020400457460311825670423438456551983202200477824746954606715454777572171338072595463648319687279859427306787306509669454587505942593547068846408425666008833035029366514525328713339609172639368543886291288200447611698748441593459920236225093315001729474600911978170842383659092665509
? ##
*** last result: cpu time 283 ms, real time 284 ms.</pre>
 
=={{header|Perl}}==
{{libheader|ntheory}}
<syntaxhighlight lang="perl" line>use v5.36;
use ntheory qw(forprimes);
 
my $s;
forprimes { $s += log(1 - 1/$_)+1/$_ } 1e9;
say my $result = $s + .57721566490153286;</syntaxhighlight>
{{out}}
<pre>0.261497212871049</pre>
 
=={{header|Phix}}==
Converges very slowly, I started running out of memory (from generating so many primes) before it showed any signs of getting stuck.
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> <span style="color: #000080;font-style:italic;">-- (but perhaps a bit too slow)</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">mmc</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0.2614972128476427837554268386086958590516</span><span style="color: #0000FF;">,</span>
<span style="color: #000000;">smmc</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"0.2614972128476427837554268386086958590516"</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0.57721566490153286</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">dpa</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">p10</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">pn</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">adp</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">0</span>
<span style="color: #004080;">string</span> <span style="color: #000000;">st</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"0"</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">fmt</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"%.0f"</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Primes added M\n"</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"------------ --------------\n"</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">while</span> <span style="color: #000000;">adp</span><span style="color: #0000FF;"><=</span><span style="color: #000000;">10</span> <span style="color: #008080;">do</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">rp</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">/</span><span style="color: #7060A8;">get_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">pn</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">t</span> <span style="color: #0000FF;">+=</span> <span style="color: #7060A8;">log</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">-</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span> <span style="color: #0000FF;">+</span> <span style="color: #000000;">rp</span>
<span style="color: #008080;">if</span> <span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">-</span><span style="color: #000000;">mmc</span><span style="color: #0000FF;">)<</span><span style="color: #000000;">dpa</span> <span style="color: #008080;">then</span>
<span style="color: #004080;">string</span> <span style="color: #000000;">ft</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fmt</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">trunc</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">*</span><span style="color: #000000;">p10</span><span style="color: #0000FF;">)/</span><span style="color: #000000;">p10</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (as below)</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">ft</span><span style="color: #0000FF;">=</span><span style="color: #000000;">st</span> <span style="color: #008080;">then</span>
<span style="color: #000080;font-style:italic;">--
-- We have to artifically calculate a "truncated t", aka tt,
-- to prevent say 0.2..299[&gt;5..] being automatically rounded
-- by printf() to 0.2..300, otherwise it just "looks wrong".
--</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">tt</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">trunc</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">*</span><span style="color: #000000;">1e12</span><span style="color: #0000FF;">)/</span><span style="color: #000000;">1e12</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%,11d %0.12f (accurate to %d d.p.)\n"</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">pn</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">tt</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">adp</span><span style="color: #0000FF;">})</span>
<span style="color: #000000;">adp</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
<span style="color: #000000;">fmt</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"%%.%df"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">adp</span><span style="color: #0000FF;">)</span>
<span style="color: #000000;">st</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">smmc</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">..</span><span style="color: #000000;">2</span><span style="color: #0000FF;">+</span><span style="color: #000000;">adp</span><span style="color: #0000FF;">]</span>
<span style="color: #000000;">dpa</span> <span style="color: #0000FF;">/=</span> <span style="color: #000000;">10</span>
<span style="color: #000000;">p10</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">10</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #000000;">pn</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"(actual value 0.26149721284764278375542683860869)\n"</span><span style="color: #0000FF;">)</span>
<!--</syntaxhighlight>-->
{{out}}
<small>(Couldn't be bothered with making it an inner loop to upgrade the "accurate to 4dp" to 5dp)</small>
<pre>
Primes added M
------------ --------------
1 0.384068484341 (accurate to 0 d.p.)
3 0.288793158252 (accurate to 1 d.p.)
6 0.269978901636 (accurate to 2 d.p.)
38 0.261990332075 (accurate to 3 d.p.)
1,940 0.261499998883 (accurate to 4 d.p.)
1,941 0.261499997116 (accurate to 5 d.p.)
5,471 0.261497999951 (accurate to 6 d.p.)
34,891 0.261497299999 (accurate to 7 d.p.)
303,447 0.261497219999 (accurate to 8 d.p.)
9,246,426 0.261497212999 (accurate to 9 d.p.)
24,304,615 0.261497212899 (accurate to 10 d.p.)
(actual value 0.26149721284764278375542683860869)
</pre>
===analytical/gmp===
{{trans|PARI/GP}}
<!--<syntaxhighlight lang="phix">-->
<span style="color: #008080;">without</span> <span style="color: #008080;">js</span> <span style="color: #000080;font-style:italic;">-- no mpfr_zeta[_ui]() in mpfr.js, as yet, or mpfr_log() or mpfr_const_euler(), for that matter</span>
<span style="color: #008080;">include</span> <span style="color: #004080;">mpfr</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">moebius</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span> <span style="color: #008080;">return</span> <span style="color: #000000;">1</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #004080;">sequence</span> <span style="color: #000000;">f</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">prime_factors</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #004600;">true</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">2</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">f</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">f</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">f</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span> <span style="color: #008080;">then</span> <span style="color: #008080;">return</span> <span style="color: #000000;">0</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">return</span> <span style="color: #7060A8;">iff</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">odd</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">f</span><span style="color: #0000FF;">))?-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">:+</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">Meissel_Mertens</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_set_default_precision</span><span style="color: #0000FF;">(-</span><span style="color: #000000;">d</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (d decimal places)</span>
<span style="color: #004080;">mpfr</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">res</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">q</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpfr_inits</span><span style="color: #0000FF;">(</span><span style="color: #000000;">5</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">p</span> <span style="color: #008080;">in</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span><span style="color: #000000;">5</span><span style="color: #0000FF;">,</span><span style="color: #000000;">7</span><span style="color: #0000FF;">}</span> <span style="color: #008080;">do</span>
<span style="color: #000080;font-style:italic;">-- z += log(1-1/p)+1/p</span>
<span style="color: #7060A8;">mpfr_set_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_div_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">p</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_add</span><span style="color: #0000FF;">(</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_si_sub</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_log</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_add</span><span style="color: #0000FF;">(</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">k</span><span style="color: #0000FF;">=</span><span style="color: #000000;">2</span> <span style="color: #008080;">to</span> <span style="color: #000000;">d</span> <span style="color: #008080;">do</span> <span style="color: #000080;font-style:italic;">-- (see note)</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">m</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">moebius</span><span style="color: #0000FF;">(</span><span style="color: #000000;">k</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">if</span> <span style="color: #000000;">m</span> <span style="color: #008080;">then</span>
<span style="color: #7060A8;">mpfr_set_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">q</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">for</span> <span style="color: #000000;">p</span> <span style="color: #008080;">in</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span><span style="color: #000000;">5</span><span style="color: #0000FF;">,</span><span style="color: #000000;">7</span><span style="color: #0000FF;">}</span> <span style="color: #008080;">do</span>
<span style="color: #000080;font-style:italic;">-- q *= 1-power(p,-k)</span>
<span style="color: #7060A8;">mpfr_set_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">p</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_pow_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">k</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_si_sub</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">q</span><span style="color: #0000FF;">,</span><span style="color: #000000;">q</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #000080;font-style:italic;">-- y += moebius(k)*log(zeta(k)*q)/k</span>
<span style="color: #7060A8;">mpfr_zeta_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">k</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">q</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_log</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_div_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">k</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_mul_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">,</span><span style="color: #000000;">m</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_add</span><span style="color: #0000FF;">(</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">,</span><span style="color: #000000;">rp</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<span style="color: #000080;font-style:italic;">-- res := EULER+z+y</span>
<span style="color: #000000;">mpfr_const_euler</span><span style="color: #0000FF;">(</span><span style="color: #000000;">res</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_add</span><span style="color: #0000FF;">(</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">,</span><span style="color: #000000;">y</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">mpfr_add</span><span style="color: #0000FF;">(</span><span style="color: #000000;">res</span><span style="color: #0000FF;">,</span><span style="color: #000000;">res</span><span style="color: #0000FF;">,</span><span style="color: #000000;">z</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">res</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #004080;">mpfr</span> <span style="color: #000000;">m</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">Meissel_Mertens</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1001</span><span style="color: #0000FF;">)</span>
<span style="color: #0000FF;">?</span><span style="color: #7060A8;">shorten</span><span style="color: #0000FF;">(</span><span style="color: #000000;">mpfr_get_str</span><span style="color: #0000FF;">(</span><span style="color: #000000;">m</span><span style="color: #0000FF;">,</span><span style="color: #000000;">10</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1001</span><span style="color: #0000FF;">))</span>
<!--</syntaxhighlight>-->
{{out}}
Agrees with PARI/GP except for the very last digit being 8 instead of 9. Note that playing with precision and/or iterations (which seems to differ quite wildly) gave utterly incorrect results... not that I actually comprehend what the algorithm is doing. The 1,003 includes 2 from the "0.", fairly obviously.
<pre>
"0.261497212847642783...70842383659092665508 (1,003 digits)"
</pre>
 
=={{header|Python}}==
{{trans|FreeBASIC}}
<syntaxhighlight lang="python">#!/usr/bin/python
 
from math import log
 
def isPrime(n):
for i in range(2, int(n**0.5) + 1):
if n % i == 0:
return False
return True
 
 
if __name__ == '__main__':
Euler = 0.57721566490153286
m = 0
for x in range(2, 10_000_000):
if isPrime(x):
m += log(1-(1/x)) + (1/x)
 
print("MM =", Euler + m)</syntaxhighlight>
{{out}}
<pre>MM = 0.26149721577764207</pre>
 
=={{header|Raku}}==
{{trans|FreeBASIC}}
<syntaxhighlight lang="raku" line># 20221011 Raku programming solution
 
my $s;
.is-prime and $s += log(1-1/$_)+1/$_ for 2 .. 10**8;
say $s + .57721566490153286</syntaxhighlight>
{{out}}
<pre>0.26149721310571383</pre>
 
=={{header|RPL}}==
<code>'''BPRIM?'''</code> is defined at [[Primality by trial division#RPL]]
{{works with|Halcyon Calc|4.2.8}}
≪ → n
≪ 0.5 3 n '''FOR''' j
'''IF''' j R→B '''BPRIM? THEN''' j INV + '''END''' 2 '''STEP'''
n LN LN -
≫ ≫ ''''MMCON'''' STO
 
10000 '''MMCON'''
{{out}}
<pre>
1: 0.262733140866
</pre>
 
=={{header|Sidef}}==
<syntaxhighlight lang="ruby">var sum = 0
1e7.each_prime {|p|
with (1f/p) {|t|
sum += (log(1 - t) + t)
}
}
say sum+Num.EulerGamma</syntaxhighlight>
{{out}}
<pre>
0.26149721577767111119422410228297206467931376306
</pre>
 
=={{header|Wren}}==
===Summation method===
{{libheader|Wren-math}}
{{libheader|Wren-fmt}}
It will be seen that this is converging to the correct answer though we'd need to add a lot more primes to obtain a valid 11th digit.
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.
<syntaxhighlight lang="ecmascriptwren">import "./math" for Int
import "./fmt" for Fmt
 
Line 77 ⟶ 953:
for (p in primes) {
var rp = 1/p
sum = sum + (1-rp).log + rp + sum
c = c + 1
if ((c % 1e7) == 0 || c == pc) Fmt.print("$,11d $0.12f", c, sum + euler)
Line 86 ⟶ 962:
Primes added M
------------ --------------
10,000,000 0.261497213008261497212987
20,000,000 0.261497213008261497212912
30,000,000 0.261497213009261497212889
40,000,000 0.261497213008261497212878
50,000,000 0.261497213009261497212871
60,000,000 0.261497213009261497212867
70,000,000 0.261497213009261497212864
80,000,000 0.261497213009261497212862
90,000,000 0.261497213009261497212861
100,000,000 0.261497213009261497212859
105,097,565 0.261497213009261497212858
</pre>
===Analytic method===
{{trans|PARI/GP}}
{{libheader|Wren-gmp}}
This agrees with the Phix entry that the 1,000th digit after the decimal point is '8' after rounding (according to OEIS it's actually '7' but the next one is '7' also).
<syntaxhighlight lang="wren">import "./gmp" for Mpf
import "./math" for Int
import "./fmt" for Fmt
 
var isSquareFree = Fn.new { |n|
var i = 2
while (i * i <= n) {
if (n%(i*i) == 0) return false
i = (i > 2) ? i + 2 : i + 1
}
return true
}
 
var mu = Fn.new { |n|
if (n < 1) Fiber.abort("Argument must be a positive integer")
if (n == 1) return 1
var sqFree = isSquareFree.call(n)
var factors = Int.primeFactors(n)
if (sqFree && factors.count % 2 == 0) return 1
if (sqFree) return -1
return 0
}
 
var meisselMertens = Fn.new { |d|
Mpf.defaultPrec = d
var z = Mpf.zero
var y = Mpf.zero
var r = Mpf.new()
var q = Mpf.new()
var t = Mpf.new()
var m = Mpf.new()
for (p in [2, 3, 5, 7]) {
r.setUi(p).inv
t.uiSub(1, r).log.add(r)
z.add(t, z)
}
for (k in 2..d) {
q.setUi(1)
for (p in [2, 3, 5, 7]) {
r.setUi(p).inv
t.uiSub(1, r.pow(k))
q.mul(t)
}
m.setSi(mu.call(k))
t.zetaUi(k).mul(q).log.mul(m).div(k)
y.add(t, y)
}
return Mpf.euler.add(z).add(y)
}
 
Fmt.print("$20a", meisselMertens.call(3300).toString(1001))</syntaxhighlight>
 
{{out}}
<pre>
0.261497212847642783...70842383659092665508
</pre>
 
=={{header|XPL0}}==
Runs in 98 seconds on Pi4.
<syntaxhighlight lang "XPL0">func IsPrime(N); \Return 'true' if N is prime
int N, D;
[if N < 2 then return false;
if (N&1) = 0 then return N = 2;
if rem(N/3) = 0 then return N = 3;
D:= 5;
while D*D <= N do
[if rem(N/D) = 0 then return false;
D:= D+2;
if rem(N/D) = 0 then return false;
D:= D+4;
];
return true;
];
 
def Euler = 0.57721566490153286;
real M; int P;
[M:= 0.;
for P:= 2 to 100_000_000 do
if IsPrime(P) then
M:= M + Ln(1. - 1./float(P)) + 1./float(P);
Format(1, 16);
Text(0, "MM = "); RlOut(0, Euler + M); CrLf(0);
]</syntaxhighlight>
{{out}}
<pre>
MM = 0.2614972131104150
</pre>
 
=={{header|Mathematica}} / {{header|Wolfram Language}}==
<syntaxhighlight lang="mathematica">PrimeNumbers = Select[Range[100000000], PrimeQ[#] &]; (*all primes in the first 100 000 000 numbers, this takes a toll on my computer's CPU and RAM*)
MM = N[Total[Log[1 - 1/PrimeNumbers] + 1/PrimeNumbers] + EulerGamma, 10] (*Calculating it up to a precision of 10, this is correct up to 8 digits*)
AnalyticMMto305 = N[EulerGamma + Sum[MoebiusMu[n]/n Log[Zeta[n]], {n, 2, 1000}], 1000] (*Precise up to 305 digits*)
AnalyticMM = N[EulerGamma + Sum[MoebiusMu[n]/n Log[Zeta[n]], {n, 2, 10000}], 1001] (*Precise up to at least 1000 digits*)</syntaxhighlight>
 
{{out}}<pre>0.2614972131
0.2614972128476427837554268386086958590515666482611992061920642139249245108973682097141426314342466510516177288764860219977833903242700444245434874019723864066619495570939258171277477421198525880726627206414446423259002354310517723217392566322998031476383162375814905929038228475826597236342201597145878545286546864785868143588282690675811212349253390275898516648963722862532834592270193787663814788872688899729598789511243611212105643183719483028350694889140412496558830385604700978671612478906659270558671778566494523172201279234383779159082832574101983245266507469788929642499636163109342230499137582782542115136325204939772977424363327903550519767079810438569344262666067369673777630199473914034464437524475856770642386781282897163202273071362276688321001599406964213181471899136557333922267566346325514106895289284293459552062339943989806700523705145085201774716630458161117596603698964857884739475619061322196278674275162487770072794776455644351959504021182511495417580620912301758280597008866059
0.26149721284764278375542683860869585905156664826119920619206421392492451089736820971414263143424665105161772887648602199778339032427004442454348740197238640666194955709392581712774774211985258807266272064144464232590023543105177232173925663229980314763831623758149059290382284758265972363422015971458785446941586825460538918007031787714156680620570605257601785334398970354507934530971953511716888598019955346947142883673537117910619342522616975101911159537244599605203558051780574237201332999961769676911386909654186249097435916294862238555389898241954857937738258646582212506260380084370067541379219020626760709633535981989783010762417792511961619355361391684002933280522289185167238258837930443067100391254985761418536020400457460311825670423438456551983202200477824746954606715454777572171338072595463648319687279859427306787306509669454587505942593547068846408425666008833035029366514525328713339609172639368543886291288200447611698748441593459920236225093315001729474600911978170842383659092665508
</pre>
Chkas
1,985

edits

Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/331171...358257"