# Primes numbers with prime index

Write a program or function that outputs/returns the first 10000 prime-indexed prime numbers.

If we call the nth prime p(n), this list is

3, 5, 11, 17, 31, 41, 59 ... 1366661


because

p(p(1)) = p(2) = 3
p(p(2)) = p(3) = 5
p(p(3)) = p(5) = 11
p(p(4)) = p(7) = 17
...
p(p(10000)) = p(104729) = 1366661


Standard loopholes are forbidden, and standard output methods are allowed. You may answer with a full program, a named function, or an anonymous function.

• You should generally try to post challenges in the sandbox first (see the link on the right side) in order to work out the issues. – aditsu quit because SE is EVIL Oct 23 '15 at 17:43
• Optimizing for runtime is not what we do in a code-golf challenge; the shortest program always wins. – lirtosiast Oct 23 '15 at 18:53
• Primes with prime subscripts: A006450. – user12166 Oct 23 '15 at 21:55
• @bilbo Answers for code golf are usually accepted after a week, and should be accepted as the shortest successful code. If you wanted code speed, there is a tag for that. See this page about the tag code-golf. – Addison Crump Oct 24 '15 at 12:04
• All contests need an objective winning criterion; they are off topic otherwise. If you're going to judge answers by size and speed, you need to disclose a way to combine both. This should be done when the contest is posted, not 14 hours and 10 answers later. I have undone all speed-related edits, since the only other option would be to close this post for being off topic. – Dennis Oct 25 '15 at 14:10

# MATLAB/Octave, 25 bytes

p=primes(2e6)
p(p(1:1e4))


It doesn't get much more straightforward than this.

## Python, 72 bytes

P=p=1;l=[]
while p<82e5/6:l+=P%p*[p];P*=p*p;p+=1
for x in l:print l[x-1]


This terminates with a "list index out of range error" after printing the 10000 numbers, which is allowed by default.

Uses the Wilson's Theorem method to generate a list l of the primes up to the 10000th prime. Then, prints the primes with the positions in the list, shifted by 1 for zero-indexing, until we run out of bounds after the 10000th prime-th prime.

Conveniently, the upper bound of 1366661 can be estimated as 82e5/6 which is 1366666.6666666667, saving a char.

I'd like a single-loop method, printing prime-indexed primes as we add them, but it seems to be longer.

P=p=1;l=[]
while p<104730:
l+=P%p*[p]
if len(l)in P%p*l:print p
P*=p*p;p+=1

• This is way better than the garbage I was writing. +1 – user45941 Oct 23 '15 at 17:51
• This only prints 1229 numbers – aditsu quit because SE is EVIL Oct 23 '15 at 19:03
• @aditsu I think I see my mistake. Are you able to run this code with the bigger bound? – xnor Oct 23 '15 at 19:07
• It will probably take a long time :p – aditsu quit because SE is EVIL Oct 23 '15 at 19:16
• I think it finished ＼（＠；◇；＠）／ , it seems correct – aditsu quit because SE is EVIL Oct 24 '15 at 15:59

## J, 11 bytes

p:<:p:i.1e4


Outputs the primes in the format

3 5 11 17 31 41 59 67 83 109 127 ...


## Explanation

        1e4  Fancy name for 10000
i.     Integers from 0 to 9999
p:       Index into primes: this gives 2 3 5 7 11 ...
<:         Decrement each prime (J arrays are 0-based)
p:           Index into primes again


# Mathematica, 2625 23 bytes

Prime@Prime@Range@1*^4&


Pure function returning the list.

• Prime is Listable so a simple Prime@Prime@Range@1*^4& will do – user46060 Oct 23 '15 at 22:34
• I know the feeling ... In any case, I think this is the prettiest Mathematica solution I have seen on here! – user46060 Oct 23 '15 at 22:36
• Let me guess, the @ operator has higher precedence than ^ when writing Range@10^4? That's classic Mathematica messing up your game of golf. Good trick! – user46060 Oct 23 '15 at 23:14

