Méthode tau

En analyse numérique la méthode tau est, avec la méthode de Galerkine ou celle de collocation, l'une des méthodes utilisée dans le cadre des méthodes spectrales pour la résolution des équations différentielles ou aux dérivées partielles^[1],[2],[3].

On considère l'équation différentielle suivante établissant la fonction ${\displaystyle u(x)}$ définie sur l'intervalle ${\displaystyle x\in [0,1]}$ :

${\displaystyle u'(x)-u(x)=0\,,\quad u(0)=1}$

Sa solution exacte est ${\displaystyle u=\exp(x)}$. On cherche une solution approchée en écrivant l'équation :

${\displaystyle u'(x)-u(x)+\tau P(x)=0\,,\quad u(0)=1}$

où ${\displaystyle P(x)}$ est un polynôme de degré N. Il existe un polynôme tel que la solution soit également un polynôme de degré N et ${\displaystyle \tau }$ un paramètre que l'on espère petit. Par exemple, si on prend ${\displaystyle P(x)=x^{2}}$ la solution s'écrit :

${\displaystyle u(x)={\frac {x^{2}}{2}}+x+1}$ pour ${\displaystyle \tau =1/2}$.

Le résidu par rapport au développement de la solution exacte est en ${\displaystyle x^{2}/2}$. On obtient ainsi une approximation heuristique raisonnable au regard de la simplicité de la méthode. En augmentant l'ordre à la valeur N, l'erreur est en ${\displaystyle {\frac {x^{N}}{N!}}}$ et décroît donc avec N.

On peut généraliser cette méthode à un système d'équations aux dérivées partielles linéaire^[4].

Cette méthode a été généralisée, en particulier par Steven Orszag. On considère l'équation différentielle :

${\displaystyle Lu(x)-f(x)=0\,,\quad \alpha <x<\beta }$

où ${\displaystyle L}$ est un opérateur linéaire.

Les conditions aux bords sont données par les opérateurs ${\displaystyle B_{-}}$ et ${\displaystyle B_{+}}$ (linéaires, de type Dirichlet, Neumann ou Robin) :

${\displaystyle {\begin{array}{rcl}B_{-}u(\alpha )&=&g_{-}\\[0.3em]B_{+}u(\beta )&=&g_{+}\end{array}}}$

On recherche la solution sous la forme d'une série de fonctions de base :

${\displaystyle u_{N}(x)=\sum _{n=0}^{N}a_{n}\phi _{n}(x)}$

Ces fonctions de base sont orthogonales avec un poids donné w positif, adapté au choix de la fonction de base :

${\displaystyle (\phi _{i},\phi _{j})_{w}=\int \phi _{i}{\overline {\phi _{j}}}wdx=c_{i}\delta _{ij}}$

où ${\displaystyle \delta }$ est le symbole de Kronecker.

Elles vérifient les conditions aux bords :

${\displaystyle {\begin{array}{rcl}B_{-}u_{N}(\alpha )&=&g_{-}\\[0.3em]B_{+}u_{N}(\beta )&=&g_{+}\end{array}}}$

On appelle résidu la quantité :

${\displaystyle R_{N}(x)=Lu_{N}(x)-f(x)}$

Le problème est donc de résoudre le système :

${\displaystyle (R_{N},\phi _{m})_{w}=0\,,\quad \forall m=0,N-2}$

soit :

${\displaystyle \sum _{n=0}^{N}a_{n}(L\phi _{n},\phi _{m})=(f,\phi _{m})_{w}\,,\quad \forall m=0,N-2}$

Avec les deux conditions aux bords, on obtient ainsi un système linéaire de ${\displaystyle N+1}$ équations à résoudre.

Les polynômes couramment utilisés sont ceux de Tchebychev, de Legendre ou de Laguerre^[2].

On considère l'équation différentielle suivante^[5] portant sur la fonction ${\displaystyle u(x)}$ définie sur l'intervalle ${\displaystyle x\in [-1,1]}$ :

${\displaystyle u''(x)=1\,,\quad u(-1)=0\,,\quad u'(1)=0}$

Sa solution exacte est ${\displaystyle u={\frac {x^{2}}{2}}-x-{\frac {3}{2}}}$. On cherche une solution approchée sous forme d'un développement en polynômes de Tchebychev :

${\displaystyle {\begin{array}{rcl}u(x)&=&\sum _{n=0}^{N}a_{n}T_{n}(x)\\[0.3em]{\frac {du}{dx}}&=&\sum _{n=0}^{N-1}a_{n}^{(1)}T_{n}(x)\\[0.3em]{\frac {d^{2}u}{dx^{2}}}&=&\sum _{n=0}^{N-2}a_{n}^{(2)}T_{n}(x)\end{array}}}$

En substituant dans l'équation différentielle et les conditions aux limites; il vient :

${\displaystyle {\begin{array}{rcl}\sum _{n=0}^{N-2}a_{n}^{(2)}T_{n}(x)&=&1\\[0.3em]\sum _{n=0}^{N}a_{n}T_{n}(-1)&=&0\\[0.3em]\sum _{n=0}^{N-1}a_{n}^{(1)}T_{n}(1)&=&0\end{array}}}$

On utilise les relations de récurrence :

${\displaystyle \sum _{n=1}^{N}a_{n}{\frac {\mathrm {d} T}{\mathrm {d} x}}(1)=\sum _{n=1}^{N}a_{n}n^{2}}$

Les relations d'orthogonalité permettent de calculer les coefficients ${\displaystyle a_{n}^{(2)}}$ :

${\displaystyle a_{n}^{(2)}=\delta _{0n}}$

Par substitution dans les conditions aux limites et en tenant compte de ${\displaystyle T_{n}(1)=1}$ et ${\displaystyle T_{n}(-1)=(-1)^{n}}$ :

${\displaystyle \sum _{n=0}^{N}(-1)^{n}a_{n}=0}$

${\displaystyle \sum _{n=1}^{N}n^{2}a_{n}=0}$

On prend ${\displaystyle N=5}$. Dans ce cas, les termes d'ordre 2 deviennent :

${\displaystyle {\begin{array}{rclll}a_{0}^{(2)}&=&4a_{2}+32a_{4}&=&1\\[0.2em]a_{1}^{(2)}&=&24a_{3}+120a_{5}&=&0\\[0.2em]a_{2}^{(2)}&=&48a_{4}&=&0\\[0.2em]a_{3}^{(2)}&=&48a_{4}&=&0\end{array}}}$

et les conditions aux limites s'écrivent :

${\displaystyle {\begin{array}{rcl}a_{0}-a_{1}+a_{2}-a_{3}+a_{4}-a_{5}&=&0\\[0.2em]a_{1}+4a_{2}+9a_{3}+16a_{4}+25a_{5}&=&0\end{array}}}$

La solution du système matriciel des 6 équations est :

${\displaystyle a=\left(-{\frac {5}{4}},-1,{\frac {1}{4}},0,0,0\right)^{T}}$

La solution du problème est donc :

${\displaystyle {\begin{array}{rcl}u(x)&=&-{\frac {5}{4}}T_{0}(x)-T_{1}(x)+{\frac {1}{4}}T_{2}(x)\\[0.3em]&=&{\frac {x^{2}}{2}}-x-{\frac {3}{2}}\end{array}}}$

C'est la solution exacte du problème.

• Portail de l'analyse