INTEGRATION NUMERIQUE

Objectif : calculer numériquement une intégrale

Soit f : [a , b] → R une fonction continue

b
L’objectif est de calculer l’intégrale suivante : I = ∫ f ( x )dx
a

On partitionne l’intervalle [a,b] en n intervalles

[xi,xi+1], i = 0,…,n-1 avec a = x0 < x1 <…< xn = b

On a alors :
n −1
I = ∫ f ( x )dx = ∑ ∫
b x i +1
f ( x )dx
a xi
i =0

1
On doit donc déterminer l’intégrale suivante :
x i +1
∫xi
f ( x )dx i =0,…,n-1

En se ramenant sur l’intervalle [-1,1] via le changement de variable

x − xi x i +1 − x i
t=2 −1 ou encore x= ( t + 1) + x i
x i +1 − x i 2

x i +1 x − xi 1
∫xi
f ( x )dx = i +1
2 ∫−1
g i ( t )dt

 x − xi 
avec g i ( t ) = f  i +1 ( t + 1) + x i  i = 0,…,n-1
 2 

2
Formule de quadrature

Soit g une fonction continue sur [-1,1]. On appelle formule de

quadrature l’application J définie par :
p
J :a ∑ ω j g (t j )
j =0

− 1 ≤ t 0 < t 1 < ... < t p ≤ 1 points d’intégration de la formule de quadrature J

ω 0 , ω1 ,..., ω p poids de la formule de quadrature J.

Les poids et les points d’intégration sont choisis de telle manière que :

p
g ( t )dt ≅ ∑ ω j g ( t j )
1
∫−1
j= 0

3
Formule de quadrature exacte

Soit Q ∈ P q, l’ensemble des polynômes de degré inférieur ou égal à q.

La formule de quadrature J est dite exacte sur Pq, si

1
∀Q ∈ Pq ∫−1
Q( t )dt = J(Q)

c’est à dire p
Q( t )dt = ∑ ω j Q( t j )
1
∀Q ∈ Pq ∫ −1
j= 0

Formule de quadrature exacte : base de Lagrange

Théorème.
Soit L0, L1,…, Lp la base de Lagrange de Pp relativement aux points

− 1 ≤ t 0 < t 1 < ... < t p ≤ 1

On a l’équivalence suivante :
p
J (g) = ∑ ωi g( t j ) exacte sur Pp ⇔ ω j = ∫ L j ( t )dt
1
j = 0, 1,…,p 4
−1
j= 0
Formule de Newton-Cotes
On appelle méthode de Newton-Cotes les méthodes précédentes pour
lesquelles on choisit des points (xi), i = 0,…,n équidistants :

x i = x 0 + i.h b−a
h=
n
x i +1 x i +1 − x i 1
Dès lors : ∫xi
f ( x )dx =
2 ∫−1
g i ( t )dt

h
≈ J (g i )
2
h p
≈ ∑ ω jg i (t j )
2 j= 0
h p h 
≈ ∑ ω j f  ( t j + 1) + x i 
2 j= 0 2 
b
On peut à présent approcher I = ∫
a
f ( x)dx par la formule composite: 5
 n −1
I h = ∫ f ( x )dx = ∑ ∫
b x i +1
f ( x )dx
a xi
i =0

h n −1 p
≈ ∑∑ ω j f  (t j + 1) + x i 
h 
2 i =0 j=0 2 

Formules de Newton classiques ainsi que le degré p de Pp.

Degré p Appellation
0 Méthode des rectangles
1 Méthode des trapèzes
2 Méthode de Simpson
3 Méthode de Simpson 3/8
4 Méthode de Bode
6
Méthode des rectangles

La méthode des rectangles est une

méthode à 1 point. On a p = 0 et
t0 = 0. La base de Lagrange
associée à P0 est

L0(t) = 1

Par suite :
1 1 A
ω0 = ∫ L 0 ( t )dt = ∫ dt = 2
−1 −1
0 0.5 1 1.5 2 2.5 3

La formule de quadrature est alors définie par : J (g) = ω 0 g ( t 0 ) = 2g (0)

On a alors la formule composite des rectangles :
 x + x i +1 
n −1
I h = h∑ f  i  7
i =0  2 
Méthode des trapèzes

La méthode des trapèzes est une

méthode à 2 points. On a

p =1 et t0 = -1, t1 = 1. La base de
Lagrange associée à P1 est :

1− t 1+ t
L 0 (t) = L1 ( t ) =
2 2
Par suite :
1 1
A
ω 0 = ∫ L 0 ( t )dt = 1 ω1 = ∫ L1 ( t )dt = 1
−1 −1
0 0.5 1 1.5 2 2.5 3

La formule de quadrature est alors définie par :

J (g) = ω 0 g ( t 0 ) + ω l g( t 1 ) = g (−1) + g (1)

h n −1
On a alors la formule composite des trapèzes : I h = ∑ (f ( x i ) + f ( x i +1 ) )
2 i=0
8
Méthode des trapèzes améliorée

La formule des trapèzes peut être améliorée en procédant comme suit,

on développe f(x) en série de Taylor au voisinage de xi :
1 1
f ( x) = f ( xi ) + ( x − xi ) f ′( xi ) +
( x − xi ) 2 f ′′( xi ) + ( x − xi ) 3 f ( 3) ( xi ) + ...
2! 3!
en intégrant de xi à xi+1, on obtient :

xi +1 h2 h3 h 4 ( 3)
∫ f ( x)dx = hf ( xi ) + f ′( xi ) + f ′′( xi ) + f ( xi ) + ...
xi 2! 3! 4!
2
h
comme : f i +1 = f ( x i +1 ) = f i + hf ′( x i ) + f ′′( x i ) + ...
2!

