OFFSET
1,1
LINKS
Amiram Eldar, Table of n, a(n) for n = 1..10000 (terms 1..1000 from Robert Israel)
Eric Weisstein's World of Mathematics, Smith Brothers.
MAPLE
issmith:= proc(n)
if isprime(n) then return false fi;
convert(convert(n, base, 10), `+`) = add(t[2]*convert(convert(t[1], base, 10), `+`), t=ifactors(n)[2])
end proc:
S:= select(issmith, {$4..10^5}):
sort(convert(S intersect map(`-`, S, 1), list)); # Robert Israel, Jan 15 2018
MATHEMATICA
smithQ[n_] := !PrimeQ[n] && Total[Flatten[IntegerDigits[Table[#[[1]], {#[[2]]}]& /@ FactorInteger[n]]]] == Total[IntegerDigits[n]];
Select[Range[10^5], smithQ[#] && smithQ[#+1]&] (* Jean-François Alcover, Jun 07 2020 *)
PROG
(PARI) isone(n) = {if (!isprime(n), f = factor(n); sumdigits(n) == sum(k=1, #f~, f[k, 2]*sumdigits(f[k, 1])); ); }
isok(n) = isone(n) && isone(n+1); \\ Michel Marcus, Jul 17 2015
(Python)
from sympy import factorint
from itertools import count, islice
def sd(n): return sum(map(int, str(n)))
def smith():
for k in count(1):
f = factorint(k)
if sum(f[p] for p in f) > 1 and sd(k) == sum(sd(p)*f[p] for p in f):
yield k
def agen():
prev = -1
for s in smith():
if s == prev + 1: yield prev
prev = s
print(list(islice(agen(), 38))) # Michael S. Branicky, Dec 23 2022
CROSSREFS
KEYWORD
nonn,base
AUTHOR
EXTENSIONS
Offset corrected by Arkadiusz Wesolowski, May 08 2012
STATUS
approved