The Question

A Sophie Germain prime is a prime \$p\$ such that \$2p+1\$ is prime as well. For example, 11 is a Sophie Germain prime because 23 is prime as well. Write the shortest program to calculate Sophie Germain primes in ascending order


  • The Sophie Germain primes must be generated by your program, not from an external source.
  • Your program must calculate all Sophie Germain primes under \$2^{32}-1\$
  • You must print each distinct Sophie Germain prime your program finds.
  • The person with the lowest score wins


  • 2 points per byte of your code
  • -10 if you can show a prime generated by your program greater than \$2^{32}-1\$
  • \$\begingroup\$ Comments are not for extended discussion; this conversation has been moved to chat. \$\endgroup\$ – Martin Ender May 11 '15 at 23:41

16 Answers 16


Pyth, 19 bytes * 2 - 10 = 28

Note that the online compiler/executor doesn't show output because it's an infinite loop.



K1                      K=1
  #                     While true:
   ~K1                  K+=1
      I                 If
       &                logical AND
        !tPK            K is prime
            !tPhyK      2*K+1 is prime (y is double, h is +1)
                  K     Print K
  • \$\begingroup\$ PZ doesn't return a truthy or falsy value. It returns the prime factorization of Z. Testing for prime is !tPZ, which checks if the prime factorization only contains one factor. \$\endgroup\$ – Jakube May 11 '15 at 21:39
  • \$\begingroup\$ Yes. Now it works. !tP mistakes 0 and 1 to be prime though, since their prime factorization only contains 1 factor. Easy fix is to replace all Z by K and assign K2 at the beginning. \$\endgroup\$ – Jakube May 11 '15 at 21:45
  • \$\begingroup\$ Some other golfs: assign K1 instead of K2, and swap the if and the increment. This way you can remove the ). And +1*K2 is the same thing as hyK. \$\endgroup\$ – Jakube May 11 '15 at 21:49
  • \$\begingroup\$ Ah, I had just read about those on the tutorial page. Does it work for you on pyth.herokuapp.com/?code=K2%23I%26!tPK!tPhyKK)~K1&debug=0 \$\endgroup\$ – mbomb007 May 11 '15 at 21:51
  • \$\begingroup\$ The online compiler doesn't show a result, because the program is stuck in an infinite loop. And the website shows only output, after the program finishes. I've tested the code using the offline compiler. It works. \$\endgroup\$ – Jakube May 11 '15 at 21:57


For 17 chars we get full enumeration up to 2^32:


For 4 chars more, we get a range just large enough to include an SG prime greater than 2^32:


since 4294967681 = 2^32 + 385 < 2^32 + 400.

Of course, we could equally extend the range for free as

  • \$\begingroup\$ This means you can submit it without the bonus for 17 characters or with the bonus for 21 characters \$\endgroup\$ – Meow Mix May 11 '15 at 21:02
  • \$\begingroup\$ @user3502615, or with the bonus for 17 characters. Although it's debatable whether the SG prime I list was actually generated "by my program", since I don't have a powerful enough computer to run it that far. \$\endgroup\$ – Peter Taylor May 11 '15 at 21:04
  • \$\begingroup\$ I, treats I as a signed 32-bit integer, so the maximum value for I is 2 ** 31 - 1. \$\endgroup\$ – Dennis May 11 '15 at 21:07
  • 2
    \$\begingroup\$ @Dennis, is that a documented property of the language or an implementation quirk of the Java implementation? \$\endgroup\$ – Peter Taylor May 11 '15 at 21:11
  • \$\begingroup\$ It's not documented, but the behavior is consistent for both the Java and the online interpreter. \$\endgroup\$ – Dennis May 11 '15 at 21:13

Husk, 16 8 7 * 2 - 10 = 4

-8 12 thanks to @Jo King


Try it online!

A list of all natural numbers that are Sophie Germain primes. I hope it's valid, because it doesn't terminate when you try to show the entire list, but Husk lazily evaluates it, so you can just take the required amount.

f         Filter
     İp   the prime numbers
    D     that, when doubled
   →      and incremented by 1
 ȯṗ       are prime

Another way to do it, also 7 bytes:

  • 1
    \$\begingroup\$ I had this exact same program, but it didn't display anything(like you have said), so I wrote an APL answer. \$\endgroup\$ – Razetime Oct 11 '20 at 4:00

05AB1E, 7 bytes, score = 4

Thanks to @ovs for -1 byte!


(Technically) Prints an infinite list of Sophie Germain primes! (Hence, score is subtracted by 10)

Try it online!


∞ʒ        # For each natural number 'p'...
  x>‚     # Pair p with 2p + 1 (The ‚ is not a comma, it is a weird sort of quotation)
     p    # Are they both prime?
      P   # Then multiply the boolean values! (And implicitly print if the product is 1)
  • \$\begingroup\$ ∞ʒx>‚pP for 7 bytes. \$\endgroup\$ – ovs Oct 30 '20 at 9:14
  • \$\begingroup\$ ∞<Øʒ·>p for an alternative 7-byter. Too bad 2 is truthy, otherwise the < could have been removed. \$\endgroup\$ – Kevin Cruijssen Oct 30 '20 at 10:04
  • \$\begingroup\$ This is ‚ SINGLE LOW-9 QUOTATION MARK, not , COMMA. This can be really confusing, but the first pairs up two numbers while the latter prints to STDOUT \$\endgroup\$ – ovs Oct 30 '20 at 13:18
  • \$\begingroup\$ @ovs Just figured that out... :P \$\endgroup\$ – SunnyMoon Oct 30 '20 at 13:18

Pyth - 2 * 16 bytes - 10 = 22

