|
%I #58 Feb 24 2024 01:03:35
%S 1,4,16,27,64,108,256,432,729,1024,1728,2916,3125,4096,6912,11664,
%T 12500,16384,19683,27648,46656,50000,65536,78732,84375,110592,186624,
%U 200000,262144,314928,337500,442368,531441,746496,800000,823543
%N Numbers k such that Sum_i ( e(i)/p(i) ) is an integer, where the prime factorization of n is Product_i ( p(i)^e(i) ).
%C Also, numbers k such that k divides k', the arithmetic derivative of k. As shown by Ufnarovski and Ahlander, all terms in this sequence have the form Product_{j=1..r} (pj^pj)^ej, where the pj are primes. The quotient k'/k equals Sum_{j=1..r} ej. - _T. D. Noe_, Jan 04 2006
%C Multiplicative closure of A051674. - _Reinhard Zumkeller_, Jan 21 2012
%C The number of terms < 10^k: 2, 5, 9, 15, 25, 36, 52, 73, 98, 128, 167, 213, 270, 338, 421, 517, 632, 768, 920, 1101, ..., . - _Robert G. Wilson v_, Jan 19 2016
%D See A003415.
%H Robert G. Wilson v, <a href="/A072873/b072873.txt">Table of n, a(n) for n = 1..10000</a> (terms 1..2500 from Nathaniel Johnston)
%F A124010(a(n),k) mod A027748(a(n),k) = 0 for k = 1 .. A001221(a(n)). - _Reinhard Zumkeller_, Jan 21 2012
%F Sum_{n>=1} 1/a(n) = Product_{p prime} p^p/(p^p-1) = 1.38506028520448917638... - _Amiram Eldar_, Sep 27 2020
%e 108 is in the sequence because 108 = 2^2*3^3 and 2/2 + 3/3 = 2 is an integer.
%t Select[Range[1000000],IntegerQ[Total[#[[2]]/#[[1]]&/@FactorInteger[#]]]&] (* _Harvey P. Dale_, Jul 04 2014 *)
%t lst = {}; Do[n = 2^e2*3^e3*5^e5*7^e7; If[n < 10^11, AppendTo[lst, n]], {e2, 0, 36, 2}, {e3, 0, 23, 3}, {e5, 0, 15, 5}, {e7, 0, 13, 7}]; Take[ Sort@ lst, 40] (* _Robert G. Wilson v_, Jan 19 2016 *)
%o (Haskell)
%o import Data.Set (empty, fromList, deleteFindMin, union)
%o import qualified Data.Set as Set (null)
%o a072873 n = a072873_list !! (n-1)
%o a072873_list = 1 : h empty [1] a051674_list where
%o h s mcs xs'@(x:xs)
%o | Set.null s || x < m = h (s `union` fromList (map (* x) mcs)) mcs xs
%o | otherwise = m : h (s' `union` fromList (map (* m) $ init (m:mcs)))
%o (m:mcs) xs'
%o where (m, s') = deleteFindMin s
%o -- _Reinhard Zumkeller_, Jan 21 2012
%o (PARI) is(n)=my(f=factor(n)); for(i=1,#f~,if(f[i,2]%f[i,1], return(0))); 1 \\ _Charles R Greathouse IV_, Oct 28 2014
%o (Python)
%o from itertools import count, islice
%o from sympy import factorint
%o def A072873_gen(startvalue=1): # generator of terms >= startvalue
%o return (k for k in count(max(startvalue,1)) if not any(e%p for p, e in factorint(k).items()))
%o A072873_list = list(islice(A072873_gen(),20)) # _Chai Wah Wu_, Sep 15 2023
%Y Cf. A003415, A051674, A027748, A085731, A048102, A124010.
%K nonn
%O 1,2
%A _Benoit Cloitre_, Jul 28 2002
%E More terms from _T. D. Noe_, Jan 04 2006
|
