Why does GMP only run Miller-Rabin test twice when generating a prime?

In mpz_nextprime(), after some sieving with small primes, an MR test function is called, with the number of trials set to 25 (https://github.com/alisw/GMP/blob/master/mpz/nextprime.c#L118):

      if (mpz_millerrabin (p, 25))
        goto done;

But then in mpz_millerrabin(), for large enough candidates (candidates bigger than $31\times2^{46}$), the number of trials is suddenly reduced to 1 (https://github.com/alisw/GMP/blob/master/mpz/millerrabin.c#L142):

      reps -= 24;

with no explanations.

As a result, for large candidates only three tests are run:

• MR with base 2
• Lucas
• MR with random base

But the comment at the top of the millerrabin.c says that "Knuth indicates that 25 passes are reasonable." What is the reason the number of trials is unconditionally slashed this way?

The article, "Strengthening the Baillie-PSW primality test" referred to above, suggests adding a third test to the standard BPSW (MR test base 2 combined with Lucas).

But there's a third congruence that is true for primes, and rarely true for composite $N$: $V_{N+1} = 2 Q \pmod N$. A composite that satisfies this congruence is called a $V$-pseudoprime.

Up to 1E15, there are about 2 million base-2 psp's, and about 2 million Lucas pseudoprimes (and no overlap between these two sets).

However, there are only five $V$-pseudoprimes up to 1E15. None of these is either a base-2 psp or a Lucas pseudoprime.

Therefore, we recommended that a primality test include a check on the value of $V_{N+1}$. The $V$'s are typically computed along with the $U$'s, so almost no additional computation would be required.

Further thoughts: Doing dozens of MR tests (that is, strong Fermat tests) tests may not be the way to go.

1. It takes roughly 4 times as long to do a BPSW test (Fermat + Lucas) as to do a Fermat test. It is not clear that the additional time needed to do one or two dozen Fermat tests produces more trustworthy results than BPSW.

2. Most importantly: Fermat tests to different bases are not independent. The Pomerance/Selfridge/Wagstaff paper The Pseudoprimes to 25⋅10⁹, has data to back this up. For example, Table 6 shows that there are 21853 psp(2) below 25*10^9, but that 4709 of these (21%) are also psp(3).

   The Baillie/Wagstaff paper Lucas Pseudoprimes, explains why: If N is a pseudoprime to some base a, this is probably because N is one of those few numbers that is pseudoprime to many bases. Therefore, N is more likely than most numbers of that size to be pseudoprime to another base, b. Theorem 1 in that paper gives a formula for the number of such bases mod N. For example, N = 341 is psp(2), but it is also psp to 99 other bases between 2 and N - 2.

By contrast, there are hints that Fermat and Lucas tests are independent. For example, Fermat and Lucas psp's tend (with exceptions, of course) to fall into residue classes +1 and -1 (mod m) for small m.