A374291 Squares of powerful numbers.

1, 16, 64, 81, 256, 625, 729, 1024, 1296, 2401, 4096, 5184, 6561, 10000, 11664, 14641, 15625, 16384, 20736, 28561, 38416, 40000, 46656, 50625, 59049, 65536, 82944, 83521, 104976, 117649, 130321, 153664, 160000, 186624, 194481, 234256, 250000, 262144, 279841, 331776

Offset: 1

Author

Amiram Eldar, Jul 02 2024

Keywords

• nonn
• easy

Comments

First differs from A340588 at n = 12.

4-full (or 3-full) squares.

Numbers whose exponents in their prime factorization are all even numbers >= 4.

This sequence is closed under multiplication.

The sequence {A000290(n)*A078615(A000290(n)), n>=1} is a permutation of this sequence, and the sequence {a(n)/A078615(a(n)), n>=1} is a permutation of {A000290(n), n>=1}.

The sequence {A335988(n)*A007947(A335988(n)), n>=1} is a permutation of this sequence, and the sequence {a(n)/A007947(a(n)), n>=1} is a permutation of A335988.

Links

• Amiram Eldar, Table of n, a(n) for n = 1..10000
• Index entries for sequences related to powerful numbers.

Crossrefs

Intersection of A000290 and A036967 (or A036966).

Intersection of A000290 and A337050.

Subsequence of A322449.

Cf. A000290, A001694, A007947, A078615, A335988, A340588.

Programs

Mathematica

powQ[n_] := n==1 || AllTrue[FactorInteger[n][[;; , 2]], # > 1 &]; Select[Range[600], powQ]^2

PARI

is(k) = issquare(k) && ispowerful(sqrtint(k));

Python

from math import isqrt
from sympy import mobius, integer_nthroot
def A374291(n):
    def squarefreepi(n):
        return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))
    def bisection(f, kmin=0, kmax=1):
        while f(kmax) > kmax: kmax <<= 1
        while kmax-kmin > 1:
            kmid = kmax+kmin>>1
            if f(kmid) <= kmid:
                kmax = kmid
            else:
                kmin = kmid
        return kmax
    def f(x):
        c, l = n+x, 0
        j = isqrt(x)
        while j>1:
            k2 = integer_nthroot(x//j**2, 3)[0]+1
            w = squarefreepi(k2-1)
            c -= j*(w-l)
            l, j = w, isqrt(x//k2**3)
        c -= squarefreepi(integer_nthroot(x, 3)[0])-l
        return c
    return bisection(f,n,n)**2 # Chai Wah Wu, Sep 10 2024

Formula

a(n) = A000290(A001694(n)) = A001694(n)^2.

Sum_{n>=1} 1/a(n) = Product_{p prime} (1 + 1/(p^2*(p^2-1))) = zeta(4)*zeta(6)/zeta(12) = 15015/(1382*Pi^2) = 1.10082313486953808844... .

Sum_{n>=1} 1/a(n)^s = Product_{p prime} (1 + 1/(p^(2*s)*(p^(2*s)-1))) = zeta(4*s)*zeta(6*s)/zeta(12*s), for s > 1/4.