OFFSET
1,1
COMMENTS
If m is in the sequence and d|m, then m^d is also a term. Note that this sequence contains all infinite subsequences of the form p^(p^k) for k>0, where p is a prime. - Amiram Eldar and Thomas Ordowski, Feb 06 2019
If one selects some composite k, k >= 8, and decomposes (k - sopfr(k)) into an additive partition having only prime parts, then those parts, when taken as a product with k, yield an element of this sequence. - Christopher Hohl, Jul 30 2019
LINKS
François Huppé, Table of n, a(n) for n = 1..50000 (terms 1..1000 from T. D. Noe)
K. Alladi and P. Erdős, On an additive arithmetic function, Pacific J. Math., Volume 71, Number 2 (1977), 275-294. See "special numbers" on page 287.
EXAMPLE
a(38) = 884 = 2 * 2 * 13 * 17 -> 2 + 2 + 13 + 17 = 34 so 884 / 34 = 26.
MAPLE
isA046346 := proc(n)
if isprime(n) then
false;
elif modp(n, A001414(n)) = 0 then
true;
else
false;
end if;
end proc:
for n from 2 to 1000 do
if isA046346(n) then
printf("%d, ", n);
end if;
end do: # R. J. Mathar, Jan 12 2016
MATHEMATICA
Select[Range[2, 1170], !PrimeQ[#]&&IntegerQ[#/Total[Times@@@FactorInteger[#]]]&] (* Jayanta Basu, Jun 02 2013 *)
PROG
(PARI) sopfr(n) = {my(f=factor(n)); sum(k=1, #f~, f[k, 1]*f[k, 2]); }
lista(nn) = forcomposite(n=2, nn, if (! (n % sopfr(n)), print1(n, ", ")); ); \\ Michel Marcus, Jan 06 2016
(MATLAB) m=1; for u=2:1200 if and(isprime(u)==0, mod(u, sum(factor(u)))==0); sol(m)=u; m=m+1; end; end; sol % Marius A. Burtea, Jul 31 2019
(Magma) [k:k in [2..1200]| not IsPrime(k) and k mod (&+[m[1]*m[2]: m in Factorization(k)]) eq 0]; // Marius A. Burtea, Jul 31 2019
(Python)
from sympy import factorint
def ok(n):
f = factorint(n)
return sum(f[p] for p in f) > 1 and n % sum(p*f[p] for p in f) == 0
print(list(filter(ok, range(1250)))) # Michael S. Branicky, Apr 16 2021
CROSSREFS
KEYWORD
nonn
AUTHOR
Patrick De Geest, Jun 15 1998
EXTENSIONS
Description corrected by Robert A. Stump (bee_ess107(AT)yahoo.com), Jan 09 2002
STATUS
approved