|
%I #46 Jul 10 2022 03:56:01
%S 4,9,16,25,49,64,81,121,169,289,361,529,625,729,841,961,1024,1369,
%T 1681,1849,2209,2401,2809,3481,3721,4096,4489,5041,5329,6241,6889,
%U 7921,9409,10201,10609,11449,11881,12769,14641,15625,16129,17161,18769,19321
%N Prime powers with special exponents: q^(p-1) where p > 2 and q are prime numbers.
%C Composite numbers with a prime number of divisors.
%H T. D. Noe, <a href="/A036454/b036454.txt">Table of n, a(n) for n = 1..1000</a>
%F d(d(a(n))) = 2, where d(x) = tau(x) = sigma_0(x) is the number of divisors of x.
%F a(n) = (n log n)^2 + 2n^2 log n log log n + O(n^2 log n). - _Charles R Greathouse IV_, Apr 26 2012
%F (1 - A010051(a(n))) * A010055(a(n)) * A010051(A100995(a(n))+1) = 1. - _Reinhard Zumkeller_, Jun 05 2013
%F A036459(a(n)) = 2. - _Ivan Neretin_, Jan 25 2016
%F a(n) = A283262(n)^2. - _Amiram Eldar_, Jul 04 2022
%F Sum_{n>=1} 1/a(n) = Sum_{k>=2} P(prime(k)-1) = 0.54756961912815344341..., where P is the prime zeta function. - _Amiram Eldar_, Jul 10 2022
%e From powers of 2: 4,16,64,1024,4096,65536 are in the sequence since exponent+1 is also prime. The same powers of any prime base are also included.
%p N:= 10^5:
%p P1:= select(isprime,[2,seq(2*i+1,i=1..floor((sqrt(N)-1)/2))]):
%p P2:= select(`<=`,P1,1+ilog2(N))[2..-1]:
%p S:= {seq(seq(p^(q-1), q = select(`<=`,P2,1+floor(log[p](N)))),p=P1)}:
%p sort(convert(S,list)); # _Robert Israel_, May 18 2015
%t specialPrimePowerQ[n_] := With[{f = FactorInteger[n]}, Length[f] == 1 && PrimeQ[f[[1, 1]]] && f[[1, 2]] > 1 && PrimeQ[f[[1, 2]] + 1]]; Select[Range[20000], specialPrimePowerQ] (* _Jean-François Alcover_, Jul 02 2013 *)
%t Select[Range[20000], ! PrimeQ[#] && PrimeQ[DivisorSigma[0, #]] &] (* _Carlos Eduardo Olivieri_, May 18 2015 *)
%o (PARI) for(n=1,34000, if(isprime(n), n++,x=numdiv(n); if(isprime(x),print(n))))
%o (PARI) list(lim)=my(v=List(),t);lim=lim\1+.5;forprime(p=3,log(lim)\log(2) +1, t=p-1; forprime(q=2,lim^(1/t),listput(v,q^t))); vecsort(Vec(v))
%o \\ _Charles R Greathouse IV_, Apr 26 2012
%o (Haskell)
%o a009087 n = a009087_list !! (n-1)
%o a009087_list = filter ((== 1) . a010051 . (+ 1) . a100995) a000961_list
%o -- _Reinhard Zumkeller_, Jun 05 2013
%o (Magma) [n: n in [1..20000] | not IsPrime(n) and IsPrime(DivisorSigma(0, n))]; // _Vincenzo Librandi_, May 19 2015
%Y Intersection of A002808 and A009087.
%Y Cf. A000005, A010051, A010553, A010055, A100995, A036450, A036452, A036459, A283262.
%K nonn,easy,nice
%O 1,1
%A _Labos Elemer_
|
