cp's OEIS Frontend

This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.

Showing 1-8 of 8 results.

A001694 Powerful numbers, definition (1): if a prime p divides n then p^2 must also divide n (also called squareful, square full, square-full or 2-powerful numbers).

Original entry on oeis.org

1, 4, 8, 9, 16, 25, 27, 32, 36, 49, 64, 72, 81, 100, 108, 121, 125, 128, 144, 169, 196, 200, 216, 225, 243, 256, 288, 289, 324, 343, 361, 392, 400, 432, 441, 484, 500, 512, 529, 576, 625, 648, 675, 676, 729, 784, 800, 841, 864, 900, 961, 968, 972, 1000
Offset: 1

Views

Author

Keywords

Comments

Numbers of the form a^2*b^3, a >= 1, b >= 1.
In other words, if the prime factorization of n is Product_k p_k^e_k then all e_k are greater than 1.
Numbers n such that Sum_{d|n} phi(d)*phi(n/d)*mu(d) > 0; places of nonzero A300717. - Benoit Cloitre, Nov 30 2002
This sequence is closed under multiplication. The primitive elements are A168363. - Franklin T. Adams-Watters, May 30 2011
Complement of A052485. - Reinhard Zumkeller, Sep 16 2011
The number of terms less than or equal to 10^k beginning with k = 0: 1, 4, 14, 54, 185, 619, 2027, 6553, 21044, ...: A118896. - Robert G. Wilson v, Aug 11 2014
a(10^n): 1, 49, 3136, 253472, 23002083, 2200079025, 215523459072, 21348015504200, 2125390162618116, ... . - Robert G. Wilson v, Aug 15 2014
a(m) mod prime(n) > 0 for m < A258599(n); a(A258599(n)) = A001248(n) = prime(n)^2. - Reinhard Zumkeller, Jun 06 2015
From Des MacHale, Mar 07 2021: (Start)
A number m is powerful if and only if |R/Z(R)| = m, for some finite non-commutative ring R.
A number m is powerful if and only if |G/Z(G)| = m, for some finite nilpotent class two group G (Reference Aine Nishe). (End)
Numbers n such that Sum_{k=1..n} phi(gcd(n,k))*mu(gcd(n,k)) > 0. - Richard L. Ollerton, May 09 2021

Examples

			1 is a term because for every prime p that divides 1, p^2 also divides 1.
2 is not a term since 2 divides 2 but 2^2 does not.
4 is a term because 2 is the only prime that divides 4 and 2^2 does divide 4. - _N. J. A. Sloane_, Jan 16 2022
		

References

  • G. E. Hardy and M. V. Subbarao, Highly powerful numbers, Congress. Numer. 37 (1983), 277-307.
  • Aleksandar Ivić, The Riemann Zeta-Function, Wiley, NY, 1985, see p. 407.
  • Richard A. Mollin, Quadratics, CRC Press, 1996, Section 1.6.
  • Aine NiShe, Commutativity and Generalisations in Finite Groups, Ph.D. Thesis, University College Cork, 2000.
  • Paulo Ribenboim, Meine Zahlen, meine Freunde, 2009, Springer, 9.1 Potente Zahlen, pp. 241-247.
  • N. J. A. Sloane, A Handbook of Integer Sequences, Academic Press, 1973 (includes this sequence).
  • N. J. A. Sloane and Simon Plouffe, The Encyclopedia of Integer Sequences, Academic Press, 1995 (includes this sequence).
  • Gérald Tenenbaum, Introduction to analytic and probabilistic number theory, Cambridge University Press, 1995, p. 54, exercise 10 (in the third edition 2015, p. 63, exercise 70).

Crossrefs

Disjoint union of A062503 and A320966.
Cf. A007532 (Powerful numbers, definition (2)), A005934, A005188, A003321, A014576, A023052 (Powerful numbers, definition (3)), A046074, A013929, A076871, A258599, A001248, A112526, A168363, A224866, A261883, A300717.
Cf. A052485 (complement), A076446 (first differences), A376361, A376362.

