This site is supported by donations to The OEIS Foundation.
User:Peter Luschny/GeneralizedBernoulliNumbers
From OeisWiki
Generalized Bernoulli
Numbers and Polynomials
KEYWORDS: Bernoulli numbers, Bernoulli polynomials, Stirling numbers, Eulerian numbers, StirlingFrobenius numbers.
Concerned with sequences: A225480, A225481, A226038, A226039, A226040, A226156, A226157
This is the third part of a series of posts which build on each other.
It is recommended that the reader skims through the previous parts first to become familiar with the definitions and notations given there which we will not repeat here.
Contents 
You can also read this post much better rendered (with MathJax) on my homepage
The generalized Bernoulli numbers
The Bernoulli numbers (Concrete Mathematics, section 6.5), which are rational numbers, can be generalized via the StirlingFrobenius subset numbers.
We will call these numbers the generalized Bernoulli numbers. In the case m = 1 we will suppress the index m and write
These are the classical Bernoulli numbers by a theorem due to Frobenius relating Eulerian numbers to Stirling set numbers . We use notations introduced by D. Knuth.


To simplify things we will use the Clausen numbers as the common denominators for all . The Clausen numbers are defined as and for as where p is prime (A027760, see also A160014).
This means, for instance, that we will not reduce to lowest terms but will write the numbers in the form (see A225480 for the numerators):
Of course it requires a proof that this is possible. For m = 1 this is the famous von StaudtClausen theorem and for m > 1 this follows from the fact that the denominator of divides . In other words we have a generalization of the von StaudtClausen theorem from "For " to "For and divides "
Clausennormalized numerators of B_{n}^{(m)} 
0  1  2  3  4  5  6  7 
m = 1, A027641  1  1  1  0  1  0  1  0 
m = 2, A225480  1  0  2  0  14  0  62  0 
m = 3  1  1  3  2  39  10  363  98 
m = 4  1  2  2  6  14  50  62  854 
m = 5  1  3  1  12  145  148  4537  3892 
m = 6  1  4  6  20  546  340  22506  12740 
m = 7  1  5  13  30  1321  670  71533  33910 
The table above can be computed with Sage and this script:
def Clausen(n): #  The Clausen numbers  if n == 0: return 1 return mul(filter(lambda s: is_prime(s), map(lambda i: i+1, divisors(n)))) def ClausenB(n): #  Clausen numbers, implementation better to read but much slower  return mul(filter(lambda p: n%(p1) == 0, primes(n+2))) @CachedFunction def EulerianNumber(n, k, m) : #  The generalized Eulerian numbers  if n == 0: return 1 if k == 0 else 0 return (m*(nk)+m1)*EulerianNumber(n1,k1,m) + (m*k+1)*EulerianNumber(n1,k,m) def gen_bernoulli_number(n, m): #  The generalized Bernoulli numbers  return add(add(EulerianNumber(n, j, m)*binomial(j, n  k) for j in (0..n))*(1)^k/(k+1) for k in (0..n)) for m in (1..7): [gen_bernoulli_number(n,m)*Clausen(n) for n in (0..7)]
Plotting abs(gen_bernoulli_number(2n, m)) for and with the vertical axis in logarithmic scale we see how fast the numbers grow in size (the graph of the classical Bernoulli numbers is red):
Note that for . Moreover the are integers. Thus is the integer sequence
This is up to sign A009843 (Egf. x/cos(x)), but more exciting these are the numerators of the lost Bernoulli numbers as I called them in my blog post on the EulerBernoulli diamond, the numerators of the Bernoulli secant numbers A009843 / A193476.
The observation that are integers is just a special case of a much more general fact: All are integers for and . Thus the Bernoulli numbers are only semirationals (pardon the expression), a fact which is of no great importance in the classical case because all classical Bernoulli numbers with an odd index are zero.
m \ n  1  2  3  4  5  6  7 
3  1  5  49  809  20317  722813  34607305 
4  3  25  427  12465  555731  35135945  2990414715 
5  6  74  1946  88434  6154786  607884394  80834386026 
6  10  170  6370  415826  41649850  5922729722  1134081384850 
7  15  335  16955  1503963  204957775  39666688711  10337889346275 
Fortunately all these sequences can be computed by a single formula:
These sequences can be found (modulo sign, offset, suppressed zeros and other OEIStypical peculiarities) in A002111, A009843, A069852, A069994, A083011, A083012, A083013, A083014. But much more is true.
With regard to the StirlingFrobenius numbers and generalized Eulerian numbers we thus have the identities
 (I) 
