How do I calculate chances of finding a prime, while testing multiple exponents?

[drkirkby]

First, my background is engineering, not mathematics, and I'm not too hot on maths.



I'm trying to understand how mprime calculates the chances of finding a prime in one of the exponents I'm testing at the minute. I'm testing 32 exponents between M110274583 and M110904847. This is what mprime shows.

Below is a report on the work you have queued and any expected completion
dates.
<snip>
The chance that one of the 32 exponents you are testing will yield a
Mersenne prime is about 1 in 26429.

So that's a range of numbers between
2^110274583 -1 = 1.74x10^33195957
2^110904847 -1 = 4.09x10^33385685

The mean of those is around 2x10^33385685
The probability of any one of them being prime is then

probability_of_one_being_prime=1/log(n)=1/log(2x10^33385685)=1.3x10^-8

The probability of one being composite is then

probability_of_one_being_composite=1-probability_of_one_being_prime
=0.9999999869915960129554356230753

Testing 32, I would expect the probability of them all being composite to be
probability_of_all_composite=probability_of_one_being_composite^32
=0.99999958373115634697586958993

Probability of at least one being prime =

1-probability_of_all_composite
1-0.99999958373115634697586958993
=4.2x10^-7

So I would expect the chance of finding a prime to be one in
1/4.2x10^-7 = 2.4x10^6.

Yet mprime says the chance of finding a Mersenne prime is 1 in 26429, which is about 90 times higher than I calculate at 1 in 2400000.

Clearly I'm doing something wrong, as I calculate the chances of finding a prime are about 90x worse than mprime. I realise that the fact that the exponents have been trial factored increases the chances of finding a prime compared to a random number, but does that account for all the difference?

[Uncwilly]

I think that the error is in the mprime display. Your value feels right.

[R. Gerbicz]

Your error is in the assumptation that mp is prime by 1/log(mp) probability. The true is that it is prime by a "much" larger probability: c*log(log(mp))/log(mp) probability [what is constant*log(p)/p], for c see https://primes.utm.edu/notes/faq/NextMersenne.html . Here in the search for Mersenne primes you even have a higher chance to find a prime, because these Mp values survived trial factoring (TF), and most of them some P-1 factoring works.

ps. it doesn't mean that it is a golden treasure in all possible way, because Mp's prime divisors are in the q=2*k*p+1 form, so you basically lost that log(log(mp)) factor in TF.

[drkirkby]

Originally Posted by R. Gerbicz

Your error is in the assumptation that mp is prime by 1/log(mp) probability. .

Yes, I also rather stupidly ignored that for Mersenne Prime testing, one is only testing odd numbers, so the probability of a random odd number being prime is 2/log(mp), not 1/log(mp). But that only changes my original estimate by a factor of 2, which is still a factor of 45 different to mprime.

Originally Posted by R. Gerbicz

The true is that it is prime by a "much" larger probability: c*log(log(mp))/log(mp) probability [what is constant*log(p)/p], for c see https://primes.utm.edu/notes/faq/NextMersenne.html . Here in the search for Mersenne primes you even have a higher chance to find a prime, because these Mp values survived trial factoring (TF), and most of them some P-1 factoring works.

ps. it doesn't mean that it is a golden treasure in all possible way, because Mp's prime divisors are in the q=2*k*p+1 form, so you basically lost that log(log(mp)) factor in TF.

I'm not following all that, but I get that the probability 2^p -1 being prime is higher than predicted by the Prime Number Theorem as various tests have been done to eliminate easy cases - there's also the fact that we only test odd numbers. .

Dave

[R. Gerbicz]

Originally Posted by drkirkby

Yes, I also rather stupidly ignored that for Mersenne Prime testing, one is only testing odd numbers, so the probability of a random odd number being prime is 2/log(mp), not 1/log(mp). But that only changes my original estimate by a factor of 2, which is still a factor of 45 different to mprime.

I've written a quite different formula. At least one PRP/LL test has performed on each prime for p<10^8.
Lets check out what your original formula says about the number of Mersenne primes for p in [2,10^8]:

? s=0;forprime(p=2,10^8,s+=1/(p*log(2)));s
%1 = 4.5805210191516384202989894082293513116

Even if you multiple it by two, you are still very far from the true number.

Setting the exact conjectured formula for p>2:

? prob(p)=exp(0.57721)*log(if(p%4==1,6,2)*p)/(p*log(2))
%2 = (p)->exp(0.57721)*log(if(p%4==1,6,2)*p)/(p*log(2))
? 
? 
? s=0;forprime(p=3,10^8,s+=prob(p));s
%3 = 51.085172813343511293322996421505390399
?

It is quite good, because (as we know M2 is prime) this says that the expected number of Mersenne primes is 52.08, and we know
51 [a double check still can find one or more]. That is impressive. Returning to your problem:

? p=1.1e7;
? 1-(1-prob(p)*log(p)/log(2^76))^32
%5 = 3.8894933977891824505859557809830582239 E-5
? 
? 1/%5
%6 = 25710.289174636665666484139748462123138

This is how Maths works. Not far from the original number, used the same setup p=1.1e8 and 2^76 TF limit, but not
included the P-1 works. Still the probability is close to the program's number.