p=[x|x<-[2..],all((>0).mod x)[2..x-1]]
f=take 10000$map((0:p)!!)p  Outputs: [3,5,11,17,31,41,59,67,83,109,127.....<five hours later>...,1366661] Not very fast. How it works: p is the infinite list of primes (naively checking all mod x ys for y in [2..x-1]) . Take the first 10000 elements of the list you get when 0:p!! (get nth element of p) is mapped on p. I have to adjust the list of primes where I take the elements from by prepending one number (-> 0:), because the index function (!!) is zero based. # AWK - 129 bytes ...oookay... too long to win points for compactness... but maybe it can gain some honor for the speed? The x file: BEGIN{n=2;i=0;while(n<1366662){if(n in L){p=L[n];del L[n]}else{P[p=n]=++i;if(i in P)print n}j=n+p;while(j in L)j=j+p;L[j]=p;n++}}  Running: $ awk -f x | nl | tail
9991  1365913
9992  1365983
9993  1366019
9994  1366187
9995  1366327
9996  1366433
9997  1366483
9998  1366531
9999  1366609
10000  1366661


BEGIN {
n=2
i=0
while( n<1366662 ) {
if( n in L ) {
p=L[n]
del L[n]
} else {
P[p=n]=++i
if( i in P ) print n
}
j=n+p
while( j in L ) j=j+p
L[j]=p
n++
}
}


The program computes a stream of primes using L as "tape of numbers" holding found primes jumping around on L to flag the nearby numbers already known to have a divisor. These jumping primes will advance while the "tape of numbers" L is chopped off number by number from its beginning.

While chopping off the tape head L[n] being empty means there is no known (prime) divisor.

L[n] holding a value means, this value is a prime and known to divide n.

So either we have found a prime divisor or a new prime. Then ths prime will be advanced to the next L[n+m*p] on the tape found being empty.

This is like the Sieve of Eratosthenes "pulled thru a Klein's bottle". You always act on the tape start. Instead of firing multiples of primes thru the tape, you use the primes being found already as cursors jumping away from the tape start by multiple distances of their own value until a free position is found.

While the outer loop generates one prime or not prime decission per loop, the primes found get counted and stored in P as key, the value of this (key,value) pair is not relevant for the program flow.

If their key i happens to be in P already (i in P), we have a prime of the p(p(i)) breed.

Running:

$time awk -f x.awk | wc -l 10000 real 0m3.675s user 0m3.612s sys 0m0.052s  Take into account that this code does not use external precalculated prime tables. Time taken on my good old Thinkpad T60, so I think it deserves to be called fast. Tested with mawk and gawk on Debian8/AMD64 • good 129 bytes in gawk: now with Debian10/AMD64 on my corei7-i870@3.6Ghz: real 0m2,417s user 0m2,205s sys 0m0,042s – JeanClaudeDaudin Jan 25 '20 at 14:42 • you can save one byte with: BEGIN{n=2;i=0;while(n<1366662){if(n in L){p=L[n];del L[n]}else{P[p=n]=++i;if(i in P)print n}j=n+p;while(j in L)j+=p;L[j]=p;n++}} – JeanClaudeDaudin Jan 25 '20 at 15:08 # PARI/GP, 25 bytes apply(prime,primes(10^4))  # CJam, 19 3D#{mp},_1e4<:(\f=p  You can try it online, but you'll need a little patience :p For the record, the last number is 1366661. # Husk, 11 10 bytes (or 6 5 bytes to list the first 10,000 prime-indexed prime numbers, followed by [the infinite number of] all the rest of them) Thanks to user and Razetime for suggesting the 6 5-byte version! ↑!4İ⁰Ṡm!İp  Try it online! (unfortunately times-out on TIO before printing all 10,000 prime-indexed primes, but at least outputs what it's found so far...) ↑!4İ⁰Ṡm!İp ↑ # print the first !4 # 4th index of İ⁰ # powers of 10 (so, the first 10,000), of Ṡ # hook: apply function to its own argument m! # map index function to İp # prime numbers # (implied by hook) to prime numbers  • Nice answer. The 10000 primes restriction is a shame - too many bytes are wasted on that. – user Nov 18 '20 at 15:22 • If you output an infinite list, it will still have the first 10000 items ;) – Razetime Nov 18 '20 at 15:27 • @user & Razetime - thanks for the suggestion from both of you - I've incorporated it (& credited you both)! – Dominic van Essen Nov 18 '20 at 15:45 • It didn't actually occur to me to print an infinite list - that's all Razetime, but thanks anyway. – user Nov 18 '20 at 16:51 # Perl, 55 bytes use ntheory':all';forprimes{print nth_prime$_,\$/}104729


