login
a(n) is the number of binomial coefficients C(n,k), with 0<=k<=n, that are not practical.
3

%I #22 Aug 11 2023 13:47:47

%S 0,0,2,0,4,2,6,1,2,4,6,2,8,10,14,0,6,2,8,5,10,10,20,2,8,14,14,6,16,14,

%T 30,0,2,4,6,2,6,10,16,4,10,10,24,16,14,20,36,4,8,10,16,12,30,14,30,6,

%U 14,18,34,14,34,38,62,0,2,2,6,10,12,10,24,2,8,14,14

%N a(n) is the number of binomial coefficients C(n,k), with 0<=k<=n, that are not practical.

%H Paolo Leonetti and Carlo Sanna, <a href="https://arxiv.org/abs/1905.12023">Practical numbers among the binomial coefficients</a>, arXiv:1905.12023 [math.NT], 2019.

%p A334082 := proc(n)

%p local a,k ;

%p a := 0 ;

%p for k from 0 to n do

%p if not isA005153(binomial(n,k)) then # reuses code of A005153

%p a := a+1 ;

%p end if;

%p end do:

%p a ;

%p end proc:

%p seq(A334082(n),n=1..100) ; # _R. J. Mathar_, Jul 07 2023

%t PracticalQ[n_] := Module[{f, p, e, prod = 1, ok = True}, If[n < 1 || (n > 1 && OddQ[n]), False, If[n == 1, True, f = FactorInteger[n]; {p, e} = Transpose[f]; Do[If[p[[i]] > 1 + DivisorSigma[1, prod], ok = False; Break[]]; prod = prod*p[[i]]^e[[i]], {i, Length[p]}]; ok]]];

%t a[n_] := Select[Table[Binomial[n, k], {k, 0, n}], !PracticalQ[#]&] // Length;

%t Table[a[n], {n, 1, 100}] (* _Jean-François Alcover_, Aug 11 2023, after _T. D. Noe_ in A005153 *)

%o (PARI) a(n) = sum(k=0, n, !is_A005153(binomial(n,k)));

%o (Python)

%o from math import comb

%o from sympy import factorint

%o def A334082(n):

%o m = 0

%o for k in range(1,n):

%o c = comb(n,k)

%o l = (~c & c-1).bit_length()

%o if l>0:

%o P = (1<<l+1)-1

%o for p, e in factorint(c>>l).items():

%o if p > 1+P:

%o break

%o P *= (p**(e+1)-1)//(p-1)

%o else:

%o continue

%o m += 1

%o else:

%o m += 1

%o return m # _Chai Wah Wu_, Jul 05 2023

%Y Cf. A005153 (practical numbers), A007318 (binomial coefficients).

%Y Cf. A334083, A334084.

%K nonn

%O 1,3

%A _Michel Marcus_, Apr 14 2020