|
%I #31 Jun 17 2021 09:13:46
%S 1,2,3,4,5,0,6,7,8,0,9,0,10,0,0,11,12,0,13,0,0,0,14,0,15,0,16,0,17,0,
%T 18,19,0,0,0,0,20,0,0,0,21,0,22,0,0,0,23,0,24,0,0,0,25,0,0,0,0,0,26,0,
%U 27,0,0,28,0,0,29,0,0,0,30,0,31,0,0,0,0,0,32,0,33,0,34,0,0,0,0,0,35,0
%N a(n) = k if n = A000961(k) (powers of primes), a(n) = 0 if n is not in A000961.
%C The name has been edited to clarify that the indices k refer to A000961 ("powers of primes" = {1} U A246655) and not to the list A246655 of proper prime powers. - _M. F. Hasler_, Jun 16 2021
%H Reinhard Zumkeller, <a href="/A095874/b095874.txt">Table of n, a(n) for n = 1..10000</a>
%F a(n) = Sum_{1 <= k <= n} A010055(k); [corrected by _M. F. Hasler_, Jun 15 2021]
%F a(n) = A065515(n)*(A065515(n)-A065515(n-1)).
%F a(n) = A065515(n)*A069513(n). - _M. F. Hasler_, Jun 16 2021
%t Join[{1},Module[{k=2},Table[If[PrimePowerQ[n],k;k++,0],{n,2,100}]]] (* _Harvey P. Dale_, Aug 15 2020 *)
%o (Haskell)
%o a095874 n | y == n = length xs + 1
%o | otherwise = 0
%o where (xs, y:ys) = span (< n) a000961_list
%o -- _Reinhard Zumkeller_, Feb 16 2012, Jun 26 2011
%o (PARI) a(n)=if(isprimepower(n), sum(i=1,logint(n,2), primepi(sqrtnint(n,i)))+1, n==1) \\ _Charles R Greathouse IV_, Apr 29 2015
%o (PARI) {M95874=Map(); A095874(n,k)=if(mapisdefined(M95874,n,&k),k, isprimepower(n), mapput(M95874,n, k=sum(i=1,exponent(n), primepi(sqrtnint(n,i)))+1); k,n==1)} \\ Variant with memoization, possibly useful to compute A097621, A344826 and related. One may omit "isprimepower(n)," (possibly requiring factorization) and ",n==1" if n is known to be a power of a prime, i.e., to get a left inverse for A000961. - _M. F. Hasler_, Jun 15 2021
%Y Cf. A000961 (right inverse), A049084, A097621.
%K nonn
%O 1,2
%A _Reinhard Zumkeller_, Jun 10 2004
%E Edited by _M. F. Hasler_, Jun 15 2021
|