This proves that our choice of the parameter m in generalizing the Eulerian numbers and the Stirling numbers was a sensible choice.
But wait a minute! Aren't the polynomials the Bernoulli polynomials? Despite the appearance they are not (DLMF). Anyway, let us experiment a little bit: What happens if we replace the Bernoulli numbers in the formula by the Bernoulli polynomials?
Generalized Bernoulli polynomials
More precisely we define the generalized Bernoulli polynomials as
 (II) 
The reader may convince himself that are the classical Bernoulli polynomials. Moreover
In particular we now can state the generalized Bernoulli numbers in a form independent from the StirlingFrobenius numbers and generalized Eulerian numbers in the same way as the classical Bernoulli numbers are associated with the classical Bernoulli polynomials:
Many special cases of the formulas above could be added to the sequences in OEIS. We only give one for the purpose of illustration: the Glaisher G numbers A002111 (J.W.L. Glaisher, On a set of coefficients analogous to the Eulerian numbers) written as generalized Bernoulli numbers.
# Generalized Bernoulli polynomials def gen_bernoulli_polynomial(n, m, x): p = add(add(add(((1)^(nv)/(j+1))*binomial(n,k)*binomial(j,v)*(m*(vx))^k for v in (0..j)) for j in (0..k)) for k in (0..n)) return expand(p) # Alternative description of the generalized Bernoulli polynomials # based on standard Bernoulli polynomials def gen2_bernoulli_polynomial(n, m, x): return (1)^n*sum(binomial(n,k)*bernoulli_polynomial(x,k)*m^k for k in (0..n)) # Generalized Bernoulli numbers # based on generalized Bernoulli polynomials def gen_bernoulli_number(n, m): return gen_bernoulli_polynomial(n, m, 0)
The scaled generalized Bernoulli numbers
The definition of the StirlingFrobenius numbers
suggests a second way to generalize Bernoulli numbers in the spirit of Frobenius:
We will call these numbers the scaled generalized Bernoulli numbers. Clearly .


