cp's OEIS Frontend

A358173 First differences of A286708.

%I A358173 #14 Sep 11 2024 00:34:35
%S A358173 36,28,8,36,52,4,16,9,63,36,68,8,32,9,43,16,76,72,27,1,108,16,64,36,
%T A358173 68,4,28,89,36,27,4,69,71,27,29,20,72,77,47,32,128,36,36,136,8,56,25,
%U A358173 91,188,8,188,92,9,99,4,40,144,28,109,62,49,64,49,18,97,11,81
%N A358173 First differences of A286708.
%C A358173 Consider the sequence of powerful numbers A001694, superset of A246547, the sequence of composite prime powers. Let s = A001694(k) such that omega(s) > 1 be followed by t = A001694(k+1) such that omega(t) = 1.
%C A358173 Since A286708 = A001694 \ A246547, prime powers t are missing in A286708. We consider s = A286708(j) and note that the difference A286708(j+1) - A286708(j) > A001694(k+1) - A001694(k).
%C A358173 Therefore we see a subset S containing s in A286708 that plots "out of place" with respect to the complementary subset R = A286708 \ S; some of this subset S exceeds the maxima of R in the scatterplot of this sequence. The plot of the R resembles the scatterplot of A001694.
%H A358173 Michael De Vlieger, <a href="/A358173/b358173.txt">Table of n, a(n) for n = 1..10000</a>
%H A358173 Michael De Vlieger, <a href="/A358173/a358173.png">Scatterplot of a(n)</a>, n = 1..2^16, highlighting in red terms s that in A001694 are succeeded by a prime power t.
%e A358173 The number 36 is the smallest powerful number that is not a prime power; the next powerful number that is not a prime power is 72, and their difference is 36, hence a(1) = 36.
%t A358173 With[{nn = 2^25}, Differences@ Select[Rest@ Union@ Flatten@ Table[a^2*b^3, {b, nn^(1/3)}, {a, Sqrt[nn/b^3]}], ! PrimePowerQ[#] &]]
%o A358173 (Python)
%o A358173 from math import isqrt
%o A358173 from sympy import integer_nthroot, primepi, mobius
%o A358173 def A358173(n):
%o A358173     def squarefreepi(n): return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))
%o A358173     def bisection(f, kmin=0, kmax=1):
%o A358173         while f(kmax) > kmax: kmax <<= 1
%o A358173         while kmax-kmin > 1:
%o A358173             kmid = kmax+kmin>>1
%o A358173             if f(kmid) <= kmid:
%o A358173                 kmax = kmid
%o A358173             else:
%o A358173                 kmin = kmid
%o A358173         return kmax
%o A358173     def f(x):
%o A358173         c, l = n+x+1+sum(primepi(integer_nthroot(x, k)[0]) for k in range(2, x.bit_length())), 0
%o A358173         j = isqrt(x)
%o A358173         while j>1:
%o A358173             k2 = integer_nthroot(x//j**2,3)[0]+1
%o A358173             w = squarefreepi(k2-1)
%o A358173             c -= j*(w-l)
%o A358173             l, j = w, isqrt(x//k2**3)
%o A358173         c -= squarefreepi(integer_nthroot(x,3)[0])-l
%o A358173         return c
%o A358173     return -(a:=bisection(f,n,n))+bisection(lambda x:f(x)+1,a,a) # _Chai Wah Wu_, Sep 10 2024
%Y A358173 Cf. A001694, A053707, A076446, A246547, A286708.
%K A358173 nonn
%O A358173 1,1
%A A358173 _Michael De Vlieger_, Nov 01 2022