cp's OEIS Frontend

A035166 Let d(m) = denominator of Sum_{k=1..m} 1/k^2 and consider f(m) = product of primes which appear to odd powers in d(m); sequence lists m such that f(m) is different from f(m-1).

%I A035166 #25 Sep 06 2023 14:58:46
%S A035166 1,10,15,20,25,42,49,50,55,66,75,78,91,100,110,121,125,136,153,156,
%T A035166 164,169,171,182,189,190,205,250,253,272,276,289,294,342,343,354,361,
%U A035166 375,406,413,435,465,473,496,500,506,516,529,555,592,605,625
%N A035166 Let d(m) = denominator of Sum_{k=1..m} 1/k^2 and consider f(m) = product of primes which appear to odd powers in d(m); sequence lists m such that f(m) is different from f(m-1).
%C A035166 The prime 479 first appears in f(m) at m = 2395, ahead of 71, which first appears in f(2485).
%C A035166 The first occurrence of four distinct primes is at m = 2500, with 5^7, 17^3, 71 and 479.
%C A035166 For 1890 < m < 2006, d(m) is a square (f(m)=1). The lone prime in 1875 .. 1890 is 61 and in 2006 .. 2027 it is 59.
%C A035166 It appears that f(m) can differ from f(m-1) in at most one prime.
%C A035166 (f from definition) = A007913, squarefree part. - _Reinhard Zumkeller_, Jul 06 2012
%H A035166 Reinhard Zumkeller, <a href="/A035166/b035166.txt">Table of n, a(n) for n = 1..250</a>
%e A035166 f(10) = 5 is the first time f(m) > 1. The 5 persists until it disappears at m = 15.
%t A035166 d[n_] := Denominator[ HarmonicNumber[n, 2]]; f[n_] := Times @@ Select[ FactorInteger[d[n]], OddQ[#[[2]]]&][[All, 1]]; A035166 = Join[{1}, Select[ Range[1000], f[#] != f[#-1]&]] (* _Jean-François Alcover_, Feb 26 2016 *)
%o A035166 (Macsyma) for k:1 do (subset(factor_number(denom(harmonic(k,2))), lambda([x],oddp(second(x)))), if old#old:%% then print(k,%%))
%o A035166 (Haskell)
%o A035166 import Data.List (findIndices)
%o A035166 a035166 n = a035166_list !! (n-1)
%o A035166 a035166_list = map (+ 1) $ findIndices (/= 0) $ zipWith (-) (tail gs) gs
%o A035166    where gs = 0 : map a007913 a007407_list
%o A035166 -- _Reinhard Zumkeller_, Jul 06 2012
%o A035166 (PARI) d(m) = denominator(sum(k=1, m, 1/k^2));
%o A035166 f(m) = my(f=factor(d(m))); for (k=1, #f~, if (!(f[k,2] % 2), f[k,2] = 0)); factorback(f);
%o A035166 isok(m) = if (m==1, 1, f(m) != f(m-1)); \\ _Michel Marcus_, Sep 06 2023
%Y A035166 Cf. A075326, A075327, A007407.
%K A035166 nonn,nice
%O A035166 1,2
%A A035166 _Bill Gosper_, Sep 04 2002