It is convenient to use as common denominators for all the numbers which are defined and for as where p is prime. In other words is the product over all primes p ≤ n + 1 such that p divides n + 1 or p − 1 divides n. We call this the weak Clausen condition because it relaxes the Clausen condition by logical disjunction with . This results in sequence A225481:
Of course the claim that for all n and m the denominator of divides requires a proof. Moreover we claim that the terms of all other sequences with this property are multiples of the terms of this sequence. The two claims can be conflated to
The right hand side makes sense as, for n fixed, there are in fact only a finite number of different (where the are reduced to lowest terms). The additional '2' is necessary because the Clausen condition is true for the prime even if 2 does not divide n.
C^{~}lausennormalized numerators of 
0  1  2  3  4  5  6  7 
m = 1, A226156  1  1  1  0  1  0  1  0 
m = 2, A226157  1  1  2  2  14  33  62  132 
m = 3  1  3  7  6  411  26  10767  1178 
m = 4  1  5  28  0  1546  1191  27868  24202 
m = 5  1  7  61  28  2941  5580  68123  96212 
m = 6  1  9  106  90  3426  15175  594474  172520 
m = 7  1  11  163  198  1111  31362  2119277  3682 
The table above can be computed with this Sage script:
#  The generalized scaled Bernoulli numbers  def genscale_bernoulli_number(n, m): return add(add(EulerianNumber(n, j, m)*binomial(j, n  k) for j in (0..n))/((m)^k*(k+1)) for k in (0..n)) def Clausen_tilde(n): #  The weak Clausen condition  if n == 0: return 1 return mul(filter(lambda p: ((n+1)%p == 0) or (n%(p1) == 0), primes(n+2))) def Clausen_tildeA(n): #  Same as Clausen_tilde but much faster  if n == 0: return 1 F = filter(lambda s: is_prime(s), map(lambda i: i+1, divisors(n))) G = [x for x in prime_divisors(n+1) if not x in F] return mul(F)*mul(G) def divides(a,b): return b%a == 0 #  Best to read but slow  def A225481(n): #  A225481 is the OEIS name for Clausen_tilde  return mul(filter(lambda p: divides(p,n+1) or divides(p1,n), primes(n+2))) for m in (1..7): [genscale_bernoulli_number(n,m)*Clausen_tilde(n) for n in (0..7)]
Appendix: The sequence C^{~}_{n} / C_{n}
From the definition and for
it follows that divides . The quotients are A226040:
Let's spell out the definition: is the product over all primes p such that p divides n + 1 and p − 1 does not divide n.
Analyzing this sequence one finds two other interesting sequences: The positions n such that and their complement, the positions n such that .
=/= : 5, 9, 11, 13, 14, 17, 19, 20, 21, 23, 25, 27, 29, 32, 33, 34, 35, = : 0, 1, 2, 3, 4, 6, 7, 8, 10, 12, 15, 16, 18, 22, 24, 26, 28, 30,
I expected to find well known numbertheoretic sequences and so I threw the above number slices into the searchmouth of OEIS. And indeed there were two hits. A080765 for and called "m such that divides ". And A181062 for named "prime powers minus 1". Simple relations like these were exactly what I was hoping for. But do they fit our definitions?
The answer is 'No'. According to our definition the case says that every prime which divides also divides And says at least one prime divides which does not divide
The first spoiler is . The only primes which divide 45 are 3 and 5. Since 2 and 4 divide 44 the primes 3 and 5 also divide . However 45 is not the power of a prime. So I decided to submit two new sequences: A226038 and A226039.
def is_A080765(n): return not is_prime_power(n+1) def A080765_list(n): return filter(is_A080765, (0..n)) def is_A181062(n): return is_prime_power(n+1) def A181062_list(n): return filter(is_A181062, (0..n)) def F(n): return filter(lambda p: ((n+1)%p == 0) and (n%(p1) <> 0), primes(n)) def A226040(n): return mul(F(n)) def A226038_list(n): return filter(lambda n: []==F(n), (0..n)) def A226039_list(n): return filter(lambda n: []<>F(n), (0..n))
As if things are not already confusing enough there is a third sequence closely related to and : the denominators of the Bernoulli polynomials. We will use them in the next section.
0  1  2  3  4  5  6  7  8  9  10  11  12  13  14  
A141056 Clausen  1  2  6  2  30  2  42  2  30  2  66  2  2730  2  6 
A225481 Weak Clausen  1  2  6  2  30  6  42  2  30  10  66  6  2730  14  30 
A144845 Denom Bernoulli polynomials  1  2  6  2  30  6  42  6  30  10  66  6  2730  210  30 
Appendix: The generalized Bernoulli polynomials B_{n,2}(x).
Again it is convenient to list the generalized Bernoulli polynomials not with coefficients which are reduced to lowest terms, rather with coefficients relative to the denominators of the classical Bernoulli polynomials (which are in A144845). We will call this the normalized form of the generalized Bernoulli polynomials. Thus the beautiful polynomials in the case m = 2 can be written as
1  1  
4x  2  
24x^{2} − 2  6  
16x^{3} − 4x  2  
480x^{4} − 240x^{2} + 14  30  
192x^{5} − 160x^{3} + 28x  6  
2688x^{6} − 3360x^{4 }+ 1176x^{2} − 62  42  
768x^{7} − 1344x^{5} + 784x^{3} − 124x  6 
Summary
Many generalizations of the Eulerian, the Stirling and the Bernoulli numbers have been proposed and studied in the literature. We aimed at generalizing these three kind of numbers simultaneously using the most basic relations between these numbers.
Inspired by the connection of the cardinal Bsplines (Euler splines) with the Eulerian polynomials we first generalized the Eulerian numbers to
with the special value Next we introduced generalized Stirling numbers via the connection between the classical Stirling subset numbers and the Eulerian numbers due to Frobenius.
Based on that we introduced the generalized Bernoulli numbers
They can be seen as the value of the generalized Bernoulli polynomials at the origin.