2 3
h h h
en multipliant par h/2 , on a : f ′( x i ) = (f i +1 − f i ) − f ′′( x i ) + ...
2! 2 4

9
l’intégrale recherchée (sur [a,b]) vaut alors :

h h 2 n −1
∑ hf ′′( xk )
xn
∫ x0
f ( x)dx ≈ [ f 0 + 2 f1 + 2 f 2 + ... + 2 f n −1 + f n ] −
2 12 k = 0
xn =b
le dernier terme (somme) est une approximation de
∫ x0 = a
f ′′( x)dx

en le remplaçant par f’(b)-f’(a), on obtient la formule des trapèzes améliorée :

h h2
f ( x)dx ≈ [ f 0 + 2 f1 + 2 f 2 + ... + 2 f n −1 + f n ] − [ f ′(b) − f ′(a)]
b
∫ a
2 12
Erreur commise en h4 , la méthode est donc d’ordre 4.
Remarque :
La première formule permet de déterminer l’erreur de la méthode
des trapèzes, soit :
h2
E (h) = − (b − a) f ′′(ξ ) , ξ ∈ [a,b].
12 10
Méthode de Simpson

La méthode de Simpson est une méthode à 3 points. On a p = 2

et t0 = -1, t1 = 0, t2 = 1. La base de Lagrange associée à P2 est

t2 − t t2 + t
L 0 (t) = L1 ( t ) = 1 − t 2
L 2 (t ) =
2 2

Par suite :
1 1 1 4 1
ω 0 = ∫ L 0 ( t )dt =
1
−1
ω1 = ∫ L1 ( t )dt = ω 2 = ∫ L 2 ( t )dt =
3 −1 3 −1 3

La formule de quadrature est alors définie par :

1 4 1
J (g ) = ω 0 g ( t 0 ) + ω1g ( t 1 ) + ω 2 g ( t 2 ) = g (−1) + g (0) + g (1)
3 3 3
On a alors la formule composite de Simpson :

h   x i + x i +1  
Ih = ∑
6 

 f ( x i ) + 4f 



+ f ( x i +1 
) 
2  11
Méthode de Romberg

Soit R0(n) la somme des aires des n trapèzes de la méthode des trapèzes.

Lorsque n → ∞ , la suite R0(n) converge en général assez lentement vers la

valeur exacte de l’intégrale

La méthode de Romberg permet d’accélérer considérablement la convergence.

formule des trapèzes améliorée

b 2 4 b−a
R 0 ( n ) = ∫ f ( x )dx + C.h + O( h ) , h = C : constante indépendante de h
a n

En doublant le nombre d’intervalles on a :

b−a
2
b h 4
R 0 ( 2n ) = ∫ f ( x )dx + C + O(h ) , h =
a 4 n
12
 Par élimination de l’inconnue C on obtient :

b 4R 0 ( 2 n ) − R 0 ( n ) 4
∫a f (x )dx = 3
+ O( h )

La précision en h2 (formule des trapèzes) passe à une précision en h4.

Construisons une nouvelle suite :

4R 0 (2n ) − R 0 (n )
R 1 ( 2n ) = , n = 1,2,4,8,...
3
Cette suite converge plus vite vers la valeur de l’intégrale que la suite R0(n)
(si C ≠ 0).

Exercice : Montrer que R1(2n) est égale à l’approximation fournie par la

méthode de Simpson pour une subdivision de [a,b] en 2n intervalles égaux.
13
Ayant :
b 4 6 b−a
R 1 (n ) = ∫ f ( x )dx + C * h + O(h ) , h =
a n

Où C* est une constante indépendante de h.

En éliminant encore la constante C*, on obtient :

b 16R 1 ( 2n ) − R 1 ( n ) 6
∫a
f ( x )dx =
15
+ O( h )

La précision passe d’un ordre en h4 à un ordre en h6, la suite

16R 1 (2n ) − R 1 (n )
R 2 ( 2n ) = , n = 2,4,8,...
15
converge donc en général plus vite vers la valeur exacte de l’intégrale.
14
Algorithme de Romberg :

Pour n = 1, 2, 4, 8, … calculer :

h b−a
R 0 (n ) = [f 0 + 2f1 + 2f 2 + ... + 2f n −1 + f1 ] , h=
2 n

Pour n = 2k-1, 2k, 2k+1, … ; m = 1, 2, 3, …

calculer :
k
4 R k −1 (2n ) − R k −1 (n )
R k ( 2n ) = k
4 −1

15
 b−a n
Quadrature de Gauss : I= ∑ ω jf (x j )
2 j=1

b+a b−a
On pose : x= + ξ
2 2

Méthode de Gauss : interpolation de la fonction par un

polynôme de degré 2n-1

les valeurs des ωj ainsi que celles de ξ sont tabulée.

16
 Racines et poids pour la méthode de Gauss - Legendre
n
± ξk ωk
2 0,5773502692 1,0000000000
3 0,0000000000 0,8888888889
0,7745966692 0,5555555556
4 0,3399810436 0,6521451549
0,8611363116 0,3478548451
5 0,0000000000 0,5688888889
0,5384693101 0,4786286705
0,9061798459 0,2369268850
6 0,2386191861 0,4679139346
0,6612093865 0,3607615730
0,9324695142 0,1713244924
7 0,0000000000 0,4179591837
0,4058451514 0,3818300505
0,7415311856 0,2797053915
0,949079123 0,1294849662
17