La méthode des éléments finis (MÉF) est une manière numérique de résoudre certains des problèmes de physique. C'est une méthode qui permet de déterminer une solution approchée sur un domaine spatial, c'est-à-dire qui permet de calculer un champ (de scalaires, de vecteurs, de tenseurs) qui correspond à certaines équations et à certaines conditions imposées.
La méthode consiste à découper le domaine spatial en petits éléments, également appelés mailles, et à rechercher une formulation simplifiée du problème sur chaque élément, c'est-à-dire à transformer le système d'équations quelconque en un système d'équations linéaires. Chaque système d'équations linéaires peut se représenter par une matrice. Les systèmes d'équations pour tous les éléments sont ensuite rassemblés, ce qui forme une grande matrice ; la résolution de ce système global donne la solution approchée au problème.
La quantité de calculs nécessite l'emploi d'un ordinateur.
Théorie et pratique.
Les logiciels modernes utilisant la méthode des éléments finis bénéficient d'une interface graphique rendant leur utilisation relativement simple. Par ailleurs, un certain nombre de tâches sont automatisables. On peut donc quasiment lancer un calcul sur ordinateur sans connaître la méthode.
Cependant, le modèle utilisé risque d'être inadapté au problème, on aura donc un résultat très éloigné de la réalité. L'utilisateur doit avoir des connaissances suffisantes pour être en mesure de
D'autre part, les premières versions de ces logiciels utilisaient un langage de programmation (langage DMAP pour Nastran, MAPDL pour Ansys, Code_Aster). Si les fonctions principales — les plus couramment utilisées — ont été intégrées dans l'interface graphique, certaines opérations ne sont possible qu'en entrant des instructions. Dans le cas d'une utilisation avancée du logiciel, il faut donc s'intéresser au langage de programmation.
L'utilisation d'un logiciel de résolution par la méthode des éléments finis est donc faussement simple, ce qui n'est pas sans poser problème :
On peut essayer de déterminer des compétences minimales requises en fonction de l'utilisation. Définissons trois niveaux d'utilisation en mécanique (le niveau 2 inclut le niveau 1, le niveau 3 inclut le niveau 2) :
Les compétences requises recommandées sont :
Principe de la méthode.
Équation aux dérivées partielles.
Un certain nombre de problèmes physiques sont décrits par des équations aux dérivées partielles (ÉDP) sur un domaine spatial, un volume. Il s'agit d'une généralisation des équations différentielles aux fonctions de plusieurs variables. Par exemple, si l'on a une fonction de trois variables ƒ("x"1, "x"2, "x"3), l'équation suivante : 
est une équation aux dérivées partielles. Cette équation est assortie de conditions aux limites : valeurs de la fonction ou de ses dérivées partielles en certains points.
Notons que
La résolution exacte, analytique, de telles équations devient vite impossible manuellement. Par contre, si l'on découpe le domaine spatial en petites cellules, appelées « éléments finis » (ÉF), on peut résoudre simplement l'ÉDP sur chaque élément.
Démarche générale.
La méthode des éléments finis (MÉF) consiste donc à :
On retrouve la démarche générale analyse-synthèse.
D'un point de vue pratique, la mise en œuvre de la méthode comporte les étapes suivantes :
Maillage.
Le maillage "(meshing)" consiste donc à découper l'espace en petits domaines appelés éléments finis. Dans le cas général, on utilise des éléments volumiques (3D), mais on peut :
On utilise plusieurs types d'éléments finis. Dans un premier temps, retenons qu'il y a principalement deux types d'éléments pour les coques et les volumes :
Un élément fini est une maille ; il est défini par ses nœuds, c'est-à-dire les sommets pour un polygone ou un polyèdre.
La MÉF consiste à calculer les valeurs de la fonction ƒ — solution des ÉDP — aux nœuds ; on ne recherche pas la valeur en tout point de l'espace, mais uniquement en certains points.
La fonction ƒ est une fonction vectorielle de dimension "d". En chaque nœud, il faut donc déterminer "d" valeurs. On dit que le nœud a « "d" degrés de liberté ».
Interpolation des fonctions.
Donc, on s'intéresse aux valeurs de la fonction ƒ uniquement aux nœud du maillage. La valeur aux autres points est ensuite déduite par interpolation polynomiale.
Sauf cas exceptionnel, dans une interpolation, plus on s'éloigne des points connus, plus l'écart avec la valeur théorique est important. Donc, plus le maillage est grossier, plus l'erreur commise en dehors des nœuds est grande.
Considérons un élément donné dont les nœuds sont numérotés de 1 à "n" ; on note ƒ"i" la valeur de la fonction au nœud "i". Ce sont ces valeurs ƒ"i" que l'on veut déterminer par la MÉF.
Au sein de l'élément, on remplace donc la fonction par son interpolation, sous la forme
où les N"i" sont des fonction qui ne dépendent que de la forme de l'élément fini.
On utilise différents types d'éléments finis ; pour chaque type, on a un élément de référence, qui est l'élément fini dont les coordonnées des nœuds valent 0, 1 ou -1. On peut dire en quelque sorte que l'élément fini réel est l'élément fini de référence déformé et placé dans l'espace. Les fonctions N"i" sont tabulées pour les éléments de référence ; on applique ensuite un changement de repère.
On écrit souvent cette opération d'interpolation sous la forme d'un produit matriciel :
Formulations intégrales.
Le système d'ÉDP peut s'écrire
les fonctions R"i" étant appelées « résidus ». Pour une fonction ψ"i" arbitraire, dite « fonction de pondération », on a ψ"i"×R"i"(ƒ) = 0 et donc "a fortiori"
Ceci constitue la forme intégrale forte. Les fonctions sont notées W par analogie avec les travaux virtuels. On choisit les fonctions de pondération de manière judicieuse afin de faciliter la résolution ; ce ne sont pas nécessairement des fonctions connues, elles peuvent contenir la fonction ƒ.
Les ÉDP sont en général d'ordre 2. Les fonctions résidus contiennent donc des termes du type ∂2ƒ/∂"x"2. À une dimension, on peut faire une intégration par parties :
Cette opération permet :
Pour des intégrales de volume, l'équivalent de l'intégration par partie est le théorème de Green :
où
Les conditions limites apparaissent dans l'intégrale sur la surface.
Lorsque l'on applique le théorème de Green aux fonctions W"i", on obtient la « forme intégrale faible ».
Matrice de rigidité.
La méthode consiste donc à résoudre le système d'équations
les fonctions W"i" étant les formulations intégrales faibles des fonctions résidu pondérées (les fonctions résidu étant les fonctions originales du système d'ÉDP).
Sur un élément fini "j" donné, on ne travaille pas sur la fonction ƒ elle-même mais sur son interpolation linéaire
La formulation intégrale faible fait de ce fait intervenir les intégrales des fonctions N"k" qui sont connues. Ainsi, la fonction W"i" est une fonction faisant intervenir (pour simplifier, on reprend la forme à une dimension) :
On obtient ainsi une équation matricielle de la forme :
où
Si l'élément fini a "n" degrés de liberté, alors l'équation s'écrit avec les composantes :
La méthode des éléments finis a d'abord été développée pour la résistance des matériaux. Dans ce domaine, F"i" est une composante d'une force qui s'exerce sur un nœud "i" et U"i" est une composante du déplacement d'un nœud. On a un rapport évident avec la loi des ressorts
où "k" est l raideur du ressort. De fait, la matrice K est appelée « matrice de rigidité. »
Comme indiqué précédemment, on ne connaît pas les conditions limites pour un élément fini quelconque. On connaît donc la matrice K, qui est caractéristique de l'ÉF, mais rien sur F ni sur U.
Lors de la synthèse, les vecteurs F locaux sont mis bout-à-bout pour former un vecteur global Fg, contenant toutes les composantes de tous les ÉF. De même, les vecteurs locaux U viennent former un grand vecteur global Ug, et les matrices de rigidité locales viennent former une matrice globale Kg. Cette opération est appelée « expansion-assemblage ».
On a ainsi un système linéaire de la forme
qui contient les conditions aux limites (valeurs imposées de Fg ou Ug).
Habituellement, dans une équation matricielle comme celle-ci, soit l'on connaît entièrement Ug et l'on cherche Fg, soit l'on connaît entièrement Fg et l'on cherche Ug. Ici, c'est différent, puisque certaines composantes de Ug sont connues et d'autres pas ; et certaines composantes de Fg sont connues et d'autres sont inconnues.
Mais on se retrouve globalement à résoudre un système linéaire.
Exemples d'application.
Les quatre principaux domaines pour lesquels on a recours à la MÉF sont :
Résistance des matériaux.
Historiquement, c'est le domaine pour lequel la MÉF a été inventée.
Le système d'équations aux dérivées partielles est fourni par les équations de la statique (les efforts extérieurs s'exerçant sur un élément de matière s'annulent) :
où
Les conditions aux limites sont de deux ordres :
Le déplacement d'un point sous l'effet de la charge est représenté par le vecteur
On impose donc un champ de déplacement formula_18 sur une surface Su. On relie ce déplacement aux contraintes (grandeurs de l'ÉDP) par le passage aux déformations :
et l'utilisation de la loi de Hooke :
où H est une matrice de coefficients d'élasticité, caractéristique du matériau ; σ et ε sont des vecteurs dont les composantes sont respectivement les contraintes et les déformations.
On impose des efforts surfaciques
formula_22
(des pressions) sur une surface Sf. Si le vecteur normal unitaire à la surface est notée
alors l'équilibre de la surface impose
Notons bien que
C'est ce que calcule l'ordinateur. Ce champ est transformé en champ de déformation ε (avec la définition fournie ci-dessus), puis le champ des contraintes σ est calculé avec la loi de Hooke.
Transfert de chaleur.
C'est le domaine pour lequel les équations sont les plus simples.
L'équation aux dérivées partielles est l'équation de la chaleur (équation de Fourier) :
où
En régime stationaire, on a
Les conditions limite consistent à imposer :
Mécanique des fluides.
La simulation numérique des écoulement est en générale appelée MFN (mécanique des fluides numérique) ou bien CFD "(computational fluids dynamics)".
On utilise parfois la méthode des éléments finis, mais la plupart du temps, c'est la méthode des volumes finis qui est utilisée : les équations différentielles ne sont pas résolues de la même manière. Dans le cas des logiciels munis d'une interface graphique, on voit peu de différence avec un logiciel d'éléments finis, mais dans tous les cas, le traitement de la géométrie du système est très différent : on ne peut pas simplifier les détails de la même manière, et la stratégie de maillage obéit à des impératifs différents.
Algorithmes de résolution.
La méthode étant posée, les logiciels ont recours à des méthodes numériques — calculs approchés — pour effectuer les calculs. Nous présentons ici quelques méthodes utilisées classiquement.
Un modèle d'éléments finis peut contenir dix mille, cent mille voire parfois plus d'un million d'éléments. Le nombre de calcul est très important, la résolution est donc extrêmement gourmande en ressources et en temps : plusieurs heures voire plusieurs jours de calcul, et peut nécessiter le lancement d'un calcul parallélisé sur un superordinateur. L'utilisation d'algorithmes de calcul performants n'est donc pas inhérent à la méthode des éléments finis, mais est toutefois indispensable pour des raisons de performance.
Intégration.
De manière intuitive, pour résoudre une équation différentielle, il faut inverser le processus de dérivation c'est-à-dire calculer une intégrale. Dans le MÉF, les intégrales sur les éléments finis sont calculées de manière approchée, numérique. On utilise en général la :
Les points "xi" sont appelés « points de Gauss ». Les coefficients ω"i" sont des constantes. Les valeurs des "xi" et des ω"i" dépendent du nombre de points "n" que l'on considère. Plus ce nombre est élevé, plus le calcul est précis, mais plus il est long.
Cette méthode se généralise aux fonctions de plusieurs variables. La figure ci-contre montre un exemple de carré avec quatre points de Gauss (il y a plusieurs possibilité) ; les coefficients de pondération ω"i" sont ici tous égaux à 1. Soit ƒ("x", "y" ) une fonction définie sur [-1 ; 1]2 ; on peut approcher son intégrale par
ce calcul est exact si ƒ est un polynôme de dergé inférieur ou égal à trois, c'est-à-dire comportant des termes de la forme "xi"×"yj" avec "i" + "j" ≤ 3.
Ceci se généralise pour d'autres formes (en particulier triangle) et pour des intégrations sur des volumes (cubes, tétraèdres, prismes droits, …).
Bien évidemment, il n'y a aucune raison pour qu'un élément fini ait ses frontières aux coordonnées 1 ou -1. On fait donc un changement de repère pour relier l'élément fini du maillage à un élément fini de référence. Ce changement d'espace de coordonnées dans le contexte dérivée/intégrale fait classiquement intervenir les matrices jacobiennes.
Rappelons que l'on ne travaille pas directement avec la fonction ƒ mais avec son interpolation linéaire ∑"j" = 1"m" N"j"·ƒ"j" (si l'élément fini a "m" nœuds). Le calcul de l'intégrale — à une dimension "x" — devient donc
Rappelons que tout comme les "xi" et ω"i", les N"j" sont tabulées pour les éléments de référence. On peut faire ressortir les inconnues :
les coefficients α"j" étant connus.
Notes et références.
Introduction < > Rappels de mécanique