# Efficient algorithm for finding Sophie Germain primes

What's the industry standard for an efficient finding large Sophie Germain primes?

As a part of request handling in my application, I need to generate Paillier key.

My current approach is to generate a pseudorandom probably prime number q backed up by Miller-Rabin test, then setting p = 2q + 1 and checking if p is also prime with Miller-Rabin and Baillie-PSW. This is repeated until p is prime.

loop {
q <- randomPrime(nbits - 1) // generates probably prime number,
// performs Miller-Rabin tests on the result
p = 2q + 1
return if isPrime(p)  //applying the Miller-Rabin & Baillie-PSW tests
}


For nbits = 1024 one loop step takes about 100ms on my laptop. Not that bad, but in order to find desired p and q I need to execute several loop steps. The best result so far was 156 retries (~16s), the worst one was 1616 retries (~160s).

It's a way above required latency so I wonder what's considered the correct approach nowadays.

• See math.stackexchange.com/questions/870626/… for a short cut which uses only a base-2 Fermat test for the second number. And there is a similar question by @fgrieu at mathoverflow.net/questions/295516/… Jun 26, 2018 at 16:01
• So it's enough to check 2^{p-1} ≡ 1 (mod n) instead of running expensive isPrime(p) ? Jun 26, 2018 at 21:20
• Yes, if you know that $(p-1)/2$ is prime, then $2^{p-1} \equiv 1 \pmod p$ implies that $p$ is prime as well. On the other hand, that observation doesn't save you that much time (essentially, a few MR iterations at the very end when you do find a prime $p$) Jun 26, 2018 at 21:57

I don't know of any industry standards, however it's obvious that your code is fairly suboptimal.

The initial randomPrime() will effectively pick a candidate number $q$, and then spend a lot of time ensuring that number is prime. Now, half the time it'll pick a $q \equiv 1 \pmod 3$; if so, p will be a multiple of 3, and so you've just wasted all that time testing if $q$ is prime.

Now, what all random prime routines I've seen do is first ensure that their candidate number doesn't have any small factors (so we don't need to run our expensive primality tests on those values); this can be done either by computing $q \bmod s$ for small primes $s$, or by creating a sieve (a bitmap where each bit $i$ corresponds a value $a + bi$; the bit gets marked if $a + bi$ has a small factor).

Both methods can be extended to check either $q$ or $p = 2q+1$ has a small factor:

• For the explicit $q \bmod s$ computation, we skip it if that is either 0 or $(s-1)/2$

• For the sieve case, we mark bits where $a + bi \bmod s$ is either 0 or $(s-1)/2$

Doing this will significantly reduce the number of iterations you'll need to perform.

Now, if your randomPrime() logic is in a library you can't touch, this advice doesn't help you. However, if you can do it, this will be of help.

• That's a very good suggestion. I dig into randomPrime and isPrime implementation and found it's already doing an optimization you mentioned. That is, checks mod for 15 small primes before running expensive tests. I wonder if it's possible to select q in some special way to increase the possibility it's a Sophie Germain prime but did not find anything about it. Jun 26, 2018 at 21:17
• Two things: 1) you could do better than checking the first 15 small primes. If performing a $q \bmod s$ for single precision $s$ is 1000 times faster than a Miller-Rabin iteration (and it should be considerably faster than that), then it's a win to test all primes $< 1000$. Doing sieving is an even bigger when (as it makes sense to remove more small primes), but it's more work to implement 2) if you could modify the randomPrime code, you could insert the test I suggested directly (which would ensure that $p$ doesn't have tiny factors either) Jun 26, 2018 at 21:50
• These comments should move into the answer. Also, see this answer Jan 24, 2020 at 17:10