login

Year-end appeal: Please make a donation to the OEIS Foundation to support ongoing development and maintenance of the OEIS. We are now in our 61st year, we have over 378,000 sequences, and we’ve reached 11,000 citations (which often say “discovered thanks to the OEIS”).

Prime numbers that have the same base-10 digits as their prime index (with multiplicity), disregarding zero digits.
3

%I #48 Jul 07 2022 08:03:42

%S 5701,27091,50417,55049,57089,60601,61051,63607,66301,72019,132607,

%T 270913,284057,574031,654301,936007,936907,950647,1250609,1461001,

%U 1504417,1580921,1682069,1703287,1750631,1810553,1810573,1837601,2050241,2145089,2485001,2641109,2709169,2800333,2805703

%N Prime numbers that have the same base-10 digits as their prime index (with multiplicity), disregarding zero digits.

%e 15097067 is a term since its index primepi(15097067) = 976571 shares the same multiset of nonzero digits {1,5,6,7,7,9}.

%t a[max_] := Module[{l}, Select[{#, Prime[#]} & /@ Range[max], (l = IntegerDigits[#[[2]]]; SortBy[Tally[l], First] === SortBy[Tally[PadLeft[IntegerDigits[#[[1]]], Length[l]]], First]) &]]; a[10^6][[All, 2]] (* Gives the first 108 terms *)

%o (Python)

%o from sympy import nextprime

%o from itertools import islice

%o def b10s(n): return "".join(sorted(str(n))).lstrip("0")

%o def agen():

%o k, pk = 1, 2

%o while True:

%o if b10s(k) == b10s(pk): yield pk

%o k, pk = k+1, nextprime(pk)

%o print(list(islice(agen(), 32))) # _Michael S. Branicky_, Jun 28 2022

%o (PARI) digs(k) = vecsort(select(x->(x>0), digits(k)));

%o isok(p) = if (isprime(p), digs(p) == digs(primepi(p))); \\ _Michel Marcus_, Jul 05 2022

%Y Cf. A355318 (the prime indices), A355539.

%Y Cf. A355418.

%K base,easy,nonn

%O 1,1

%A _Xiaofeng Wang_, Jun 28 2022