A195441 a(n) = denominator(Bernoulli_{n+1}(x) - Bernoulli_{n+1}).

1, 1, 2, 1, 6, 2, 6, 3, 10, 2, 6, 2, 210, 30, 6, 3, 30, 10, 210, 42, 330, 30, 30, 30, 546, 42, 14, 2, 30, 2, 462, 231, 3570, 210, 6, 2, 51870, 2730, 210, 42, 2310, 330, 2310, 210, 4830, 210, 210, 210, 6630, 1326, 858, 66, 330, 110, 798, 114, 870, 30, 30, 6

Offset: 0

Author

Peter Luschny, Sep 18 2011

Keywords

• nonn

Comments

If s(n) is the smallest number such that s(n)*(1^n + 2^n + ... + x^n) is a polynomial in x with integer coefficients then a(n)=s(n)/(n+1) (see A064538).

a(n) is squarefree, by the von Staudt-Clausen theorem on the denominators of Bernoulli numbers. - Kieren MacMillan and Jonathan Sondow, Nov 20 2015

Kellner and Sondow give a detailed analysis of this sequence and provide a simple way to compute the terms without using Bernoulli polynomials and numbers. They prove that a(n) is the product of the primes less than or equal to (n+2)/(2+(n mod 2)) such that the sum of digits of n+1 in base p is at least p. - Peter Luschny, May 14 2017

The equation a(n-1) = denominator(Bernoulli_n(x) - Bernoulli_n) = rad(n+1) has only finitely many solutions, where rad(n) = A007947(n) is the radical of n. It is conjectured that S = {3, 5, 8, 9, 11, 27, 29, 35, 59} is the full set of all such solutions. Note that (S\{8})+1 joined with {1,2} equals A094960. More precisely, the set S implies the finite sequence of A094960. See Kellner 2023. - Bernd C. Kellner, Oct 18 2023

As was observed in the example section of A318256: denominator(B_n(x)) = rad(n+1) if n is in {0, 1, 3, 5, 9, 11, 27, 29, 35, 59} = {A094960(n) - 1: 1 <= n <= 10}. - Peter Luschny, Oct 18 2023

Links

• Peter Luschny, Table of n, a(n) for n = 0..10000 (terms 0..1000 from G. C. Greubel)
• Olivier Bordellès, Florian Luca, Pieter Moree, and Igor E. Shparlinski, Denominators of Bernoulli polynomials, Mathematika 64 (2018), 519-541.
• Harald Hofstätter, Denominators of coefficients of the Baker-Campbell-Hausdorff series, arXiv:2010.03440 [math.NT], 2020. Mentions this sequence.
• Bernd C. Kellner, On a product of certain primes, J. Number Theory, 179 (2017), 126-141; arXiv:1705.04303 [math.NT], 2017.
• Bernd C. Kellner, On the finiteness of Bernoulli polynomials whose derivative has only integral coefficients, J. Integer Seq. 27 (2024), Article 24.2.8, 11 pp.; arXiv:2310.01325 [math.NT], 2023.
• Bernd C. Kellner and Jonathan Sondow, Power-Sum Denominators, Amer. Math. Monthly, 124 (2017), 695-709; arXiv:1705.03857 [math.NT], 2017.
• Bernd C. Kellner and Jonathan Sondow, The denominators of power sums of arithmetic progressions, Integers 18 (2018), #A95, 17 pp.; arXiv:1705.05331 [math.NT], 2017.
• Bernd C. Kellner and Jonathan Sondow, On Carmichael and polygonal numbers, Bernoulli polynomials, and sums of base-p digits, Integers 21 (2021), #A52, 21 pp.; arXiv:1902.10672 [math.NT], 2019.

Crossrefs

Cf. A002997, A064538, A094960, A144845, A286516, A286762, A286763, A318256, A324369, A324370, A324371.

Programs

Julia

using Nemo, Primes
function A195441(n::Int)
    n < 4 && return ZZ([1,1,2,1][n+1])
    P = primes(2, div(n+2, 2+n%2))
    prod([ZZ(p) for p in P if p <= sum(digits(n+1, base=p))])
end
println([A195441(n) for n in 0:59]) # Peter Luschny, May 14 2017

Maple

A195441 := n -> denom(bernoulli(n+1, x)-bernoulli(n+1)):
seq(A195441(i),i=0..59);
# Formula of Kellner and Sondow:
a := proc(n) local s; s := (p,n) -> add(i,i=convert(n,base,p));
select(isprime,[$2..(n+2)/(2+irem(n,2))]); mul(i,i=select(p->s(p,n+1)>=p,%)) end: seq(a(n), n=0..59); # Peter Luschny, May 14 2017

Mathematica

a[n_] := Denominator[Together[(BernoulliB[n + 1, x] - BernoulliB[n + 1])]]; Table[a[n], {n, 0, 59}] (* Jonathan Sondow, Nov 20 2015 *)
SD[n_, p_] := If[n < 1 || p < 2, 0, Plus @@ IntegerDigits[n, p]]; DD[n_] := Times @@ Select[Prime[Range[PrimePi[(n+2)/(2+Mod[n, 2])]]], SD[n+1, #] >= # &]; Table[DD[n], {n, 0, 59}] (* Bernd C. Kellner, Oct 18 2023 *)

PARI

a(n) = {my(vp = Vec(bernpol(n+1, x)-bernfrac(n+1))); lcm(vector(#vp, k, denominator(vp[k])));} \\ Michel Marcus, Feb 08 2016

Python

from math import prod
from sympy.ntheory.factor_ import primerange, digits
def A195441(n): return prod(p for p in primerange((n+2)//(2|n&1)+1) if sum(digits(n+1,p)[1:])>=p) # Chai Wah Wu, Oct 04 2023

Sage

A195441 = lambda n: mul([p for p in (2..(n+2)//(2+n%2)) if is_prime(p) and sum((n+1).digits(base=p))>=p])
print([A195441(n) for n in (0..59)]) # Peter Luschny, May 14 2017

Formula

a(n) = A064538(n)/(n+1). - Jonathan Sondow, Nov 12 2015

A001221(a(n)) = A001222(a(n)). - Kieren MacMillan and Jonathan Sondow, Nov 20 2015

a(2*n)/a(2*n+1) = A286516(n+1). - Bernd C. Kellner and Jonathan Sondow, May 24 2017

a(n) = A007947(A338025(n+1)). - Harald Hofstätter, Oct 10 2020

From Bernd C. Kellner, Oct 18 2023: (Start)

Note that the formulas here are shifted in index by 1 due to the definition of a(n) using index n+1!

a(n) = A324369(n+1) * A324370(n+1).

a(n) = A144845(n) / A324371(n+1).

a(n-1) = lcm(a(n), rad(n+1)), if n >= 3 is odd.

If n+1 is composite, then rad(n+1) divides a(n-1).

If m is a Carmichael number (A002997), then m divides both a(m-1) and a(m-2).

See papers of Kellner and Kellner & Sondow. (End)

Extensions

Definition simplified by Jonathan Sondow, Nov 20 2015