|
%I #31 Dec 16 2023 13:11:04
%S 1,1,1,2,3,1,2,8,7,1,4,12,26,15,1,2,24,56,80,31,1,6,24,124,240,242,63,
%T 1,4,48,182,624,992,728,127,1,6,48,342,1200,3124,4032,2186,255,1,4,72,
%U 448,2400,7502,15624,16256,6560,511,1,10,72,702,3840
%N Array of values of Jordan function J_k(n) read by antidiagonals (version 1).
%D L. Comtet, Advanced Combinatorics, Reidel, 1974, p. 199, #3.
%D R. Sivaramakrishnan, "The many facets of Euler's totient. II. Generalizations and analogues", Nieuw Arch. Wisk. (4) 8 (1990), no. 2, 169-187.
%H Enrique Pérez Herrero, <a href="/A059379/b059379.txt">Table of n, a(n) for n = 1..10000</a>
%F J_k(n) = sum( d divides n, d^k*mu(n/d)) - _Benoit Cloitre_ and Michael Orrison (orrison(AT)math.hmc.edu), Jun 07 2002
%e Array begins:
%e 1, 1, 2, 2, 4, 2, 6, 4, 6, 4, 10, 4, ...
%e 1, 3, 8, 12, 24, 24, 48, 48, 72, 72, ...
%e 1, 7, 26, 56, 124, 182, 342, 448, 702, ...
%e 1, 15, 80, 240, 624, 1200, 2400, 3840, ...
%p J := proc(n,k) local i,p,t1,t2; t1 := n^k; for p from 1 to n do if isprime(p) and n mod p = 0 then t1 := t1*(1-p^(-k)); fi; od; t1; end;
%p #alternative
%p A059379 := proc(n,k)
%p add(d^k*numtheory[mobius](n/d),d=numtheory[divisors](n)) ;
%p end proc:
%p seq(seq(A059379(d-k,k),k=1..d-1),d=2..12) ; # _R. J. Mathar_, Nov 23 2018
%t JordanTotient[n_,k_:1]:=DivisorSum[n,#^k*MoebiusMu[n/#]&]/;(n>0)&&IntegerQ[n];
%t A004736[n_]:=Binomial[Floor[3/2+Sqrt[2*n]],2]-n+1;
%t A002260[n_]:=n-Binomial[Floor[1/2+Sqrt[2*n]],2];
%t A059379[n_]:=JordanTotient[A004736[n],A002260[n]]; (* _Enrique Pérez Herrero_, Dec 19 2010 *)
%o (PARI)
%o jordantot(n,k)=sumdiv(n,d,d^k*moebius(n/d));
%o A002260(n)=n-binomial(floor(1/2+sqrt(2*n)),2);
%o A004736(n)=binomial(floor(3/2+sqrt(2*n)),2)-n+1;
%o A059379(n)=jordantot(A004736(n),A002260(n)); \\ _Enrique Pérez Herrero_, Jan 08 2011
%o (Python)
%o from functools import cache
%o def MoebiusTrans(a, i):
%o @cache
%o def mb(n, d = 1):
%o return d % n and -mb(d, n % d < 1) + mb(n, d + 1) or 1 // n
%o def mob(m, n): return mb(m // n) if m % n == 0 else 0
%o return sum(mob(i, d) * a(d) for d in range(1, i + 1))
%o def Jrow(n, size):
%o return [MoebiusTrans(lambda m: m ** n, k) for k in range(1, size)]
%o for n in range(1, 8): print(Jrow(n, 13))
%o # Alternatively:
%o from sympy import primefactors as prime_divisors
%o from fractions import Fraction as QQ
%o from math import prod as product
%o def J(n: int, k: int) -> int:
%o t = QQ(pow(k, n), 1)
%o s = product(1 - QQ(1, pow(p, n)) for p in prime_divisors(k))
%o return (t * s).numerator # the denominator is always 1
%o for n in range(1, 8): print([J(n, k) for k in range(1, 13)])
%o # _Peter Luschny_, Dec 16 2023
%Y See A059379 and A059380 (triangle of values of J_k(n)), A000010 (J_1), A059376 (J_3), A059377 (J_4), A059378 (J_5). Columns give A000225, A024023, A020522, A024049, A059387, etc.
%Y Main diagonal gives A067858.
%K nonn,tabl
%O 1,4
%A _N. J. A. Sloane_, Jan 28 2001
|
