%I #39 Aug 05 2024 16:03:38
%S 4,9,12,18,20,25,28,36,44,45,49,50,52,60,63,68,75,76,84,90,92,98,99,
%T 100,116,117,121,124,126,132,140,147,148,150,153,156,164,169,171,172,
%U 175,180,188,196,198,204,207,212,220,225,228
%N Cubefree numbers which are not squarefree.
%C a(n)=m iff A051903(m)=2.
%C Let us introduce a function D(n)=sigma_0(n)/(2^(alpha(1)+...+alpha(r)), sigma_0(n) number of divisors of n (A000005), prime factorization of n=p(1)^alpha(1) * ... * p(r)^alpha(r), alpha(1)+...+alpha(r) is sequence (A086436). This function splits the set of positive integers into subsets, according to the value of D(n). Squarefree numbers (005117) has D(n)=1, other numbers are "deviated" from the squarefree ideal and have 0 < D(n) < 1. So for D(n)=1/2 we have A048109, D(n)=3/4 we have A067295. - _Ctibor O. Zizka_, Sep 21 2008
%C The asymptotic density of this sequence is 1/zeta(3) - 1/zeta(2) = A088453 - A059956 = 0.22398... - _Amiram Eldar_, Jul 09 2020
%H Reinhard Zumkeller, <a href="/A067259/b067259.txt">Table of n, a(n) for n = 1..10000</a>
%H Eric Weisstein's World of Mathematics, <a href="http://mathworld.wolfram.com/Cubefree.html">Cubefree</a>
%H Eric Weisstein's World of Mathematics, <a href="http://mathworld.wolfram.com/Squarefree.html">Squarefree</a>
%F A212793(a(n)) * (1 - A008966(a(n))) = 1. - _Reinhard Zumkeller_, May 27 2012
%t f[n_]:=Union[Last/@FactorInteger[n]][[ -1]]; lst={}; Do[If[f[n]==2,AppendTo[lst,n]],{n,2,6!}]; lst (* _Vladimir Joseph Stephan Orlovsky_, Feb 12 2010 *)
%t Select[Range[500], Not[SquareFreeQ[#]] && FreeQ[FactorInteger[#], {_,k_ /;k>2}]&] (* _Vaclav Kotesovec_, Jul 09 2020 *)
%o (Haskell)
%o a067259 n = a067259_list !! (n-1)
%o a067259_list = filter ((== 2) . a051903) [1..]
%o -- _Reinhard Zumkeller_, May 27 2012
%o (PARI) is(n)=n>3 && vecmax(factor(n)[,2])==2 \\ _Charles R Greathouse IV_, Oct 15 2015
%o (Python)
%o from math import isqrt
%o from sympy import mobius, integer_nthroot
%o def A067259(n):
%o def f(x): return n+x+sum(mobius(k)*(x//k**2-x//k**3) for k in range(1, integer_nthroot(x,3)[0]+1))+sum(mobius(k)*(x//k**2) for k in range(integer_nthroot(x,3)[0]+1,isqrt(x)+1))
%o m, k = n, f(n)
%o while m != k:
%o m, k = k, f(k)
%o return m # _Chai Wah Wu_, Aug 05 2024
%Y Cf. A004709, A005117, A059956, A088453.
%K nonn
%O 1,1
%A _Reinhard Zumkeller_, Feb 20 2002