|
%I #16 Feb 24 2020 15:38:04
%S 1,1,1,1,1,3,1,1,1,5,1,3,1,7,5,1,1,3,1,5,7,11,1,3,1,13,1,7,1,15,1,1,
%T 11,17,35,3,1,19,13,5,1,21,1,11,1,23,1,3,1,5,17,13,1,3,55,7,19,29,1,
%U 15,1,31,7,1,13,33,1,17,23,35,1,3,1,37,5,19,77,39
%N a(n) = product{ p prime such that p divides n + 1 and p - 1 does not divide n }.
%H Charles R Greathouse IV, <a href="/A226040/b226040.txt">Table of n, a(n) for n = 0..10000</a>
%H Peter Luschny, <a href="http://www.luschny.de/math/euler/GeneralizedBernoulliNumbers.html">Generalized Bernoulli numbers</a>.
%F a(n) = A225481(n) / A141056(n).
%e a(41) = 21 = 3*7 = product({2,3,7} setminus {2}).
%p s:= (p, n) -> ((n+1) mod p = 0) and (n mod (p-1) <> 0);
%p A226040 := n -> mul(z, z = select(p->s(p,n), select('isprime', [$2..n])));
%p seq(A226040(n), n=0..77);
%t a[n_] := Times @@ Select[ FactorInteger[n+1][[All, 1]], !Divisible[n, #-1] &]; a[0] = 1; Table[a[n], {n, 0, 77}] (* _Jean-François Alcover_, Jun 27 2013, after Maple *)
%o (Sage)
%o def A226040(n):
%o F = filter(lambda p: ((n+1) % p == 0) and (n % (p-1)), primes(n))
%o return mul(F)
%o [A226040(n) for n in (0..77)]
%o (PARI) a(n)=my(f=factor(n+1)[,1],s=1);prod(i=1,#f,if(n%(f[i]-1),f[i],1)) \\ _Charles R Greathouse IV_, Jun 27 2013
%Y Cf. A160014, A226038, A226039, A225481, A141056.
%K nonn
%O 0,6
%A _Peter Luschny_, May 26 2013
|
