Pour la théorie des nombres de Bernoulli, voir entre autres ici: http://fr.wikipedia.org/wiki/Nombre_de_Bernoulli

Je m'y suis intéressé parce que ces nombres interviennent dans certains calculs, en particulier dans le calcul de la tangente par série entière.

Je suis parti de l'équation récurrente suivante:

• <m> sum{k=0}{n}c_n_1^{k}}*{B_k}} = 0 </m>

• avec <m>B_0 = 1</m>

• et <m>B_i = 0</m> si i est un nombre impair supérieur à 1

Dans cette formule:

• <m>{C_{n+1}}^k</m> est le nombre de combinaison de n+1 objets pris k à k (voir analyse combinatoire sur ce site)

• <m>{B_k}</m> est le nombre de Bernoulli pour la valeur k

Mais en quoi cette équation est-elle “récurrente”?

Pour le savoir, développons pour n=3:

• <m>c_4^0}*{B_0}} ~+~ c_4^1}*{B_1}} ~+~ c_4^2}*{B_2}} ~+~ c_4^3}*{B_3}} ~=~ 0</m>

Il faudrait donc que l'on puisse calculer <m>B_3</m>.

Mais là, on voit bien le problème: on connait <m>{B_0}=1</m> par définition, mais pour calculer <m>B_3</m>, il faut calculer avant: <m>B_1</m> et <m>B_2</m>.

Comment? En appliquant la formule pour n=1 et n=2 et, à chaque fois, en calculant le dernier B.

Voilà le code qui fait cela (Python version 2.6):

#!/usr/bin/python
# -*- coding: utf-8 -*-
from __future__ import division
 
from decimal import *
 
##############################################################################
def combin(n, k):
    """Nombre de combinaisons de n objets pris k a k"""
    if k > n//2:
        k = n-k
    x = 1
    y = 1
    i = n-k+1
    while i <= n:
        x = (x*i)//y
        y += 1
        i += 1
    return x
 
##############################################################################
def bernoulli(n):
    """calcul du nombre de Bernoulli Bn """ 
    prec = getcontext().prec
    getcontext().prec = 50   
    B = [Decimal("1")]
    for m in xrange(1,n+1):
        if m>1 and m%2!=0:
            B.append(Decimal("0"))
        else:
            c = Decimal("0")
            for k in xrange(0,m):
                c += Decimal(str(combin(m+1,k)))*B[k]
            B.append(-c/Decimal(str(combin(m+1,m))))
    getcontext().prec = prec 
    return float(B[n])