Uses @DanaJ's Math::Prime::Util module for perl (loaded with the pragma ntheory). Get it with:

cpan install Math::Prime::Util
cpan install Math::Prime::Util::GMP


# Japt, 11 bytes

Èj ©°Tj}jL²


Test it

# Python 3, 76 109 bytes - disqualified - requalified?

edit: I am sorry, I was going for primes but missed the putting, didn´t reach the goal of "primes with primeindex". Your are free to grab and enhance. I'll be glad to be learning from you.

Based on Wilson Theorem

def p(x):
c=n=2
for i in range(3,x):
c*=i-1
if c%i!=0:
print(n,i)
n+=1


## without index, 62 bytes

def p(x):
c=2
for i in range(3,x):
c*=i-1
if c%i!=0:print(i)


Anyone can give me an example or good link, how to make a lambda function from it?

## edit: just to have it done right: 109bytes

def p(x):
c=n=2
p=[2]
for i in range(2,x):
c*=i-1
if c%i!=0:
p.append(i)
if n in p: print(p[-1])
n+=1

• I think you've misunderstood the question. Your functions seem to output all prime numbers (and the first one precedes them with an index). The question defines 'prime-indexed prime numbers' as those in the list of all primes that have an index value that is also a prime. For instance, 7 - the 4th prime - isn't prime-indexed because 4 isn't a prime... – Dominic van Essen Nov 21 '20 at 16:25
• I am sorry. You are right, I misunderstood it :( Thank you very much for spending your time on clearing this up for me! – Bambeleme Nov 22 '20 at 0:59

# VyxalH, 8 7 6 5 bytes

Saved 1 + 1 bytes thanks to Aaron Miller

²ʀǎ‹ǎ


Try it Online!

H flag - Push 100 to the stack
²ʀǎƛ‹ǎ
²        Squared (10,000)
ʀ       Range [0, 10,000)
ǎ      Map each a in that range to the ath prime
‹     Decrement each, because the challenge uses 1-indexing
ǎ    Get the nth prime for each of those

• 7 bytes – Aaron Miller May 3 at 2:26
• @AaronMiller Thanks, forgot about auto-closing! – user May 3 at 2:27
• I completely forgot about the H tag, which presets the stack to 100, allowing for 5 bytes – Aaron Miller May 3 at 15:31

# 05AB1E, 7 bytes (non-competing)

Code:

4°L<Ø<Ø


Try it online!, note that I have changed the 4 into a 2. If you have a lot of time, you can change the 2 back to 4, but this will take a lot of time. I need to fasten the algorithm for this.

Explanation:

4°       # Push 10000 (10 ^ 4)
L      # Create the list [1 ... 10000]
<     # Decrement on every element, [0 ... 9999]
Ø    # Compute the nth prime
<   # Decrement on every element
Ø  # Compute the nth prime


# Jelly, 6 bytes

ȷ4RÆN⁺


Try it online!

Obviously, this times out on TIO. Here is a version that prints as it goes, which reaches $$\57493\$$ (the $$\765\$$th term) before timing out

## How it works

ȷ4RÆN⁺ - Main link. Takes no arguments
ȷ4     - 10000
R    - [1, 2, 3, ..., 9999, 10000]
⁺ - Do twice:
ÆN  -   n'th prime of each


# 05AB1E, 6 bytes

тnp<Ø


Try it online!

     Ø  # prime numbers with...
<   # 1-based...
Ø  # indices in...
Åp    # first...
т       # 100...
n      # ^ 2...
Åp    # primes
# implicit output