OFFSET
1,1
COMMENTS
From Bernard Schott, Oct 20 2021: (Start)
Except for 11, all terms have an odd number of digits.
Except for terms of the form 10^k+1, k >= 2, the middle digit is always odd; the unique known term of the form 10^k+1 for 2 <= k <= 100000 is 101 (see comment in A000533). (End)
REFERENCES
Charles W. Trigg, Journal of Recreational Mathematics, Vol. 20(2), 1988.
LINKS
Michel Marcus, Table of n, a(n) for n = 1..1215
Mike Mudge, Morph code, Hands On Numbers Count, Personal Computer World, May 1997, p. 290.
EXAMPLE
373 -> 3^2 + 7^2 + 3^2 = 67, which is prime.
MATHEMATICA
Select[Prime@ Range[2, 10^4], And[PalindromeQ@ #, PrimeQ@ Total[IntegerDigits[#]^2]] &] (* Michael De Vlieger, Oct 20 2021 *)
PROG
(PARI) isok(p) = my(d=digits(p)); isprime(p) && (d==Vecrev(d)) && isprime(sum(k=1, #d, d[k]^2)); \\ Michel Marcus, Oct 17 2021
(Python)
from sympy import isprime
def ok(n):
s = str(n)
return s==s[::-1] and isprime(n) and isprime(sum(int(d)**2 for d in s))
print([k for k in range(10**5) if ok(k)]) # Michael S. Branicky, Nov 23 2021
(Python) # second version for going to large terms
from sympy import isprime
from itertools import product
def ok(pal):
return isprime(pal) and isprime(sum(int(d)**2 for d in str(pal)))
def agentod(maxdigs):
yield 11
for d in range(3, maxdigs+1, 2):
pal = 10**(d-1) + 1
if ok(pal): yield pal
for first in "1379":
for left in product("0123456789", repeat=(d-3)//2):
left = "".join(left)
for mid in "13579":
pal = int(first + left + mid + left[::-1] + first)
if ok(pal): yield pal
print([an for an in agentod(5)]) # Michael S. Branicky, Nov 23 2021
CROSSREFS
KEYWORD
nonn,base
AUTHOR
Patrick De Geest, Dec 15 1999
STATUS
approved