Programs

  • Haskell
    a001694 n = a001694_list !! (n-1)
    a001694_list = filter ((== 1) . a112526) [1..]
    -- Reinhard Zumkeller, Nov 30 2012
    
  • Maple
    isA001694 := proc(n) for p in ifactors(n)[2] do if op(2,p) = 1 then return false; end if; end do; return true; end proc:
    A001694 := proc(n) option remember; if n = 1 then 1; else for a from procname(n-1)+1 do if isA001694(a) then return a; end if; end do; end if; end proc:
    seq(A001694(n),n=1..20) ; # R. J. Mathar, Jun 07 2011
  • Mathematica
    Join[{1}, Select[ Range@ 1250, Min@ FactorInteger[#][[All, 2]] > 1 &]]
    (* Harvey P. Dale, Sep 18 2011; modified by Robert G. Wilson v, Aug 11 2014 *)
    max = 10^3; Union@ Flatten@ Table[a^2*b^3, {b, max^(1/3)}, {a, Sqrt[max/b^3]}] (* Robert G. Wilson v, Aug 11 2014 *)
    nextPowerfulNumber[n_] := Block[{r = Range[ Floor[1 + n^(1/3)]]^3}, Min@ Select[ Sort[ r*Floor[1 + Sqrt[n/r]]^2], # > n &]]; NestList[ nextPowerfulNumber, 1, 55] (* Robert G. Wilson v, Aug 16 2014 *)
  • PARI
    isA001694(n)=n=factor(n)[,2];for(i=1,#n,if(n[i]==1,return(0)));1 \\ Charles R Greathouse IV, Feb 11 2011
    
  • PARI
    list(lim,mn=2)=my(v=List(),t); for(m=1,sqrtnint(lim\1,3), t=m^3; for(n=1,sqrtint(lim\t), listput(v,t*n^2))); Set(v) \\ Charles R Greathouse IV, Jul 31 2011; edited Sep 22 2015
    
  • PARI
    is=ispowerful \\ Charles R Greathouse IV, Nov 13 2012
    
  • Python
    from sympy import factorint
    A001694 = [1]+[n for n in range(2,10**6) if min(factorint(n).values()) > 1]
    # Chai Wah Wu, Aug 14 2014
    
  • Python
    from math import isqrt
    from sympy import mobius, integer_nthroot
    def A001694(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) # Chai Wah Wu, Sep 09 2024
    
  • Sage
    sloane.A001694.list(54) # Peter Luschny, Feb 08 2015

Formula

A112526(a(n)) = 1. - Reinhard Zumkeller, Sep 16 2011
Bateman & Grosswald prove that there are zeta(3/2)/zeta(3) x^{1/2} + zeta(2/3)/zeta(2) x^{1/3} + O(x^{1/6}) terms up to x; see section 5 for a more precise error term. - Charles R Greathouse IV, Nov 19 2012
a(n) = A224866(n) - 1. - Reinhard Zumkeller, Jul 23 2013
Sum_{n>=1} 1/a(n) = zeta(2)*zeta(3)/zeta(6). - Ivan Neretin, Aug 30 2015
Sum_{n>=1} 1/a(n)^s = zeta(2*s)*zeta(3*s)/zeta(6*s), s > 1/2 (Golomb, 1970). - Amiram Eldar, Oct 02 2022

Extensions

More terms from Henry Bottomley, Mar 16 2000
Definition expanded by Jonathan Sondow, Jan 03 2016

A003557 n divided by largest squarefree divisor of n; if n = Product p(k)^e(k) then a(n) = Product p(k)^(e(k)-1), with a(1) = 1.

Original entry on oeis.org

1, 1, 1, 2, 1, 1, 1, 4, 3, 1, 1, 2, 1, 1, 1, 8, 1, 3, 1, 2, 1, 1, 1, 4, 5, 1, 9, 2, 1, 1, 1, 16, 1, 1, 1, 6, 1, 1, 1, 4, 1, 1, 1, 2, 3, 1, 1, 8, 7, 5, 1, 2, 1, 9, 1, 4, 1, 1, 1, 2, 1, 1, 3, 32, 1, 1, 1, 2, 1, 1, 1, 12, 1, 1, 5, 2, 1, 1, 1, 8, 27, 1, 1, 2, 1, 1, 1, 4, 1, 3, 1, 2, 1, 1, 1, 16, 1, 7
Offset: 1

Views

Author

Keywords

Comments

a(n) is the size of the Frattini subgroup of the cyclic group C_n - Ahmed Fares (ahmedfares(AT)my-deja.com), Jun 07 2001.
Also of the Frattini subgroup of the dihedral group with 2*n elements. - Sharon Sela (sharonsela(AT)hotmail.com), Jan 01 2002
Number of solutions to x^m==0 (mod n) provided that n < 2^(m+1), i.e. the sequence of sequences A000188, A000189, A000190, etc. converges to this sequence. - Henry Bottomley, Sep 18 2001
a(n) is the number of nilpotent elements in the ring Z/nZ. - Laszlo Toth, May 22 2009
The sequence of partial products of a(n) is A085056(n). - Peter Luschny, Jun 29 2009
The first occurrence of n in this sequence is at A064549(n). - Franklin T. Adams-Watters, Jul 25 2014
From Hal M. Switkay, Jul 03 2025: (Start)
For n > 1, a(n) is a proper divisor of n. Thus the sequence n, a(n), a(a(n)), ... eventually becomes 1. This yields a minimal factorization of n as a product of squarefree numbers (A005117), each factor dividing all larger factors, in a factorization that is conjugate to the minimal factorization of n as a product of prime powers (A000961), as follows.
Let f(n,0) = n, and let f(n,k) = a(f(n,k-1)) for k > 0. A051903(n) is the minimal value of k such that f(n,k) = 1. A051903(n) <= log(n)/log(2). Since n/a(n) = A007947(n) is always squarefree by definition, n is a product of squarefree factors in the form Product_{i=1..A051903(n)} [f(n,i-1)/f(n,i)].
The two factorizations correspond to conjugate partitions of bigomega(n) = A001222(n). (End)

Crossrefs

Cf. A007947, A062378, A062379, A064549, A300717 (Möbius transform), A326306 (inv. Möbius transf.), A328572.
Sequences that are multiples of this sequence (the other factor of a pointwise product is given in parentheses): A000010 (A173557), A000027 (A007947), A001615 (A048250), A003415 (A342001), A007434 (A345052), A057521 (A071773).
Cf. A082695 (Dgf at s=2), A065487 (Dgf at s=3).

Programs

  • Haskell
    a003557 n = product $ zipWith (^)
                          (a027748_row n) (map (subtract 1) $ a124010_row n)
    -- Reinhard Zumkeller, Dec 20 2013
    
  • Julia
    using Nemo
    function A003557(n)
        n < 4 && return 1
        q = prod([p for (p, e) ∈ Nemo.factor(fmpz(n))])
        return n == q ? 1 : div(n, q)
    end
    [A003557(n) for n in 1:90] |> println  # Peter Luschny, Feb 07 2021
  • Magma
    [(&+[(Floor(k^n/n)-Floor((k^n-1)/n)): k in [1..n]]): n in [1..100]]; // G. C. Greubel, Nov 02 2018
    
  • Maple
    A003557 := n -> n/ilcm(op(numtheory[factorset](n))):
    seq(A003557(n), n=1..98); # Peter Luschny, Mar 23 2011
    seq(n / NumberTheory:-Radical(n), n = 1..98); # Peter Luschny, Jul 20 2021
  • Mathematica
    Prepend[ Array[ #/Times@@(First[ Transpose[ FactorInteger[ # ] ] ])&, 100, 2 ], 1 ] (* Olivier Gérard, Apr 10 1997 *)
  • PARI
    a(n)=n/factorback(factor(n)[,1]) \\ Charles R Greathouse IV, Nov 17 2014
    
  • PARI
    for(n=1, 100, print1(direuler(p=2, n, (1 - p*X + X)/(1 - p*X))[n], ", ")) \\ Vaclav Kotesovec, Jun 20 2020
    
  • Python
    from sympy.ntheory.factor_ import core
    from sympy import divisors
    def a(n): return n / max(i for i in divisors(n) if core(i) == i)
    print([a(n) for n in range(1, 101)]) # Indranil Ghosh, Apr 16 2017
    
  • Python
    from math import prod
    from sympy import primefactors
    def A003557(n): return n//prod(primefactors(n)) # Chai Wah Wu, Nov 04 2022
    
  • Sage
    def A003557(n) : return n*mul(1/p for p in prime_divisors(n))
    [A003557(n) for n in (1..98)] # Peter Luschny, Jun 10 2012
    

Formula

Multiplicative with a(p^e) = p^(e-1). - Vladeta Jovovic, Jul 23 2001
a(n) = n/rad(n) = n/A007947(n) = sqrt(J_2(n)/J_2(rad(n))), where J_2(n) is A007434. - Enrique Pérez Herrero, Aug 31 2010
a(n) = (J_k(n)/J_k(rad(n)))^(1/k), where J_k is the k-th Jordan Totient Function: (J_2 is A007434 and J_3 A059376). - Enrique Pérez Herrero, Sep 03 2010
Dirichlet convolution of A000027 and A097945. - R. J. Mathar, Dec 20 2011
a(n) = A000010(n)/|A023900(n)|. - Eric Desbiaux, Nov 15 2013
a(n) = Product_{k = 1..A001221(n)} (A027748(n,k)^(A124010(n,k)-1)). - Reinhard Zumkeller, Dec 20 2013
a(n) = Sum_{k=1..n}(floor(k^n/n)-floor((k^n-1)/n)). - Anthony Browne, May 11 2016
a(n) = e^[Sum_{k=2..n} (floor(n/k)-floor((n-1)/k))*(1-A010051(k))*Mangoldt(k)] where Mangoldt is the Mangoldt function. - Anthony Browne, Jun 16 2016
a(n) = Sum_{d|n} mu(d) * phi(d) * (n/d), where mu(d) is the Moebius function and phi(d) is the Euler totient function (rephrases formula of Dec 2011). - Daniel Suteu, Jun 19 2018
G.f.: Sum_{k>=1} mu(k)*phi(k)*x^k/(1 - x^k)^2. - Ilya Gutkovskiy, Nov 02 2018
Dirichlet g.f.: Product_{primes p} (1 + 1/(p^s - p)). - Vaclav Kotesovec, Jun 24 2020
From Richard L. Ollerton, May 07 2021: (Start)
a(n) = Sum_{k=1..n} mu(n/gcd(n,k))*gcd(n,k).
a(n) = Sum_{k=1..n} mu(gcd(n,k))*(n/gcd(n,k))*phi(gcd(n,k))/phi(n/gcd(n,k)). (End)
a(n) = A001615(n)/A048250(n) = A003415/A342001(n) = A057521(n)/A071773(n). - Antti Karttunen, Jun 08 2021

Extensions

Secondary definition added to the name by Antti Karttunen, Jun 08 2021

A346485 Möbius transform of A342001, where A342001(n) = A003415(n)/A003557(n).

Original entry on oeis.org

0, 1, 1, 1, 1, 3, 1, 1, 1, 5, 1, 2, 1, 7, 6, 1, 1, 1, 1, 4, 8, 11, 1, 2, 1, 13, 1, 6, 1, 14, 1, 1, 12, 17, 10, 0, 1, 19, 14, 4, 1, 20, 1, 10, 4, 23, 1, 2, 1, 1, 18, 12, 1, 1, 14, 6, 20, 29, 1, 8, 1, 31, 6, 1, 16, 32, 1, 16, 24, 34, 1, 0, 1, 37, 2, 18, 16, 38, 1, 4, 1, 41, 1, 12, 20, 43, 30, 10, 1, 4, 18, 22, 32, 47
Offset: 1

Views

Author

Antti Karttunen, Aug 26 2021

Keywords

Comments

Conjecture 1: After the initial zero, the positions of other zeros is given by A036785.
Conjecture 2: No negative terms. Checked up to n = 2^24.

Crossrefs

Programs

Formula

a(n) = Sum_{d|n} A008683(n/d) * A342001(d).
Dirichlet g.f.: Product_{p prime} (1+p^(1-s)-p^(-s)) * Sum_{p prime} p^s/((p^s-1)*(p^s+p-1)). - Sebastian Karlsson, May 08 2022
Sum_{k=1..n} a(k) ~ c * A065464 * n^2 / 2, where c = Sum_{j>=2} (1/2 + (-1)^j * (Fibonacci(j) - 1/2))*PrimeZetaP(j) = 0.4526952873143153104685540856936425315834753528741817723313791528384... - Vaclav Kotesovec, Mar 04 2023

A317935 Numerators of rational valued sequence whose Dirichlet convolution with itself yields A003557, n divided by largest squarefree divisor of n.

Original entry on oeis.org

1, 1, 1, 7, 1, 1, 1, 25, 11, 1, 1, 7, 1, 1, 1, 363, 1, 11, 1, 7, 1, 1, 1, 25, 19, 1, 61, 7, 1, 1, 1, 1335, 1, 1, 1, 77, 1, 1, 1, 25, 1, 1, 1, 7, 11, 1, 1, 363, 27, 19, 1, 7, 1, 61, 1, 25, 1, 1, 1, 7, 1, 1, 11, 9923, 1, 1, 1, 7, 1, 1, 1, 275, 1, 1, 19, 7, 1, 1, 1, 363, 1363, 1, 1, 7, 1, 1, 1, 25, 1, 11, 1, 7, 1, 1, 1, 1335, 1, 27, 11, 133, 1, 1, 1, 25, 1
Offset: 1

Views

Author

Antti Karttunen, Aug 12 2018

Keywords

Comments

Multiplicative because A003557 is.
No negative terms among the first 2^20 terms. Is the sequence nonnegative?

Crossrefs

Cf. A003557, A046644 (denominators).
Cf. also A300717, A300719, A318317.

Programs

Formula

a(n) = numerator of f(n), where f(1) = 1, f(n) = (1/2) * (A003557(n) - Sum_{d|n, d>1, d 1.

A300718 Möbius transform of A010848, number of numbers k <= n such that at least one prime factor of n is not a prime factor of k.

Original entry on oeis.org

0, 1, 2, 1, 4, 2, 6, 2, 4, 4, 10, 4, 12, 6, 8, 4, 16, 6, 18, 8, 12, 10, 22, 8, 16, 12, 12, 12, 28, 8, 30, 8, 20, 16, 24, 10, 36, 18, 24, 16, 40, 12, 42, 20, 24, 22, 46, 16, 36, 20, 32, 24, 52, 18, 40, 24, 36, 28, 58, 16, 60, 30, 36, 16, 48, 20, 66, 32, 44, 24, 70, 20, 72, 36, 40, 36, 60, 24, 78, 32, 36, 40, 82, 24, 64, 42, 56, 40, 88, 24, 72
Offset: 1

Views

Author

Antti Karttunen, Mar 11 2018

Keywords

Crossrefs

Programs

Formula

a(n) = Sum_{d|n} A008683(n/d)*A010848(d).
a(n) = A000010(n) - A300717(n).
a(n) = A010848(n) - A300720(n).

A300719 Difference between A003557 (n divided by largest squarefree divisor of n) and its Möbius transform.

Original entry on oeis.org

0, 1, 1, 1, 1, 1, 1, 2, 1, 1, 1, 2, 1, 1, 1, 4, 1, 3, 1, 2, 1, 1, 1, 4, 1, 1, 3, 2, 1, 1, 1, 8, 1, 1, 1, 4, 1, 1, 1, 4, 1, 1, 1, 2, 3, 1, 1, 8, 1, 5, 1, 2, 1, 9, 1, 4, 1, 1, 1, 2, 1, 1, 3, 16, 1, 1, 1, 2, 1, 1, 1, 8, 1, 1, 5, 2, 1, 1, 1, 8, 9, 1, 1, 2, 1, 1, 1, 4, 1, 3, 1, 2, 1, 1, 1, 16, 1, 7, 3, 6, 1, 1, 1, 4, 1
Offset: 1

Views

Author

Antti Karttunen, Mar 12 2018

Keywords

Crossrefs

Programs

  • PARI
    A003557(n) = { my(f=factor(n)); for (i=1, #f~, f[i, 2] = max(0,f[i, 2]-1)); factorback(f); }; \\ From A003557
    A300719(n) = -sumdiv(n,d,(dA003557(d));

Formula

a(n) = -Sum_{d|n, dA008683(n/d)*A003557(d).
a(n) = A003557(n) - A300717(n).

A326306 Dirichlet g.f.: zeta(s) * zeta(s-1) * Product_{p prime} (1 - p^(1 - s) + p^(-s)).

Original entry on oeis.org

1, 2, 2, 4, 2, 4, 2, 8, 5, 4, 2, 8, 2, 4, 4, 16, 2, 10, 2, 8, 4, 4, 2, 16, 7, 4, 14, 8, 2, 8, 2, 32, 4, 4, 4, 20, 2, 4, 4, 16, 2, 8, 2, 8, 10, 4, 2, 32, 9, 14, 4, 8, 2, 28, 4, 16, 4, 4, 2, 16, 2, 4, 10, 64, 4, 8, 2, 8, 4, 8, 2, 40, 2, 4, 14, 8, 4, 8, 2, 32, 41, 4, 2, 16, 4
Offset: 1

Views

Author

Ilya Gutkovskiy, Oct 17 2019

Keywords

Comments

Inverse Moebius transform of A003557.
Dirichlet convolution of A000203 with A097945.

Crossrefs

Cf. A000010, A000079 (fixed points), A000203, A003557, A007947, A008683, A098108 (parity of a(n)), A191750, A300717, A335032.

Programs

  • Mathematica
    Table[Sum[d/Last[Select[Divisors[d], SquareFreeQ]], {d, Divisors[n]}], {n, 1, 85}]
    Table[Sum[MoebiusMu[n/d] EulerPhi[n/d] DivisorSigma[1, d], {d, Divisors[n]}], {n, 1, 85}]
    f[p_, e_] := 1 + (p^e-1)/(p-1); a[1] = 1; a[n_] := Times @@ (f @@@ FactorInteger[n]); Array[a, 100] (* Amiram Eldar, Oct 14 2020 *)
  • PARI
    for(n=1, 100, print1(direuler(p=2, n, (1 - p*X + X)/(1 - X)/(1 - p*X))[n], ", ")) \\ Vaclav Kotesovec, Jun 14 2020

Formula

G.f.: Sum_{k>=1} (k / rad(k)) * x^k / (1 - x^k), where rad = A007947.
a(n) = Sum_{d|n} A003557(d).
a(n) = Sum_{d|n} mu(n/d) * phi(n/d) * sigma(d), where mu = A008683, phi = A000010 and sigma = A000203.
a(p) = 2, where p is prime.
From Vaclav Kotesovec, Jun 20 2020: (Start)
Dirichlet g.f.: zeta(s) * Product_{primes p} (1 + 1/(p^s - p)).
Dirichlet g.f.: zeta(s) * zeta(2*s-2) * Product_{primes p} (1 + p^(1-2*s) - p^(2-2*s) + p^(-s)). (End)
Multiplicative with a(p^e) = 1 + (p^e-1)/(p-1). - Amiram Eldar, Oct 14 2020
From Richard L. Ollerton, May 07 2021: (Start)
a(n) = Sum_{k=1..n} mu(n/gcd(n,k))*sigma(gcd(n,k)).
a(n) = Sum_{k=1..n} mu(gcd(n,k))*sigma(n/gcd(n,k))*phi(gcd(n,k))/phi(n/gcd(n,k)). (End)

A349379 Möbius transform of A057521 (powerful part of n).

Original entry on oeis.org

1, 0, 0, 3, 0, 0, 0, 4, 8, 0, 0, 0, 0, 0, 0, 8, 0, 0, 0, 0, 0, 0, 0, 0, 24, 0, 18, 0, 0, 0, 0, 16, 0, 0, 0, 24, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 48, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 32, 0, 0, 0, 0, 0, 0, 0, 32, 0, 0, 0, 0, 0, 0, 0, 0, 54, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 72, 0, 0
Offset: 1

Views

Author

Antti Karttunen, Nov 18 2021

Keywords

Comments

Multiplicative with a(p^e) = 0 if e = 1, p^2 - 1 if e = 2 and p^e - p^(e-1) otherwise. - Amiram Eldar, Nov 18 2021

Crossrefs

Programs

  • Mathematica
    f[p_, e_] := Which[e > 2, p^e - p^(e - 1), e == 2, p^2 - 1, e == 1, 0]; a[1] = 1; a[n_] := Times @@ f @@@ FactorInteger[n]; Array[a, 100] (* Amiram Eldar, Nov 18 2021 *)
  • PARI
    A057521(n) = { my(f=factor(n)); prod(i=1, #f~, if(f[i, 2]>1, f[i, 1]^f[i, 2], 1)); }; \\ From A057521
    A349379(n) = sumdiv(n,d,moebius(n/d)*A057521(d));
    
  • Python
    from math import prod
    from sympy import factorint
    def A349379(n): return prod(0 if e==1 else p**e - (1 if e==2 else p**(e-1)) for p,e in factorint(n).items()) # Chai Wah Wu, Nov 14 2022

Formula

a(n) = Sum_{d|n} A008683(n/d) * A057521(d).
a(n) = Sum_{d|n} A000010(n/d) * A349441(d).
Showing 1-8 of 8 results.