OFFSET
1,1
COMMENTS
Equivalently, numbers that are equal to the sum of their proper exponential divisors, with one of them taken with a minus sign.
From Amiram Eldar, Dec 05 2025: (Start)
If k is a term, and m is a squarefree number coprime to k, then k*m is also a term. The primitive terms in this sequence (A391283) are the powerful (A001694) terms.
The asymptotic density of this sequence is Sum_{n>=1} f(A391283(n)) = 0.00088534..., where f(n) = (6/(Pi^2*n)) * Product_{prime p|n} (p/(p+1)).
The least odd term is a(5742337) = A336681(1) = 6485886225. (End)
LINKS
Amiram Eldar, Table of n, a(n) for n = 1..10000
EXAMPLE
900 is a term since 900 = 30 + 60 + 90 + 150 - 180 + 300 + 450 is the sum of its proper exponential divisors with one of them, 180, taken with a minus sign.
MATHEMATICA
dQ[n_, m_] := (n > 0 && m > 0 && Divisible[n, m]); expDivQ[n_, d_] := Module[{ft = FactorInteger[n]}, And @@ MapThread[dQ, {ft[[;; , 2]], IntegerExponent[d, ft[[;; , 1]]]}]]; esigma[n_] := Times @@ (Sum[First[#]^d, {d, Divisors[Last[#]]}] &) /@ FactorInteger[n]; expAdmQ[n_] := (ab = esigma[n] - 2*n) > 0 && EvenQ[ab] && ab/2 < n && Divisible[n, ab/2] && expDivQ[n, ab/2]; Select[Range[50000], expAdmQ]
PROG
(PARI) fun(p, e) = sumdiv(e, d, p^d);
isexpdiv(f, d) = {for(i = 1, #f~, if(d % f[i, 1], return(0)); if(f[i, 2] % valuation(d, f[i, 1]), return(0))); 1; }
isok(k) = {my(f = factor(k), s = prod(i = 1, #f~, fun(f[i, 1], f[i, 2])), ab = s - 2*k); ab > 0 && !(ab % 2) && !(k % (ab/2)) && isexpdiv(f, ab/2); } \\ Amiram Eldar, Dec 05 2025
KEYWORD
nonn
AUTHOR
Amiram Eldar, Jul 30 2020
STATUS
approved