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.

A056671 1 + the number of unitary and squarefree divisors of n = number of divisors of reduced squarefree part of n.

This page as a plain text file.
%I A056671 #59 Apr 15 2025 08:27:52
%S A056671 1,2,2,1,2,4,2,1,1,4,2,2,2,4,4,1,2,2,2,2,4,4,2,2,1,4,1,2,2,8,2,1,4,4,
%T A056671 4,1,2,4,4,2,2,8,2,2,2,4,2,2,1,2,4,2,2,2,4,2,4,4,2,4,2,4,2,1,4,8,2,2,
%U A056671 4,8,2,1,2,4,2,2,4,8,2,2,1,4,2,4,4,4,4,2,2,4,4,2,4,4,4,2,2,2,2,1,2,8,2,2,8
%N A056671 1 + the number of unitary and squarefree divisors of n = number of divisors of reduced squarefree part of n.
%C A056671 Note that 1 is regarded as free of squares of primes and is also a square number and a unitary divisor.
%H A056671 Antti Karttunen, <a href="/A056671/b056671.txt">Table of n, a(n) for n = 1..10000</a>
%H A056671 Steven R. Finch, <a href="/A007947/a007947.pdf">Unitarism and Infinitarism</a>, February 25, 2004. [Cached copy, with permission of the author]
%H A056671 <a href="/index/Eu#epf">Index entries for sequences computed from exponents in factorization of n</a>.
%F A056671 a(n) = A000005(A055231(n)) = A000005(A007913(n)/A055229(n)).
%F A056671 Multiplicative with a(p) = 2 and a(p^e) = 1 for e > 1. a(n) = 2^A056169(n). - _Vladeta Jovovic_, Nov 01 2001
%F A056671 a(n) = A034444(n) - A056674(n). - _Antti Karttunen_, Jul 19 2017
%F A056671 From _Vaclav Kotesovec_, Feb 11 2023: (Start)
%F A056671 Dirichlet g.f.: zeta(s) * Product_{primes p} (1 + 1/p^s - 1/p^(2*s)).
%F A056671 Dirichlet g.f.: zeta(s)^2 * Product_{primes p} (1 - 2/p^(2*s) + 1/p^(3*s)), (with a product that converges for s=1).
%F A056671 Let f(s) = Product_{primes p} (1 - 2/p^(2*s) + 1/p^(3*s)), then Sum_{k=1..n} a(k) ~ n * (f(1) * (log(n) + 2*gamma - 1) + f'(1)), where f(1) = Product_{primes p} (1 - 2/p^2 + 1/p^3) = A065464 = 0.42824950567709444021876..., f'(1) = f(1) * Sum_{primes p} (4*p-3) * log(p) / (p^3 - 2*p + 1) = 0.808661108949590913395... and gamma is the Euler-Mascheroni constant A001620. (End)
%F A056671 a(n) = Sum_{d|n, gcd(d,n/d)=1} mu(d)^2. - _Wesley Ivan Hurt_, May 25 2023
%F A056671 a(n) = Sum_{d|n} A343443(d)*mu(n/d). - _Ridouane Oudra_, Dec 18 2023
%e A056671 n = 252 = 2*2*3*3*7 has 18 divisors, 8 unitary and 8 squarefree divisors of which 2 are unitary and squarefree, divisors {1,7};
%e A056671 n = 2520 = 2*2*2*3*3*5*7 has 48 divisors, 16 unitary and 16 squarefree divisors of which {1,5,7,35} are both, thus a(2520) = 4.
%e A056671 a(2520) = a(2^3*3^2*5*7) = a(2^3)*a(3^2)*a(5)*a(7) = 1*1*2*2 = 4.
%t A056671 Array[DivisorSigma[0, #] &@ Denominator[#/Apply[Times, FactorInteger[#][[All, 1]]]^2] &, 105] (* or *)
%t A056671 Table[DivisorSum[n, 1 &, And[SquareFreeQ@ #, CoprimeQ[#, n/#]] &], {n, 105}] (* _Michael De Vlieger_, Jul 19 2017 *)
%t A056671 f[p_,e_] := If[e==1, 2, 1]; a[1] = 1; a[n_] := Times @@ (f @@@ FactorInteger[n]); Array[a, 100] (* _Amiram Eldar_, May 14 2019 *)
%o A056671 (PARI)
%o A056671 A057521(n) = { my(f=factor(n)); prod(i=1, #f~, if(f[i, 2]>1, f[i, 1]^f[i, 2], 1)); } \\ _Charles R Greathouse IV_, Aug 13 2013
%o A056671 A055231(n) = n/A057521(n);
%o A056671 A056671(n) = numdiv(A055231(n));
%o A056671 \\ Or:
%o A056671 A055229(n) = { my(c=core(n)); gcd(c, n/c); }; \\ This function from _Charles R Greathouse IV_, Nov 20 2012
%o A056671 A056671(n) = numdiv(core(n)/A055229(n)); \\ _Antti Karttunen_, Jul 19 2017
%o A056671 (PARI) for(n=1, 100, print1(direuler(p=2, n, (1 + X - X^2)/(1-X))[n], ", ")) \\ _Vaclav Kotesovec_, Feb 11 2023
%o A056671 (PARI) a(n) = vecprod(apply(x -> if(x == 1, 2, 1), factor(n)[, 2])); \\ _Amiram Eldar_, Apr 15 2025
%o A056671 (Scheme) (define (A056671 n) (if (= 1 n) n (* (if (= 1 (A067029 n)) 2 1) (A056671 (A028234 n))))) ;; (After the given multiplicative formula) - _Antti Karttunen_, Jul 19 2017
%o A056671 (Python)
%o A056671 from sympy import factorint, prod
%o A056671 def a(n): return 1 if n==1 else prod([2 if e==1 else 1 for p, e in factorint(n).items()])
%o A056671 print([a(n) for n in range(1, 51)]) # _Indranil Ghosh_, Jul 19 2017
%Y A056671 Cf. A000005, A007913, A034444, A055229, A055231, A056169, A056674.
%Y A056671 Cf. A343443.
%K A056671 mult,nonn
%O A056671 1,2
%A A056671 _Labos Elemer_, Aug 10 2000