Méthodes de quadrature de Gauss

Points utilisés pour une méthode de quadrature de Gauss.

Dans le domaine mathématique de l'analyse numérique, les méthodes de quadrature sont des approximations de la valeur numérique d'une intégrale. En général, on remplace le calcul de l'intégrale par une somme pondérée prise en un certain nombre de points du domaine d'intégration (voir calcul numérique d'une intégrale pour plus d'informations). La méthode de quadrature de Gauss, du nom de Carl Friedrich Gauss[1], est une méthode de quadrature exacte pour un polynôme de degré 2n – 1 avec n points pris sur le domaine d'intégration.

Les formules de Gauss jouent un rôle fondamental dans la méthode des éléments finis.

Principe général

On souhaite évaluer numériquement l'intégrale

I = a b f ( x ) ϖ ( x ) d x . {\displaystyle I=\int _{a}^{b}f(x)\varpi (x)\,\mathrm {d} x.}

Le domaine d'intégration (a, b) couvre plusieurs cas :

  • intervalles bornés : comme [a, b], [a, b[, etc.
  • demi-droite réelle : [a, +∞[, ] -∞, b],
  • la droite réelle tout entière : ℝ.

Les méthodes sont de la forme

I = a b f ( x ) ϖ ( x ) d x i = 1 n ω i f ( x i ) {\displaystyle I=\int _{a}^{b}f(x)\varpi (x)\,\mathrm {d} x\approx \sum _{i=1}^{n}\omega _{i}f(x_{i})}

ϖ : ( a , b ) R + {\displaystyle \varpi :(a,b)\to \mathbb {R} _{+}} est une fonction de pondération continue strictement positive, qui peut assurer l'intégrabilité de f. Les ω i {\displaystyle \omega _{i}} sont appelés les coefficients de quadrature (ou poids). Les points xi sont appelés les nœuds de la quadrature.

Théorème fondamental

Pour n donné :

  • Les n nœuds xi sont réels, distincts, uniques et sont les racines du polynôme unitaire de degré n, orthogonal au sous-espace des polynômes de degré n-1 pour le produit scalaire h , g = a b h ( x ) g ( x ) ϖ ( x ) d x {\displaystyle \left\langle h,g\right\rangle =\int _{a}^{b}h(x)g(x)\varpi (x)\,\mathrm {d} x} [2].
  • Pour tout i, ω i {\displaystyle \omega _{i}} est égal à a b l i ( x ) ϖ ( x ) d x {\displaystyle \int _{a}^{b}l_{i}(x)\varpi (x)\,\mathrm {d} x} , où li est le polynôme interpolateur de Lagrange de degré n – 1, prenant la valeur 1 en xi, et 0 en les xk pour k différent de i.
  • Pour les polynômes de degré inférieur ou égal à 2n – 1, il y a égalité entre a b f ( x ) ϖ ( x ) d x {\displaystyle \int _{a}^{b}f(x)\varpi (x)\,\mathrm {d} x} et i = 1 n ω i f ( x i ) {\displaystyle \sum _{i=1}^{n}\omega _{i}f(x_{i})} .
Démonstration

On se borne à montrer que, si la formule de la quadrature de Gauss est valide, alors les nœuds et les poids sont fixés de manière unique, aux valeurs annoncées. Supposons donc que nous avons des nœuds xi et des poids ω i , i = 1 , , n {\displaystyle \omega _{i},i=1,\cdots ,n} tels que

Pour tout polynôme P de degré inférieur ou égal à 2n-1,
a b P ( x ) ϖ ( x ) d x = P , 1 = i = 1 n ω i P ( x i ) {\displaystyle \int _{a}^{b}P(x)\varpi (x)\,\mathrm {d} x\;=\;\left\langle P,1\right\rangle \;=\;\sum _{i=1}^{n}\omega _{i}P(x_{i})} .

Posons alors P n = Π i = 1 n ( x x i ) {\displaystyle P_{n}=\Pi _{i=1}^{n}(x-x_{i})} . Il vient

0 k < n , P n , x k = x k P n , 1 = i = 1 n ω i x i k P n ( x i ) = 0 {\displaystyle \forall 0\leq k<n,\;\left\langle P_{n},x^{k}\right\rangle =\left\langle x^{k}P_{n},1\right\rangle =\sum _{i=1}^{n}\omega _{i}x_{i}^{k}P_{n}(x_{i})=0} . Ainsi Pn est orthogonal à tout polynôme de degré inférieur au sien. Il est par conséquent le seul polynôme unitaire de la direction orthogonale aux polynômes de degré i = 1 , , n 1 {\displaystyle i=1,\cdots ,n-1} , dans l'espace vectoriel des polynômes de degré 1 , , n {\displaystyle 1,\cdots ,n} . L'article sur les polynômes orthogonaux montre alors que Pn a n racines distinctes, ce sont les xi.

L'égalité a b f ( x ) ϖ ( x ) d x = i = 1 n ω i f ( x i ) {\displaystyle \int _{a}^{b}f(x)\varpi (x)\,\mathrm {d} x=\sum _{i=1}^{n}\omega _{i}f(x_{i})} appliquée aux li donne les valeurs des ω i {\displaystyle \omega _{i}} annoncées.

Le domaine d'intégration et la fonction de pondération déterminent le type de la quadrature de Gauss. Le tableau suivant résume les situations les plus communes.

Principales configurations de quadrature de Gauss
Domaine d'intégration (a,b) Fonction de pondération ϖ ( x ) {\displaystyle \varpi (x)} Famille de polynômes orthogonaux
[–1, 1] 1 Legendre
]–1, 1[ ( 1 x ) α ( 1 + x ) β   ,   α , β > 1 {\displaystyle (1-x)^{\alpha }(1+x)^{\beta }\ ,\ \alpha ,\beta >-1} Jacobi
]–1, 1[ 1 1 x 2 {\displaystyle {\frac {1}{\sqrt {1-x^{2}}}}} Tchebychev (premier type)
]–1, 1[ 1 x 2 {\displaystyle {\sqrt {1-x^{2}}}} Tchebychev (second type)
+ e x {\displaystyle \mathrm {e} ^{-x}} Laguerre
e x 2 {\displaystyle \mathrm {e} ^{-x^{2}}} Hermite

Une fois le type de quadrature choisi, la formule à n points s'écrit :

I ( f ) = i = 1 n ω i f ( x i ) . {\displaystyle I(f)=\sum _{i=1}^{n}\omega _{i}f(x_{i}).}

On définit l'erreur comme E ( f ) = | I I ( f ) | {\displaystyle E(f)=|I-I(f)|} . Le degré d'exactitude d'une formule de quadrature est le degré le plus élevé de la famille des polynômes annulant E(f). On a le résultat suivant : une formule à n points admet un degré d'exactitude de 2n – 1.

Cas particulier pour un intervalle fermé

Le domaine d'intégration [a, b] peut être changé (au moyen d'un changement de variable) en [–1, 1] avant d'appliquer les méthodes de quadrature de Gauss. Le changement se déroule ainsi :

a b f ( x ) d x = b a 2 1 1 f ( b a 2 x + a + b 2 ) d x . {\displaystyle \int _{a}^{b}f(x)\,\mathrm {d} x={\frac {b-a}{2}}\int _{-1}^{1}f\left({\frac {b-a}{2}}x+{\frac {a+b}{2}}\right)\,\mathrm {d} x.}

L'approximation de la valeur de l'intégrale devient :

b a 2 i = 1 n ω i f ( b a 2 x i + a + b 2 )   . {\displaystyle {\frac {b-a}{2}}\sum _{i=1}^{n}\omega _{i}f\left({\frac {b-a}{2}}x_{i}+{\frac {a+b}{2}}\right)~.}

où les xi sont ici relatifs à l'intervalle [–1, 1].

Méthodes courantes

Méthode de Gauss-Legendre

Pour le problème d'intégration le plus classique, on utilise la méthode de Gauss-Legendre[3]. Il s'agit d'intégrer la fonction f sur le segment [–1, 1]. Les n nœuds sont les racines du n-ième polynôme de Legendre, Pn(x), et les coefficients sont donnés par l'une ou l'autre égalité :

ω i = 2 ( n + 1 ) P n ( x i ) P n + 1 ( x i ) = 2 ( 1 x i 2 ) P n ( x i ) 2 . {\displaystyle \omega _{i}={\frac {-2}{(n+1)P'_{n}(x_{i})P_{n+1}(x_{i})}}={\frac {2}{(1-x_{i}^{2})P'_{n}(x_{i})^{2}}}.}

On peut aussi remarquer que la somme des coefficients est égale à 2. Le tableau suivant donne l'ensemble des informations pour réaliser le calcul approché de I pour les formules à un, deux et trois points.

Nombre de points, n Poids ( ω i ) {\displaystyle (\omega _{i})} Points ( x i ) {\displaystyle (x_{i})} Polynôme de Legendre
1 2 0 x
2 1, 1 1/3, 1/3 (3x2 – 1)/2
3 5/9, 8/9, 5/9 3/5, 0, 3/5 (5x3 – 3x)/2

Exemple

On cherche à déterminer 1 1 ( x + 1 ) 2 d x {\displaystyle \int _{-1}^{1}(x+1)^{2}\mathrm {d} x} . On cherche à intégrer un polynôme de degré 2, 2 points suffisent pour obtenir la valeur exacte.

1 1 ( x + 1 ) 2 d x = 1 ( 1 3 + 1 ) 2 + 1 ( 1 3 + 1 ) 2 = 8 3 . {\displaystyle \int _{-1}^{1}(x+1)^{2}\mathrm {d} x=1\left({\frac {1}{\sqrt {3}}}+1\right)^{2}+1\left(-{\frac {1}{\sqrt {3}}}+1\right)^{2}={\frac {8}{3}}.}

On peut facilement vérifier ce résultat car dans cet exemple, on connaît une primitive de (x + 1)2 :

1 1 ( x + 1 ) 2 d x = [ ( x + 1 ) 3 3 ] 1 1 = 8 3 . {\displaystyle \int _{-1}^{1}(x+1)^{2}\mathrm {d} x=\left[{\frac {(x+1)^{3}}{3}}\right]_{-1}^{1}={\frac {8}{3}}.}

Cet exemple ne représente pas un cas pratique. En règle générale, on n'obtient jamais un résultat exact et bien entendu, on n'applique pas ces méthodes pour les fonctions dont on connaît une primitive.

Méthode de Gauss-Tchebychev

Cette formule est associée au poids ϖ ( x ) = ( 1 x 2 ) 1 / 2 {\displaystyle \varpi (x)=(1-x^{2})^{-1/2}} sur ]–1, 1[. Pour une formule à n points[4], les nœuds sont

x i = cos ( ( 2 i 1 ) π 2 n ) {\displaystyle x_{i}=\cos \left({\frac {(2i-1)\pi }{2n}}\right)}

et les coefficients :

w i = π n . {\displaystyle w_{i}={\frac {\pi }{n}}.}

Méthode de Gauss-Laguerre

Cette formule est associée au poids ϖ ( x ) = e x {\displaystyle \varpi (x)={\rm {e}}^{-x}} sur ]0, +∞[. Les n nœuds sont les n racines du n-ième polynôme de Laguerre Ln, et les coefficients sont

w i = 1 ( n + 1 ) L n ( x i ) L n + 1 ( x i ) . {\displaystyle w_{i}={\frac {1}{(n+1)L'_{n}(x_{i})L_{n+1}(x_{i})}}.}

Les coefficients et les nœuds ne peuvent être calculés analytiquement que pour n petit[5]. Par exemple, pour n = 2 :

n x i {\displaystyle x_{i}\,} w i {\displaystyle w_{i}\,}
2 2 ± 2 {\displaystyle 2\pm {\sqrt {2}}} 2 2 4 {\displaystyle {\frac {2\mp {\sqrt {2}}}{4}}}

Maintenant, pour intégrer une fonction f sur ℝ+, il faut remarquer que

0 + f ( x ) d x = 0 + f ( x ) ϖ ( x ) ϖ ( x ) d x . {\displaystyle \int _{0}^{+\infty }f(x)\,\mathrm {d} x=\int _{0}^{+\infty }{\frac {f(x)}{\varpi (x)}}\varpi (x)\,\mathrm {d} x.}

Il reste alors à appliquer la formule de quadrature à la fonction g ( x ) = f ( x ) / ϖ ( x ) . {\displaystyle g(x)=f(x)/\varpi (x).}

Méthode de Gauss-Hermite

Sur ℝ, la formule de Gauss-Hermite est caractérisée par la pondération ϖ ( x ) = e x 2 {\displaystyle \varpi (x)={\rm {e}}^{-x^{2}}} . Pour une formule à n points[6], les xi sont calculés comme les n racines du n-ième polynôme d'Hermite Hn ; quant aux pondérations, elles sont obtenues à partir de

w i = 2 n + 1 n ! π [ H n ( x i ) ] 2 . {\displaystyle w_{i}={\frac {2^{n+1}n!{\sqrt {\pi }}}{[H_{n}'(x_{i})]^{2}}}.}

Concernant l'intégration de f sur ℝ, il suffit d'appliquer la formule de quadrature à la fonction f ( x ) e x 2 . {\displaystyle f(x){\rm {e}}^{x^{2}}.}

Autres méthodes de quadrature de Gauss

Méthodes de Gauss-Lobatto

Cette section est vide, insuffisamment détaillée ou incomplète. Votre aide est la bienvenue ! Comment faire ?

Pour le cas des méthodes de quadrature de Gauss-Lobatto (en) sur l'intervalle [a, b], on impose parmi les r + 1 points de quadrature les deux extrémités de l'intervalle :

a = x 1 < x 2 < < x r + 1 = b . {\displaystyle a=x_{1}<x_{2}<\cdots <x_{r+1}=b.}

Pour un ordre de quadrature r, les points intérieurs de quadrature deviennent alors les zéros du polynôme dérivé du r-1e polynôme orthogonal :

i = 2 , , r   ,   P r 1 ( x i ) = 0. {\displaystyle \forall i=2,\ldots ,r\ ,\ P_{r-1}'(x_{i})=0.}

Méthodes de Gauss-Radau

Cette section est vide, insuffisamment détaillée ou incomplète. Votre aide est la bienvenue ! Comment faire ?

Pour le cas des méthodes de quadrature de Gauss-Radau sur l'intervalle [a, b], on impose parmi les r+1 points de quadrature une des extrémités :

a = x 1 < x 2 < < x r + 1 . {\displaystyle a=x_{1}<x_{2}<\cdots <x_{r+1}.}

Par symétrie, on peut également fixer b comme point.

Pour un ordre de quadrature r, les points intérieurs de quadrature deviennent alors les zéros du polynôme :

P r 1 ( x ) + P r ( x ) x a . {\displaystyle {\frac {P_{r-1}(x)+P_{r}(x)}{x-a}}.}

Calcul des points et poids de quadrature

Pour obtenir les points et poids de quadrature pour un ordre élevé, on consultera avec profit l'ouvrage d'Abramowitz et Stegun[7].

Notes et références

  1. Gauss a publié les principes de cette méthode dans Methodus nova integralium valores per approximationem inveniendi, Gœttingue, Heinrich Dietrich, .
  2. Jean-Pierre Demailly, Analyse numérique et équations différentielles, Presses universitaires de Grenoble, , 309 p. (ISBN 2-7061-0421-X), p. 74
  3. (en) Eric W. Weisstein, « Legendre-Gauss quadrature », sur MathWorld.
  4. (en) Eric W. Weisstein, « Chebyshev-Gauss quadrature », sur MathWorld.
  5. (en) Eric W. Weisstein, « Laguerre-Gauss quadrature », sur MathWorld.
  6. (en) Eric W. Weisstein, « Hermite-Gauss quadrature », sur MathWorld.
  7. (en) Abramowitz et Stegun, Handbook of Mathematical Functions (lire en ligne), p. 875 et suivantes

Voir aussi

v · m
Recherche de zéro
Transformations de matrice
Résolutions de systèmes
Intégration numérique
Équations différentielles
Interpolation numérique
v · m
Algèbre
Analyse
Théorie des nombres
Statistiques
Géométrie
Physique
  • icône décorative Portail de l'analyse