cp's OEIS Frontend

A076446 Differences of consecutive powerful numbers (definition 1).

%I A076446 #23 Feb 16 2025 08:32:47
%S A076446 3,4,1,7,9,2,5,4,13,15,8,9,19,8,13,4,3,16,25,27,4,16,9,18,13,32,1,35,
%T A076446 19,18,31,8,32,9,43,16,12,17,47,49,23,27,1,53,55,16,41,23,36,61,7,4,
%U A076446 28,24,65,36,27,4,69,71,27,8,21,17,3,72,77,47,32,81,47,36,36,49,87,8
%N A076446 Differences of consecutive powerful numbers (definition 1).
%C A076446 The term 1 appears infinitely often. Erdős conjectured that two consecutive 1's do not occur. (see Guy).
%D A076446 R. K. Guy, Unsolved Problems in Number Theory, B16
%H A076446 Reinhard Zumkeller, <a href="/A076446/b076446.txt">Table of n, a(n) for n = 1..10000</a>
%H A076446 Eric Weisstein's World of Mathematics, <a href="https://mathworld.wolfram.com/PowerfulNumber.html">Powerful numbers</a>
%F A076446 a(n) = A001694(n+1) - A001694(n).
%e A076446 The first two powerful numbers are 1 and 4, their difference is 3, so a(1)=3.
%t A076446 Differences[Join[{1},Select[Range[2000],Min[FactorInteger[#][[All, 2]]]>1&]]] (* _Harvey P. Dale_, Aug 27 2017 *)
%o A076446 (Haskell)
%o A076446 a076446 n = a076446_list !! (n-1)
%o A076446 a076446_list = zipWith (-) (tail a001694_list) a001694_list
%o A076446 -- _Reinhard Zumkeller_, Nov 30 2012
%o A076446 (Python)
%o A076446 from math import isqrt
%o A076446 from sympy import mobius, integer_nthroot
%o A076446 def A076446(n):
%o A076446     def squarefreepi(n): return int(sum(mobius(k)*(n//k**2) for k in range(1, isqrt(n)+1)))
%o A076446     def bisection(f,kmin=0,kmax=1):
%o A076446         while f(kmax) > kmax: kmax <<= 1
%o A076446         while kmax-kmin > 1:
%o A076446             kmid = kmax+kmin>>1
%o A076446             if f(kmid) <= kmid:
%o A076446                 kmax = kmid
%o A076446             else:
%o A076446                 kmin = kmid
%o A076446         return kmax
%o A076446     def f(x):
%o A076446         c, l = n+x, 0
%o A076446         j = isqrt(x)
%o A076446         while j>1:
%o A076446             k2 = integer_nthroot(x//j**2,3)[0]+1
%o A076446             w = squarefreepi(k2-1)
%o A076446             c -= j*(w-l)
%o A076446             l, j = w, isqrt(x//k2**3)
%o A076446         c -= squarefreepi(integer_nthroot(x,3)[0])-l
%o A076446         return c
%o A076446     return -(a:=bisection(f,n,n))+bisection(lambda x:f(x)+1,a,a) # _Chai Wah Wu_, Sep 10 2024
%Y A076446 Cf. A001694, A076444.
%K A076446 nonn
%O A076446 1,1
%A A076446 _Jud McCranie_, Oct 15 2002