COST (GBP)

1.00

0.00

Maths › Discrete › Number Theory ›

Bernoulli B

Calculates array of Bernoulli numbers using an infinite series

Controller: CodeCogs

Contents

Interface
Bernoulli B
Page Comments

Interface

C++

Bernoulli B

voidbernoulli_B(	int	`iMax`
	double*	`dB`	)

Bernoulli numbers $\inline B_n$ are particular values of Bernoulli polynomials $\inline B_n(x): B_n=B_n(0)$ . Polynomials $\inline B_n(x)$ are expandable in the Fourier series:

$B_{2n}(x)=2{{(-1)}^{n+1}}(2n)!\sum _{r=1}^{\infty }{{(2 \pi r)}^{-2n}}cos(2 \pi r x)$

(2)

$\begin{array}{l} {B_{2n+1}}(x)= 2{{(-1)}^{n+1}}(2n+1)!\sum _{r=1}^{\infty }{{(2 \pi r)}^{-2n-1}} sin(2 {\pi rx}) \end{array}$

(3)

Substitution $\inline x=0$ gives series expansion for Bernoulli numbers:

$\begin{array}{l} {B_{2n}}=2{{(-1)}^{n+1}}(2n)! \sum _{r=1}^{\infty }{{(2 \pi r)}^{-2n}} \\ n=1,2,3,... \\ {B_{2n+1}}=0,\\ n=1,2,3,...\\ \end{array}$

(4)

However, series like this is not quite suitable for numerical evaluation, especially at small n, in view of slow convergence and uneasy accuracy estimation. One needs to use another relation [2,23.10.21]:

${B_{n(}} \big ( \frac{1}{2}\big)=-(1-{2^{1-n}}){B_n}$

(5)

Then Fourier series expansion results in the following:

$\begin{array}{l} {B_{2n}}=\frac{{{(-1)}^{n+1}}2(2n)!}{(1-{2^{1-2n}}){{(2\pi )}^{2n}}}\sum _{r=1}^{\infty }\frac{{{(-1)}^{r-1}}}{{r^{2n}}}, \\ n=1,2,3,....\\ \end{array}$

(6)

Alternating series converges faster and gives possibility for simple accuracy estimation: error arising due to series truncation is less than first term neglected. Then N, the amount of series terms to obtain accuracy $\inline (\epsilon )$ may be estimated as follows:

$N={{\epsilon }^{-1/2n}}$

(7)

In the worst case n=2 reasonable accuracy $\inline (\epsilon ={10}^{-12} )$ results in $\inline ({10}^6)$ series terms to sum. To significantly reduce a mount of summands and speed up convergence it is suitable to use Aitken transformation of partial series sums. Let

${S_n}=\sum _{k=1}^{n}{f_k}$

(8)

then sequence of new sums

${T_n}= \frac{S_{n+1}S_{n-1} - S^2_n}{ S_{n+1} + S_{n-1} - 2{S_n}}$

(9)

may converge faster compared with $\inline (S_n)$ . To avoid subtractions of near values and losses of significant digits expression should be rewritten as follows:

${T_n}= {S_n}-\frac{f_{n+1}f_n}{f_{n+1}-f_n}$

(10)

or, using

${f_n}= {(-1)^n}/{g_n}$

(11)

${T_n}= S_n + \frac{{(-1)}^n}{g_{n+1}} + g_n$

(12)

Array dimension should be iMax+1 or greater.

References:

Higher Transcendental Functions, vol.1, (1.13) by H.Bateman and A.Erdelyi (Bateman Manuscript Project), 1953
M.Abramowitz and I.A.Stegun, Handbook of Mathematical Functions, 1964 chapt.23
Yu.Luke, Mathematical functions and their approximations, 1975 chapt.14.2

Example 1

#include <stdio.h>
#include <codecogs/maths/discrete/number_theory/bernoulli_b.h>
 
#define MAX_INDEX 16
 
int main()
{
  double dBernoulli[MAX_INDEX+1];
 
  printf( "%8s%2c%20s\n", " ", 'n', "Bn" );
  printf( "%8s", " " );
  for(int i = 0; i < 22; i++ )
    printf( "%c", '-' );
  printf( "\n" );
 
  Maths::NumberTheory::bernoulli_B( MAX_INDEX, dBernoulli );
 
  printf( "%10d%20.12f\n", 0, dBernoulli[0] );
  printf( "%10d%20.12f\n", 1, dBernoulli[1] );
  for(int i = 2; i <= MAX_INDEX; i += 2 )
    printf( "%10d%20.12f\n", i, dBernoulli[i] );
  return 0;
}

Output:

n                  Bn
----------------------
0      1.000000000000
1     -0.500000000000
2      0.166666666667
4     -0.033333333333
6      0.023809523810
8     -0.033333333333
10     0.075757575758
12    -0.253113553114
14     1.166666666667
16    -7.092156862745

Parameters

iMax	input maximal index requested
dB	output pointer on the array of numbers declared in the calling module.

Authors

Anatoly Prognimack (Mar 16, 2005)
Developed tested with Borland C++ 3.1 for DOS and Microsoft Visual C++ 5.0, 6.0
Updated by Will Bateman (March 2005)

Source Code

Source code is available when you agree to a GP Licence or buy a Commercial Licence.

Not a member, then Register with CodeCogs. Already a Member, then Login.