login
Indices of records in A018804.
0

%I #66 Aug 12 2023 01:20:00

%S 1,2,3,4,5,6,8,9,10,12,15,16,18,20,24,28,30,36,40,42,48,54,56,60,72,

%T 80,84,90,105,108,120,140,144,168,180,210,240,252,270,280,288,300,330,

%U 336,360,420,480,504,540,600,630,660,720,840,990,1008,1080,1200

%N Indices of records in A018804.

%C a(n) seems to be divisible by any positive integer as n approaches infinity. Example: There seem to be 6 terms without a 2 (i.e., 1, 3, 5, 9, 15, 105), 14 terms without a 3 (as large as 280), 22 terms without a factor of 4 (as large as 6930), and 29 terms without a 5 (as large as 3276). True for a(n) < 10^8.

%C For any term a(n) > 1, it seems that there exists at least one term a(x) such that a(x) * prime(y) = a(n) where prime(y) <= (7/2) * A053669(a(x)). True for a(n) < 10^8.

%C The ratio of prime(n)/a(n) for indices {1..8, 11, 15, 17} is almost 2 (with an error of at most 1 on the numerator). The exact ratios are 2/1, 3/2, 5/3, 7/4, 11/5, 13/6, 17/8, 19/9, 31/15, 47/24, 59/30. - _Matthew Russell Downey_, Jul 25 2023

%H Matthew Russell Downey, <a href="https://gitlab.com/planetaria/onbifurcationsandbeauty/-/blob/87a13e45fb7f28ada6449de3b2990810616c912c/OnBifurcationsAndBeauty.pdf">On Bifurcations And Beauty</a> (2023): pages 40-43.

%F Conjecture: a(n) ~ exp(4*(n-1)/21). - _Matthew Russell Downey_, Jul 25 2023

%e A018804(36) = 168 is the largest value among the first 36 terms of A018804, so 36 is a term here; since it is the 18th value that sets a new record, a(18) = 36.

%t f[p_, e_] := (e*(p - 1)/p + 1)*p^e; pil[n_] := Times @@ (f @@@ FactorInteger[n]); seq[nmax_] := Module[{p, pm = 0, s = {}}, Do[If[(p = pil[n]) > pm, pm = p; AppendTo[s, n]], {n, 1, nmax}]; s]; seq[1200] (* _Amiram Eldar_, Feb 13 2023 *)

%o (Python)

%o from sympy import factorint

%o from math import prod

%o def A018804(m): return prod(p**(e-1)*((p-1)*e+p) for p, e in factorint(m).items())

%o record = 0

%o for m in range(1, 2000):

%o value = A018804(m)

%o if value > record:

%o record = value

%o print(m, end=", ")

%o (PARI) f(n) = sumdiv(n, d, n*eulerphi(d)/d); \\ A018804

%o lista(nn) = my(r=0, list=List()); for (n=1, nn, my(m=f(n)); if (m > r, listput(list, n); r = m);); Vec(list); \\ _Michel Marcus_, Feb 26 2023

%Y Indices of records in A018804. Cf. A000040.

%K nonn

%O 1,2

%A _Matthew Russell Downey_, Feb 08 2023