Let be a natural number. Denote by the set

of coprime pairs in , and let be their number. It is well-known that ; this was first proven by Ernesto Cesàro.

We are interested in itself. Even though we know that for some , we do not know how this can be chosen.

First note that , where

and that , where is Euler's function. We want to find an explicit lower bound for , which allows us to find an explicit lower bound for itself.

As in the proof of Theorem 330 in Hardy and Wright, we have

where is the Möbius function. Hardy and Wright compare that expression to and show that the error is .

More precisely, it is

First,

Next,

Since

we finally obtain

With a simple computer program, I verified that

for all , and the above bound shows

for . The computer verification took less than one second. We therefore obtain:

We have for all . Equality holds if and only if .

The program for verification of the bound is the following C++ program:

1 #include <iostream> 2 #include <iomanip> 3 #include <cmath> 4 #include <vector> 5 6 // Simple implementation of Euclid's algorithm 7 unsigned long long gcd(unsigned long long a, unsigned long long b) 8 { 9 while (true) 10 { 11 if (a == 0) 12 return b; 13 b %= a; 14 if (b == 0) 15 return a; 16 a %= b; 17 } 18 } 19 20 // Compute all primes up to a bound 21 22 std::vector<unsigned> primes; 23 24 void create_primes(unsigned max) 25 { 26 // TODO: replace with sieve of Eratosthenes 27 primes.push_back(2); 28 for (unsigned p = 3; p <= max; p += 2) 29 { 30 bool prime = true; 31 for (unsigned i = 0; i < primes.size(); ++i) 32 if ((p % primes[i]) == 0) 33 { 34 prime = false; 35 break; 36 } 37 if (prime) 38 primes.push_back(p); 39 } 40 } 41 42 // Compute the Euler phi function (for not too large arguments) 43 unsigned phi(unsigned n) 44 { 45 if (n < 2) 46 return 1; 47 unsigned r = n; 48 bool found_divisor = false; 49 for (unsigned i = 0; i < primes.size(); ++i) 50 if ((n % primes[i]) == 0) 51 { 52 r = r / primes[i] * (primes[i] - 1); 53 found_divisor = true; 54 } 55 if (!found_divisor) 56 --r; 57 return r; 58 } 59 60 int main() 61 { 62 // Compute probabilities up to: 63 unsigned N = 932453; 64 // unsigned N = 10000000; 65 66 // Compute enough primes 67 create_primes(sqrt(N)); 68 std::cout << "Prepared " << primes.size() << " primes\n"; 69 70 // Compute probabilities 71 unsigned long long C = 1; 72 for (unsigned n = 2; n < N; ++n) 73 { 74 C += phi(n); // number of pairs (x, y) with 1 <= x <= y <= n such that gcd(x, y) == 1 75 // Compute rational number representing probability 76 unsigned long long num = 2 * C - 1; 77 unsigned long long den = n; 78 den *= n; 79 unsigned long long g = gcd(num, den); 80 num /= g; 81 den /= g; 82 // Compute floating point approximation of probability 83 long double P = (long double)num / (long double)den; 84 // Output if probability is <= 0.6079 85 if (num <= (unsigned long long)((long double)0.6079 * den)) 86 std::cout << n << ": " << C << ", " << num << "/" << den << " \\approx " << std::setprecision(20) << P << "\n"; 87 } 88 }

On my laptop, it runs in less than one second. For , it requires less than 22 seconds.

## Comments.

An alterantive approach can be found here (Proposition 3.2). There, we analyzed the probability that integers are coprime (as opposed to as here).