Uses the customary method of prime checking in Pyth with the !tP and applies it both to the number and its safe-prime, with a little trick to check both at once. Goes up to 10^10, so I'm going for the bonus.


Explanation coming soon.

f          r2^TT     Filter from 2 till 10^10
 !                   Logical not to detect empty lists
  +                  List concatenation
   tP                All but the firs element of the prime factorization
    T                The filter element
   tP                All but the firs element of the prime factorization
    hyT              2n+1

Try under 1000 online.

  • 1
    \$\begingroup\$ This takes a machine with about 40 GB of RAM memory. Quite efficient ;-) \$\endgroup\$ – Jakube May 11 '15 at 21:26
  • \$\begingroup\$ I don't think you can claim the - 10 unless you have actually successfully run the code? \$\endgroup\$ – orlp May 11 '15 at 23:53
  • \$\begingroup\$ @orlp no, I asked OP and he said making the range smaller and simulating the whole program would be sufficient: chat.stackexchange.com/transcript/message/21585393#21585393 \$\endgroup\$ – Maltysen May 12 '15 at 0:17

APL (Dyalog Extended), 10 - ⍨ 2 × 27 = 44


Try it online!

Made with some golfing help and general APL help from dzaima.


                        1+⊢ add 1 to the input (initially, 1)
             1+2×⍵          double and increment it
         1⍭⍵,               prepend the primality check(1 or 0)
       ∧\                   and scan with LCM
   ⍵⍴⍨⌽                     reverse & reshape using self
                            (if test is falsy, then this results in ⍬)
 ⎕←                         and display it with newline

CJam, 34 (2 * 22 - 10)


Prints all Sophie Germain primes under 12 ** 9, which includes 4294967681 > 2 ** 32.

I estimate that this will take roughly 8 hours on my machine. I'll run it tonight.


Haskell, 2*54-10 = 98 132

i a=all((>0).rem a)[2..a-1]
p=[n|n<-[2..],i n,i$2*n+1]

i is a prime check. p takes all numbers n where both n and 2*x+1 are prime. p is an infinite list.

Edit: better way for checking if 2*n+1 is prime.


Julia, 2*49 - 10 = 88


Prints them in list format, [2,3,5,11,...]. If that format, using the primes function, or waiting until all the computation is done to print isn't acceptable, this prints them one per line as it runs.

for i=1:2^33;f(i)&&f(2i+1)&&println(i)end

It's a little longer, 52 chars. Both compute all the Sophie Germain primes up to 2^33, so they should get the 10 point discount.


Python 3, 124 123 bytes

while 1:
 for x in range(2,round(i**.5)+1):p=min(p,i%x)
 if p:
  if s in q:print(s)

How does it work?

i=3                                 # Start at 3
q=[2]                               # Create list with first prime (2), to be list of primes.
while 1:                            # Loop forever
 p=1                                # Set p to 1 (true)
 for x in range(2,round(i**0.5)+1): # Loop from 2 to the number's square root. x is the loop value
     p=min(p,i%x)                   # Set p to the min of itself and the modulo of
                                    # the number being tested and loop value (x).
                                    # If p is 0 at the end, a modulo was 0, so it isn't prime.
 if p:                              # Check if p is 0
  q+=[i]                            # Add the current number (we know is prime) to list of primes (q)
  s=(i-1)/2                         # Generate s, the number that you would double and add 1 to make a prime.

  if s in q:print(s)                # If (i-1)/2 is a prime (in the list), then that prime satifies
                                    # the condition 2p+1 is prime because i is 2p+1, and i is prime
 i+=2                               # Increment by 2 (no even numbers are prime, except 2)

Try it online here.

My computer says it is has generated 0.023283 % of all the Sophie Germain primes below 2^32.

When it's finished, I'll post it on pastebin if there are enough lines. You can use it to check you've got them all.

  • 1
    \$\begingroup\$ .5 is shorter than 0.5 \$\endgroup\$ – mbomb007 May 11 '15 at 21:30

Perl, 2*57-10 = 104

use ntheory":all";forprimes{say if is_prime(2*$_+1)}2**33


42 seconds to 2^32, 1m26s to 2^33. Will run 50% faster if 2*$_+1 is written as 1+$_<<1 but that's one more byte.

The module also installs primes.pl which has lots of filters including one for Sophie-Germain primes. So: primes.pl --so 2**33 (20 bytes)


Ruby, 61*2 - 10 = 112

require'prime';Prime.each(1.0/0)do|n|p Prime.prime?(n*2+1)end

It would take forever to print out all values up to 2**32


Shaved off a few bytes substituting Float::INFINITY for 1.0/0


PARI/GP, 46 * 2 - 10 = 82


Python 3, 2 * 64 bytes - 10 = 118

while 1:n+=1;all((n-~n)*n%k for k in range(2,n))>0!=print(n)

Try it online!

It turns out a normal prime test is shorter.

Python 3, 2 * 67 bytes - 10 = 124

while 1:P*=k*k;k+=1;P%k*k//2in d!=print(k//2);d+=P%k*[k]

Try it online!

This uses Wilson's theorem to generate a list of primes. When k is prime and k//2 is prime as well (in d), this prints k//2.

This can be done a little more memory efficient by not collecting a list of primes but rather doing two prime tests at once at the cost of 4 bytes:

while 1:P*=k*k;k+=1;R*=m*m*-~m*-~m;m+=2;P%k*R%m>0!=print(k)

Try it online!


Jelly, 10 bytes => score 10


Generates primes up to 10^10 > 2^32.

Try it online! (TIO link is to a shorter version that generates only up to 1000, so it doesn't time out)


Stax, 14 bytes


Run and debug it


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.