cp's OEIS Frontend

This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.

A364642 a(n) is the number of iterations of psi(phi(x)) starting at x = n and terminating when psi(phi(x)) = x (n is counted), -1 otherwise.

This page as a plain text file.
%I A364642 #34 Sep 13 2023 23:13:42
%S A364642 1,2,1,2,3,2,4,3,4,3,5,3,5,4,4,4,5,4,6,4,5,5,6,4,6,5,6,5,6,4,7,5,6,5,
%T A364642 6,5,7,6,6,5,7,5,7,6,6,6,7,5,7,6,6,6,7,6,7,6,7,6,7,5,8,7,7,6,7,6,8,6,
%U A364642 7,6,8,6,8,7,7,7,8,6,8,6,8,7,8,6,7,7,7,7
%N A364642 a(n) is the number of iterations of psi(phi(x)) starting at x = n and terminating when psi(phi(x)) = x (n is counted), -1 otherwise.
%C A364642 Here phi is Euler totient function and psi is the Dedekind psi function.
%C A364642 psi(phi(1)) = 1, and psi(phi(3)) = 3. Each term of the sequence is evaluated by calling psi(phi(x)) (beginning at x = n) repeatedly until psi(phi(x)) = x. a(n) is then the number of iterations.
%C A364642 If n = 3*2^k then a(n) = k+1. Hence for any x, should x = 3*2^k then the process terminates.
%C A364642 If n = 2^k for k >= 2 then a(n) = k.
%C A364642 See A364631 for additional comments.
%F A364642 a(2^k) = A003434(2^k) = k, k >= 2.
%e A364642 a(1) = 1 since psi(phi(1)) = 1.
%e A364642 a(2) = 2 since psi(phi(2)) = 1, and psi(phi(1)) = 1.
%e A364642 a(5) = 3 since psi(phi(5)) = 6, psi(phi(6)) = 3, and psi(phi(3)) = 3.
%e A364642 a(9) = 4 since psi(phi(9)) = 12, psi(phi(12)) = 4, psi(phi(4)) = 3, and psi(phi(3)) = 3.
%t A364642 psi[n_] := n*Times @@ (1 + 1/FactorInteger[n][[;; , 1]]); psi[1] = 1; a[n_] := -1 + Length@ FixedPointList[psi[EulerPhi[#]] &, n]; Array[a, 100] (* _Amiram Eldar_, Aug 04 2023 *)
%o A364642 (Python)
%o A364642 from sympy.ntheory.factor_ import totient
%o A364642 from sympy import isprime, primefactors, prod
%o A364642 def psi(n):
%o A364642     plist = primefactors(n)
%o A364642     return n*prod(p+1 for p in plist)//prod(plist)
%o A364642 def a(n):
%o A364642     i = 1
%o A364642     r = n
%o A364642     while (True):
%o A364642         rc = psi(totient(r))
%o A364642         if (rc == r):
%o A364642             break;
%o A364642         r = rc
%o A364642         i += 1
%o A364642     return i
%o A364642 (PARI) dpsi(n) = n * sumdivmult(n, d, issquarefree(d)/d); \\ A001615
%o A364642 a(n) = my(k=0, m); while (1, m=dpsi(eulerphi(n)); k++; if (m ==n, return(k)); n=m); \\ _Michel Marcus_, Aug 14 2023
%Y A364642 Cf. A000010, A001615, A003434, A364631.
%K A364642 nonn
%O A364642 1,2
%A A364642 _Torlach Rush_, Jul 30 2023