%I #34 Dec 23 2022 17:25:37
%S 728,2964,3864,4959,5935,6187,9386,9633,11695,13764,16536,16591,20784,
%T 25428,28808,29623,32696,33632,35805,39585,43736,44733,49027,55344,
%U 56336,57663,58305,62634,65912,65974,66650,67067,67728,69279,69835,73615,73616,74168
%N Smaller of Smith brothers.
%H Amiram Eldar, <a href="/A050219/b050219.txt">Table of n, a(n) for n = 1..10000</a> (terms 1..1000 from Robert Israel)
%H Eric Weisstein's World of Mathematics, <a href="http://mathworld.wolfram.com/SmithBrothers.html">Smith Brothers.</a>
%p issmith:= proc(n)
%p if isprime(n) then return false fi;
%p convert(convert(n,base,10),`+`) = add(t[2]*convert(convert(t[1],base,10),`+`),t=ifactors(n)[2])
%p end proc:
%p S:= select(issmith, {$4..10^5}):
%p sort(convert(S intersect map(`-`,S,1), list)); # _Robert Israel_, Jan 15 2018
%t smithQ[n_] := !PrimeQ[n] && Total[Flatten[IntegerDigits[Table[#[[1]], {#[[2]]}]& /@ FactorInteger[n]]]] == Total[IntegerDigits[n]];
%t Select[Range[10^5], smithQ[#] && smithQ[#+1]&] (* _Jean-François Alcover_, Jun 07 2020 *)
%o (PARI) isone(n) = {if (!isprime(n), f = factor(n); sumdigits(n) == sum(k=1, #f~, f[k,2]*sumdigits(f[k,1])););}
%o isok(n) = isone(n) && isone(n+1); \\ _Michel Marcus_, Jul 17 2015
%o (Python)
%o from sympy import factorint
%o from itertools import count, islice
%o def sd(n): return sum(map(int, str(n)))
%o def smith():
%o for k in count(1):
%o f = factorint(k)
%o if sum(f[p] for p in f) > 1 and sd(k) == sum(sd(p)*f[p] for p in f):
%o yield k
%o def agen():
%o prev = -1
%o for s in smith():
%o if s == prev + 1: yield prev
%o prev = s
%o print(list(islice(agen(), 38))) # _Michael S. Branicky_, Dec 23 2022
%Y Cf. A006753, A050220.
%K nonn,base
%O 1,1
%A _Eric W. Weisstein_
%E Offset corrected by _Arkadiusz Wesolowski_, May 08 2012
|