II-7. Somme des puissances des premiers entiers
FOIRE AUX QUESTIONS DE fr.sci.maths CHAPITRE II: DEMONSTRATIONS
II-7. Somme des puissances des premiers entiers
SUJET PRECEDENT - INDEX - SUJET SUIVANTOn cherche à connaître la valeur de la somme de i=0 à n de i^k en fonction de n. k=1. La formule est due à Léonard Euler. Il suffit d'écrire: S = 1 + 2 +...+(n-1)+ n S = n +(n-1)+...+ 2 + 1 et en additionnant: 2S=(n+1)+(n+1)+...+(n+1)+(n+1) c'est-à-dire S = n(n+1)/2 k>1. Il faut alors définir les nombres de Bernoulli, que l'on notera B_t. On peut les définir comme suit: les B_t sont les nombres de Bernoulli définis (par exemple) par la série génératrice : +oo ----- t \ x x ) B_t - = ---------- / t! exp(x) - 1 ----- t = 1 Alors on peut écrire: n k ----- ----- \ k 1 \ t (k+1-t) ) i = --- ) B_i C (n+1) / k+1 / k+1 ----- ----- i = 1 i = 0 Remarque. Les C(t, k+1) représentent les coefficients du binôme de Newton, c'est-à-dire le nombre de combinaisons de t éléments d'un ensemble à k+1 éléments. La définition des nombres de Bernoulli n'est pas tout à fait standardisée : il y traîne chez certains auteurs des facteurs (-1)^t; chez d'autres les indices t deviennent t/2, et j'en passe. Il convient donc de toujours bien regarder la définition adoptée. Avec celle-ci, on a : B_0 = 1, B_1 = -1/2, B_2 = 1/6, B_4 = -1/30, B_6 = 1/42, B_8 = -1/30, B_10 = 5/66, B_12 = -691/2730, B_14 = 7/6,... et B_3 = B_5 = B_7 = ... = 0. On connaît bien sûr des techniques de calcul rapide des nombres de Bernoulli, la plupart récurrentes. A propos de formules explicites pour calculer rapidement ces nombres, on dispose tout de même du théorème de von Staudt - Clausen qui dit que B_(2k) + la somme des 1/p, étendue aux nombre premiers p tels que (p-1) divise (2*k) est un entier. Sachant par ailleurs que, pour k > 0, B_(2k) = (-1)^(k-1) * 2 * (2*k)! * zeta(2*k)/(2*pi)^(2*k), zeta étant bien sûr la fonction de Riemann, on peut programmer le problème par un procédé hybride sans utiliser de récurrence (d'abord approcher, en multi-précision, puis obtenir la valeur rationnelle exacte grâce à von Staudt). On a alors comme valeurs pour k variant de 2 à 9 (inclus): n ----- \ 2 ) i = 1/6 n (n + 1) (2 n + 1) / ----- i = 1 n ----- \ 3 2 2 ) i = 1/4 n (n + 1) / ----- i = 1 n ----- \ 4 2 ) i = 1/30 n (2 n + 1) (n + 1) (3 n + 3 n - 1) / ----- i = 1 n ----- \ 5 2 2 2 ) i = 1/12 n (2 n + 2 n - 1) (n + 1) / ----- i = 1 n ----- \ 6 4 3 ) i = 1/42 n (2 n + 1) (n + 1) (3 n + 6 n - 3 n + 1) / ----- i = 1 n ----- \ 7 2 4 3 2 2 ) i = 1/24 n (3 n + 6 n - n - 4 n + 2) (n + 1) / ----- i = 1 n ----- \ 8 6 5 4 3 2 ) i = 1/90 n (2 n + 1)(n + 1)(5 n +15 n +5 n - 15 n -n + 9 n- 3) / ----- i = 1 n ----- \ 9 2 2 4 3 2 2 ) i = 1/20 n (n + n - 1) (2 n + 4 n - n - 3 n + 3) (n + 1) / ----- i = 1 On peut également dévélopper une approche algébrique assez détaillée. On montre dans un premier temps que cette somme S(r,n) peut s'écrire somme des H(r,i)*n^(i+1) pour i variant de 0 à r. Le fait qu'un tel polynôme existe découle de l'observation suivante: Dans le sous espace vectoriel des polynômes de degrés inférieurs ou égaux à r, les u_i (i variant de 0 à r) définis par : u_i = (1+X)^(i+1) - X^(i+1) constituent une base. Si l'on nomme (H(r,i) ; i=0,...,r) le système des coordonnées de (1+X)^r on retrouve que la somme des H(r,i)*n^(i+1) est égale à la somme des p^r pour p variant de 1 à n. L'utilisation de cette base permet de déduire que les H(r,i) sont solution du système M*v=w. où ( C(x,y) désignant le coefficient binômial bien connu...) * M est la matrice triangulaire supérieure dont le terme m(l,c) vaut C(c,l-1) pour c entre 1 et r+1 l entre 1 et c par exemple dans le cas r=2, la matrice est : C(1,0) C(2,0) C(3,0) 0 C(2,1) C(3,1) 0 0 C(3,2) * v est le vecteur dont les composantes sont H(r,0) à H(r,r) * w est le vecteur dont les composantes sont C(r,0) à C(r,r). Le pivot de Gauss donne rapidement les premières valeurs: H(r,r) = 1/(r+1) H(r,r-1) = 1/2 (constant... c'est marrant) H(r,r-2) = r/12 H(r,r-3) = 0 H(r,r-4) = C(r,3)/120 H(r,r-5) = 0 H(r,r-6) = C(r,5)/252 H(r,r-7) = 0 il y a là de quoi traiter rapidement jusqu'au cas r=7 mais... En dérivant l'égalité polynômiale somme des u_i = (1+X)^r on trouve H(r,i) = [r/(i+1)]*H(r-1,i-1) on en déduit H(r,i) = [C(r,i)/(i+1)]*H(r-i,0). Cette relation limite les recherches aux H(k,0) mais, dans le système initial, c'est celui que le pivot de Gauss donne en dernier... La relation donnant H(r,i) en fonction de H(r,i+1),...,H(r,r) devient H(r,0) = 1 - somme_des{ [C(r+1,j)/(r+1)]*H(j,0) ; j variant de 0 à r-1 } et l'on a H(0,0) =1. Cette dernière relation permet de programmer rapidement le calcul des coefficients dans les cas où r est assez grand... Par ailleurs, on a remarqué qu'il y avait pas mal de 0... En fait, on a : pour tout k : H(2k+3,0)=0. Cela peut se montrer de plusieurs façons... l'égalité polynômiale avec X=1 puis X=-1 donne : somme des H(k,k-j) ; j pair entre 0 et k est égale à somme des H(k,k-j) ; j impair entre 0 et k est égale à 1/2 en bricolant un peu avec les termes connus on trouve que somme_des { C(k+1,j)*H(j,0) ; j impair entre 3 et k } est nulle. cette égalité est valable pour tout k (au moins 3...) et l'on sait que H(3,0) est nul ... donc, de proche en proche...