A113849 Numbers whose prime factors are raised to the fourth power.

16, 81, 625, 1296, 2401, 10000, 14641, 28561, 38416, 50625, 83521, 130321, 194481, 234256, 279841, 456976, 707281, 810000, 923521, 1185921, 1336336, 1500625, 1874161, 2085136, 2313441, 2825761, 3111696, 3418801, 4477456, 4879681, 6765201, 7890481

Offset: 1

Author

Cino Hilliard, Jan 25 2006

Keywords

• easy
• nonn

Comments

This is essentially A005117 (the squarefree numbers) raised to the fourth power. - T. D. Noe, Mar 13 2013

All positive integers have a unique factorization into powers of squarefree numbers with distinct exponents that are powers of two. So every positive number is a product of at most one squarefree number (A005117), at most one square of a squarefree number (A062503), at most one 4th power of a squarefree number (term of this sequence), at most one 8th power of a squarefree number, and so on. - Peter Munn, Mar 12 2020

Examples

			1296 = 16*81 = 2^4*3^4 so the prime factors of 1296, 2 and 3, are raised to the fourth power.
		

Links

• T. D. Noe, Table of n, a(n) for n = 1..10000

Crossrefs

Proper subset of A000583.

Other powers of squarefree numbers: A005117(1), A062503(2), A062838(3), A113850(5), A113851(6), A113852(7), A072774(all).

Cf. A000351, A225546.

Programs

Mathematica

Select[ Range@50^4, Union[Last /@ FactorInteger@# ] == {4} &] (* Robert G. Wilson v, Jan 26 2006 *)
nn = 50; t = Select[Range[2, nn], Union[Transpose[FactorInteger[#]][[2]]] == {1} &]; t^4 (* T. D. Noe, Mar 13 2013 *)
Rest[Select[Range[100], SquareFreeQ]^4] (* Vaclav Kotesovec, May 22 2020 *)

PARI

allpwrfact(n,p) = { local(x,j,ln,y,flag); for(x=4,n, y=Vec(factor(x)); ln = length(y[1]); flag=0; for(j=1,ln, if(y[2][j]==p,flag++); ); if(flag==ln,print1(x",")); ) } \\ All prime factors are raised to the power p

Python

from math import isqrt
from sympy import mobius
def A113849(n):
    def f(x): return n+x-sum(mobius(k)*(x//k**2) for k in range(1, isqrt(x)+1))
    kmin, kmax = 1,2
    while f(kmax) >= kmax:
        kmax <<= 1
    while True:
        kmid = kmax+kmin>>1
        if f(kmid) < kmid:
            kmax = kmid
        else:
            kmin = kmid
        if kmax-kmin <= 1:
            break
    return kmax**4 # Chai Wah Wu, Aug 19 2024

Formula

From Peter Munn, Oct 31 2019: (Start)

a(n) = A005117(n+1)^4.

{a(n)} = {A225546(A000351(n)) : n >= 0} \ {1}, where {a(n)} denotes the set of integers in the sequence.

(End)

Sum_{k>=1} 1/a(k) = zeta(4)/zeta(8) - 1 = 105/Pi^4 - 1. - Amiram Eldar, May 22 2020

Extensions

More terms from Robert G. Wilson v, Jan 26 2006