Objectifs du chapitre
À la fin de ce chapitre, vous serez capable de:
- décrire un maillage unidimensionnel par sa table des coordonnées et sa table de connectivité, et dire comment la numérotation des nœuds façonne la matrice de rigidité;
- écrire les fonctions de forme linéaires d'un élément, sur l'élément physique comme sur l'élément de référence , et les relier aux fonctions chapeau globales;
- démontrer l'expression de la matrice de rigidité élémentaire et du , puis démontrer que l' par la table de connectivité reproduit exactement les formes et ;
- imposer les conditions de Dirichlet par élimination et par pénalisation, comparer les deux techniques, introduire les charges de Neumann et récupérer les réactions d'appui à partir du système complet;
- calculer déformations et contraintes élément par élément, lire la table de convergence du problème modèle et expliquer pourquoi les valeurs nodales y sont exactes;
- traiter une discontinuité — saut de section, interface entre deux matériaux, force ponctuelle intérieure —, et mesurer ce que coûte un maillage qui ne place pas de nœud dessus.
Le chapitre 1 a montré l'idée de la méthode sur deux éléments et une boucle d'assemblage de quatre lignes; les chapitres 2 et 3 ont fourni la formulation faible et la méthode de Galerkin, avec ses deux garanties: l'orthogonalité de l'erreur et la meilleure approximation en énergie. Ce chapitre construit l'élément fini linéaire proprement, de la fonction de forme au post-traitement, et le programme de bout en bout. Tout ce qui viendra ensuite — treillis, poutres, triangles, quadrilatères — suivra exactement le même schéma, avec des éléments plus riches.
Le maillage et la table de connectivité
Nœuds, éléments, numérotation
Reprenons la barre de longueur du chapitre 1, de rigidité axiale , chargée par une force répartie et éventuellement par des forces concentrées. Mailler , c'est choisir des nœuds et des éléments (définition 1.5). En dimension 1, un élément est un segment entre deux nœuds, et l'on pourrait croire qu'il n'y a rien à décrire: les nœuds se suivent, les éléments aussi. C'est vrai pour un maillage fait à la main, et c'était l'hypothèse cachée de la fonction elements_finis(n) du chapitre 1, dont la ligne noeuds = [e, e + 1] disait que l'élément e relie les nœuds e et e + 1.
Un logiciel ne fait pas cette hypothèse. Son mailleur numérote les nœuds dans l'ordre où il les crée — d'abord les extrémités et les points remarquables de la géométrie, ensuite les nœuds intermédiaires —, et un modèle réel mélange des barres, des poutres et des surfaces dont les nœuds n'ont aucun ordre naturel. Le programme a donc besoin de deux tables, qui sont toute la description du maillage.
La connectivité est l'unique lien entre la numérotation locale — les nœuds 1 et 2 d'un élément, ceux de ses fonctions de forme — et la numérotation globale, celle des lignes et des colonnes de la matrice . Tout l'assemblage tient dans ce lien.
Deux conventions de numérotation coexistent dans ce cours, et ce chapitre emploie les deux. Pour le problème modèle et les problèmes scalaires, les nœuds sont et les inconnues , comme aux chapitres 1 à 3. Pour les , comme dans tous les chapitres de structures qui suivent, les nœuds et les degrés de liberté sont numérotés : la barre étagée de ce chapitre a les nœuds 1, 2 et 3. Dans le code Python, enfin, les listes commencent à 0: le nœud d'une barre est l'indice , et un commentaire le rappelle dans chaque listing.
La numérotation façonne la matrice
L'élément ne contribue qu'aux quatre termes , , et . Deux nœuds qui n'appartiennent à aucun élément commun ne sont donc : . En dimension 1 avec une numérotation qui suit la barre, chaque nœud n'a que deux voisins, et est . Mais si le mailleur numérote les nœuds dans un autre ordre, les termes non nuls s'écartent de la diagonale. On appelle le plus grand écart entre deux nœuds couplés; elle vaut 1 pour une barre numérotée dans l'ordre, et le coût de l'élimination de Gauss sur une matrice bande de taille et de demi-largeur est de l'ordre de opérations au lieu de (, chapitre 3). Une renumérotation qui réduit la bande — l'algorithme de Cuthill–McKee en est l'exemple classique — fait partie de tous les logiciels; en dimension 1, elle consiste simplement à numéroter dans l'ordre, et l'exercice 4.1 le fait à la main.
Une barre est maillée en trois éléments; la table de connectivité est: élément 1 = (1, 3), élément 2 = (3, 4), élément 3 = (4, 2). Lequel de ces termes de la matrice assemblée (4 × 4) est nécessairement nul?
Les fonctions de forme linéaires
Sur l'élément physique
Sur un élément de longueur , une fonction affine est déterminée par ses deux valeurs aux extrémités. On l'écrit
Les deux fonctions de forme et sont affines, valent respectivement et aux deux nœuds de l'élément, et leur somme vaut 1 en tout point: une translation d'ensemble est reproduite exactement par . Leurs dérivées sont constantes:
La déformation est donc constante par élément: c'est l'allongement relatif du segment, ce qu'un ingénieur aurait écrit sans aucune théorie.
Le lien avec les fonctions chapeau du chapitre 1 est immédiat (figure 4.1). La fonction chapeau du nœud vaut 1 en , 0 aux autres nœuds, et elle est affine sur chaque élément. Sur un élément qui contient le nœud , sa restriction est donc la fonction de forme de l'élément attachée à ce nœud: si est le premier nœud de , s'il est le second. Sur un élément qui ne contient pas , elle est nulle. La vue globale (, définie sur tout le domaine) et la vue locale (, , définies sur un élément) décrivent le même espace ; la première sert à raisonner, la seconde à calculer.
Sur l'élément de référence
Écrire et élément par élément oblige à transporter les coordonnées de chaque élément dans chaque formule. Il est plus économique de tout écrire une seule fois, sur un élément fixe, et de ramener chaque élément physique à celui-là par un changement de variable affine.
Deux remarques justifient cette construction. D'abord, la transformation (4.4) utilise les mêmes fonctions , pour décrire la géométrie que pour décrire l'inconnue: c'est l'idée isoparamétrique, triviale ici, qui deviendra essentielle pour les éléments courbes et les quadrilatères distordus (chapitres 7 et 10). Ensuite, le facteur de (4.4) est le jacobien de la transformation, noté . Par la règle de dérivation en chaîne,
ce qui redonne , comme (4.2). Toute intégrale sur un élément se ramène à une intégrale sur :
C'est sur cette formule que travaillent les quadratures de Gauss du chapitre 7: un programme d'éléments finis n'intègre jamais sur un élément physique, il intègre sur l'élément de référence et multiplie par le jacobien.
Matrice de rigidité et vecteur des charges élémentaires
La formulation faible, découpée en éléments
Le chapitre 2 a établi la formulation faible de la barre encastrée en , de rigidité , chargée par et par une force à son extrémité libre: trouver tel que
où est l'espace des déplacements de nuls à l'encastrement. Le chapitre 3 a remplacé par l'espace des fonctions continues affines par morceaux, et montré que la solution de Galerkin est caractérisée par
Les intégrales de (4.7) sont des sommes d'intégrales sur les éléments: et , avec et . Sur un élément, et ne dépendent que des deux valeurs nodales de l'élément, rangées dans les vecteurs élémentaires et . Par (4.1),
La matrice , de taille , est la matrice de rigidité élémentaire; le vecteur est le vecteur des charges élémentaire. Il ne reste qu'à les calculer.
Démonstration. La rigidité. Par (4.2), le produit vaut si et si , et il est constant sur l'élément. Il sort de l'intégrale (4.9):
avec le signe sur la diagonale et ailleurs. C'est (4.10), avec si la rigidité est constante.
La charge linéaire. Sur l'élément de référence, : une fonction affine est égale à son interpolée. Par (4.6),
et il suffit des trois intégrales
D'où , et de même . Pour , on trouve .
La force concentrée. Le travail de dans un déplacement virtuel est , qui est avec le vecteur annoncé.
La matrice (4.10) est celle d'un ressort de raideur , et (4.10) dit ce que la résistance des matériaux sait depuis toujours: un segment de barre se comporte comme un ressort, et deux nœuds qui s'écartent de subissent les forces . Elle est , ses lignes somment à zéro — une translation d'ensemble n'y produit aucune force —, et elle est : . Un élément isolé n'a pas d'appui, et sa translation est un mouvement de corps rigide. Le chapitre 5 retrouvera ce fait pour la structure entière.
Le mot cohérent dans (4.11) mérite une remarque. On pourrait être tenté de répartir une charge répartie «à vue», moitié sur chaque nœud, quelle que soit sa variation. Le vecteur cohérent fait mieux: il est le seul qui donne exactement le travail de la charge dans tout déplacement virtuel de , et c'est cette propriété, et elle seule, qui permet à la méthode de Galerkin de garder ses garanties. Il conserve en particulier la résultante et le moment de la charge, comme le montre l'exemple suivant.
Un élément de longueur m porte une charge axiale répartie qui croît linéairement de kN/m au premier nœud à kN/m au second. Quelle force nodale cohérente revient au second nœud, en kN?
L'assemblage
La démonstration
Il reste à vérifier que la somme des contributions élémentaires est bien le système de Galerkin (4.8). Pour l'écrire sans indices, on introduit pour chaque élément une matrice de localisation.
La matrice n'est qu'une écriture de la ligne connect[e] de la table de connectivité; aucun programme ne la construit, mais elle rend la démonstration lisible.
Démonstration. Sur l'élément , ne dépend que de ses deux valeurs nodales, , et de même . Par (4.9),
Le même calcul donne , et une force au nœud fournit le travail , qui est la composante de multipliée par . La symétrie vient de celle de chaque : . Enfin, le produit est la matrice nulle partout sauf aux quatre positions , , , , où elle porte les quatre termes de : la somme de (4.12) place donc chaque à la position (nœud global de , nœud global de ).
Le théorème dit deux choses. La première est que l'assemblage n'est pas une recette d'ingénieur posée à côté de la théorie: c'est exactement la méthode de Galerkin, écrite dans la base des fonctions chapeau. La seconde est l'algorithme: on ne forme jamais , on ajoute les quatre termes de aux quatre positions que désigne la connectivité. La figure 4.2 résume les quatre étapes d'un calcul sur un maillage de trois éléments.
L'algorithme en Python
La fonction elements_finis(n) du chapitre 1 contenait déjà la boucle d'assemblage. Il suffit d'y remplacer la connectivité implicite noeuds = [e, e + 1] par la lecture de la table, et la rigidité unité par celle de l'élément:
def assembler(coord, connect, EA):
"""Matrice de rigidite complete d'une barre d'elements lineaires.
coord[i]: abscisse du noeud i; connect[e]: les deux noeuds de l'element e."""
n = len(coord)
K = [[0.0] * n for _ in range(n)]
for e in range(len(connect)):
noeuds = connect[e] # lu dans la table de connectivite
h = coord[noeuds[1]] - coord[noeuds[0]]
C'est tout. Cette fonction ne suppose rien sur l'ordre des nœuds, ni sur l'uniformité du maillage, ni sur le nombre d'éléments; elle assemblerait aussi bien une barre de 10 000 éléments que la barre de deux éléments de l'exemple 4.2. Pour un maillage de taille industrielle, on remplacerait la liste de listes par un stockage creux qui ne garde que les termes non nuls, mais la boucle resterait la même: c'est le cœur de tout logiciel d'éléments finis.
Pour résoudre, nous reprenons la fonction resoudre(K, F) du chapitre 1, sans changement:
def resoudre(K, F):
"""Resout K d = F par elimination de Gauss avec pivot partiel."""
n = len(F)
M = [ligne[:] for ligne in K]
c = list(F)
for k in range(n - 1):
p = max(range(k, n), key=lambda i: abs(M[i][k]))
M[k], M[p] = M[p], M[k]
c[k], c[p]
Un mailleur a numéroté les nœuds de la barre étagée dans l'ordre où il les a créés: le nœud 1 à l'encastrement, le nœud 2 à l'extrémité libre (x = 2 m), le nœud 3 sur la marche (x = 1,2 m). La table de connectivité est donc (1, 3), (3, 2), soit [[0, 2], [2, 1]] en Python, où le nœud k a l'indice k - 1. La fonction ci-dessous ignore cette table: elle relie toujours l'élément e aux nœuds e et e + 1, ce qui donne ici un élément de longueur négative. Corrigez-la pour qu'elle lise la connectivité. La matrice affichée doit coupler le nœud 3 (indice 2) à ses deux voisins, et pas les nœuds 1 et 2 (indices 0 et 1).
Les conditions aux limites
La matrice assemblée sur tous les nœuds est singulière: ses lignes somment à zéro, puisque chaque a cette propriété, et le vecteur est dans son noyau. C'est le mouvement de corps rigide de la barre entière, que rien n'empêche encore. Les conditions aux limites sont ce qui le bloque, et elles se traitent de deux façons très différentes selon leur nature, comme le chapitre 2 l'a montré.
Les conditions de Neumann: rien à faire, ou presque
Une force imposée à l'extrémité libre apparaît dans sous la forme : c'est un travail virtuel. Dans le système discret, elle s'ajoute simplement à la composante de du nœud situé en , comme toute force concentrée appliquée à un nœud. Une extrémité libre et non chargée — la condition de Neumann homogène — ne demande rien du tout: on n'écrit aucune équation pour elle, et la solution discrète la vérifie approximativement d'elle-même. C'est le sens de l'expression «condition naturelle».
Les conditions de Dirichlet par élimination
Partageons les nœuds en deux groupes: les nœuds libres (indice ), dont le déplacement est inconnu, et les nœuds prescrits (indice ), dont le déplacement est imposé — nul pour un encastrement, non nul pour un tassement d'appui ou un essai à déplacement imposé. En réordonnant les lignes et les colonnes, le système complet s'écrit par blocs
où est le vecteur des réactions d'appui, inconnues, que les appuis exercent sur la barre aux nœuds prescrits. Pourquoi la réaction apparaît-elle là? Parce que le système de Galerkin (4.8) n'est écrit que pour les déplacements virtuels , nuls aux appuis: les équations des nœuds libres sont les seules que la formulation faible impose. Les équations des nœuds prescrits, qu'on obtient en assemblant sur tous les nœuds, contiennent en plus la force inconnue que l'appui doit exercer pour imposer le déplacement.
La première ligne de blocs ne contient que des inconnues du premier groupe:
C'est le système réduit. Pour un encastrement, et il suffit de supprimer les lignes et les colonnes des nœuds prescrits, ce que faisait la ligne Kr = [ligne[1:n] for ligne in K[1:n]] du chapitre 1. Pour un déplacement imposé non nul, les colonnes supprimées ne disparaissent pas: elles passent au second membre, multipliées par le déplacement imposé. La matrice est symétrique définie positive dès que les appuis bloquent les mouvements de corps rigide (chapitre 3), et le système réduit a une solution unique.
Les réactions d'appui
La seconde ligne de blocs de (4.14) n'a pas servi à calculer les déplacements. Elle sert maintenant: une fois connu, elle donne directement les réactions,
Autrement dit: on multiplie les lignes complètes des nœuds prescrits par le vecteur complet des déplacements, et l'on retranche les charges qui avaient été assemblées en ces nœuds. Ce dernier terme est celui qu'on oublie. Une charge répartie dépose sur un nœud d'appui une force nodale cohérente, pour une charge constante, qui est transmise directement à l'appui sans déformer la barre; la réaction doit l'équilibrer, et (4.16) le fait grâce au terme . Le problème guidé 4.1 et la question 4.6 le mettent en évidence.
Deux contrôles accompagnent toujours le calcul des réactions. Le premier est l'équilibre global: la somme des réactions et des charges appliquées doit être nulle, , puisque la somme des lignes de est nulle. Le second est un calcul à la main de l'ordre de grandeur: pour une barre isostatique, les réactions ne dépendent que de la statique et se calculent sans aucune matrice.
Les conditions de Dirichlet par pénalisation
L'élimination oblige à renuméroter les inconnues, ce qui est fastidieux dans un grand programme où l'on aimerait garder la structure de telle que l'assemblage l'a produite. La pénalisation évite cette renumérotation.
L'interprétation mécanique est limpide: on attache le nœud à un appui situé en par un ressort très raide de raideur . L'équation devient : le ressort exerce la force , qui est la réaction. Plus est grand, plus est proche de , mais l'égalité n'est jamais exacte: l'erreur sur vaut la réaction divisée par , et elle se propage aux autres déplacements dans un rapport de l'ordre de , où est la raideur du reste de la structure.
Le programme suivant compare les deux techniques sur la barre étagée de l'exemple 4.2 (qui en donne la solution par élimination, d), en pénalisant le nœud encastré avec :
u3_ref = d[2] # solution par elimination
for p in [2, 4, 6, 8, 12]:
alpha = 10.0**p * K[0][0]
Kp = [ligne[:] for ligne in K]
Kp[0][0] += alpha # penalite sur le noeud 1: u_1 = 0
2 u1 = 5.7e-06 m u3 = 1.085079 mm ecart 5.3e-03 R1 = -60.000 kN
4 u1 = 5.7e-08 m u3 = 1.079422 mm ecart 5.3e-05 R1 = -60.000 kN
6 u1 = 5.7e-10 m u3 = 1.079366 mm ecart 5.3e-07 R1 = -60.000 kN
8 u1 = 5.7e-12 m u3 = 1.079365 mm ecart 5.3e-09 R1 = -60.000 kN
12 u1 = 5.7e-16 m u3 = 1.079365 mm ecart 5.3e-13 R1 = -60.000 kN
L'écart relatif sur est divisé par 100 chaque fois que est multiplié par 100: il vaut exactement , parce que le ressort de pénalité, sollicité par la réaction de 60 kN, s'allonge de et que tout le reste de la barre suit. La réaction, elle, est juste dès : c'est la force dans le ressort, et elle ne dépend que de l'équilibre, que la pénalisation respecte exactement.
Où s'arrêter? La pénalisation dégrade le conditionnement de la matrice (Analyse numérique, chapitre 4): le rapport entre la plus grande et la plus petite raideur du système croît comme , et en double précision, un de fois les autres termes noie ces termes dans l'arrondi de . Sur ce petit système, l'élimination de Gauss avec pivot reste insensible même à ; sur un grand modèle résolu par une méthode itérative, dont la vitesse de convergence dépend du conditionnement, elle ne l'est pas. Le compromis usuel consiste à prendre de l'ordre de à fois le plus grand terme diagonal: six à huit chiffres de précision sur la condition imposée, et une marge de huit à dix chiffres avant l'arrondi.
| élimination | pénalisation | |
|---|---|---|
| condition imposée | exactement | à près |
| taille du système | réduite aux nœuds libres | inchangée |
| renumérotation | nécessaire | aucune |
| conditionnement | celui de | dégradé, croît avec |
| réaction | ligne mise de côté, (4.16) |
Une troisième technique, très répandue dans les programmes d'enseignement, remplace la ligne de par la ligne de l'identité et par , puis corrige les autres lignes en y transportant au second membre. Elle est exacte comme l'élimination et garde la taille du système comme la pénalisation; elle n'est que l'élimination écrite sans supprimer les lignes.
On impose l'encastrement d'une barre par pénalisation avec , et l'on obtient un déplacement en bout de barre trop grand d'environ en valeur relative. Que faut-il en conclure?
La fonction penaliser(K, F, i, valeur, alpha) doit imposer d_i = valeur par pénalisation, selon (4.17). Elle ajoute bien alpha au terme diagonal, mais oublie le second membre: elle impose donc toujours d_i = 0, quelle que soit la valeur demandée. Corrigez-la. Le programme simule un essai à déplacement imposé sur la barre étagée: encastrement au nœud 1, allongement de 1 mm imposé au nœud 3. Il doit afficher u2 = 0.4286 mm et N = 45.00 kN.
Le post-traitement: déformations et contraintes
Une déformation constante par élément
Une fois les déplacements nodaux connus, le post-traitement calcule les grandeurs dont l'ingénieur a besoin. Pour l'élément linéaire, (4.2) donne directement
toutes trois constantes sur l'élément et discontinues d'un élément à l'autre. Le déplacement est continu, sa dérivée ne l'est pas: c'est le prix d'une approximation dans seulement. Or la solution exacte, elle, a un effort normal continu partout où il n'y a pas de force concentrée. Le saut de entre deux éléments voisins est donc une erreur visible, et c'est l'indicateur le plus simple de la qualité d'un maillage: le chapitre 11 en tirera une méthode d'estimation d'erreur.
Pourquoi les valeurs nodales sont exactes
Le chapitre 1 a constaté que, pour le problème modèle avec la charge intégrée exactement, les valeurs nodales des éléments linéaires sont exactes à l'arrondi près, quel que soit , et il a promis d'en donner la raison. La voici; elle est courte, et elle explique aussi tout ce que la section suivante mesurera sur les discontinuités.
Démonstration. Fixons un nœud intérieur et définissons la fonction , continue, nulle en et en , par sa dérivée
la constante étant choisie pour que . C'est le déplacement de la barre sous une force unité au point : la du point . Pour toute fonction de nulle aux deux extrémités,
Appliquons ceci à l'erreur , qui est bien nulle aux extrémités. Si est constant sur chaque élément, ou est constant sur chaque élément — le seul changement de régime, en , est un nœud —, donc est continue et affine par morceaux sur le maillage: . L'orthogonalité de Galerkin (chapitre 3), pour tout , s'applique avec et donne
La démonstration ne s'appuie sur rien d'autre que l'orthogonalité de Galerkin et la forme de . L'exercice 3.5 du chapitre 3 établit la même exactitude nodale par une autre voie, l'interpolation: pour le problème modèle, la solution de Galerkin y coïncide avec l'interpolant de ; la fonction de Green dit pourquoi, et jusqu'où cela reste vrai. Elle s'étend sans changement à une extrémité libre (exercice 4.5), et elle vaut pour toute charge pourvu que soit calculé exactement, y compris une force concentrée entre deux nœuds. Elle dit aussi exactement quand l'exactitude nodale se perd:
- si le second membre est approché — une charge nodale au lieu de —, est la solution de Galerkin d'un autre problème, et l'on retrouve les erreurs nodales des différences finies (chapitre 1, explorateur 1.1);
- si change de valeur d'un élément, y change de pente, n'appartient plus à , et les nœuds ne sont plus exacts: c'est ce que mesurera l'exemple 4.3;
La table de convergence du problème modèle
Le problème modèle avec , de solution , traité par éléments linéaires sur un maillage uniforme de éléments avec la charge intégrée exactement, donne la table suivante, fixée pour tout le cours. L'erreur nodale est celle du programme du chapitre 1; les deux normes d'erreur sont calculées élément par élément par une quadrature de Gauss à 12 points, assez précise pour que les cinq chiffres affichés n'en dépendent pas.
| erreur nodale max | rapport |
|---|
Trois lectures. Les nœuds sont exacts, au niveau de l'arrondi, ce que le théorème 4.3 explique; l'arrondi croît légèrement avec , comme le conditionnement de (chapitre 3). L'erreur en énergie, qui mesure l'écart des pentes , est divisée par un rapport qui tend vers 2 quand double: elle est d' en . Les rapports calculés, 1,939 puis 1,985, 1,996 et 1,999, approchent 2 par en dessous; ne les arrondissez pas, ils disent que l'on entre dans le régime asymptotique. , qui mesure l'écart des valeurs, est d': ses rapports tendent vers 4. Le chapitre 11 démontrera ces deux ordres; ce chapitre se contente de les mesurer.
L'ordre 1 de l'énergie est celui des déformations, donc des contraintes: une contrainte constante par élément approche une contrainte variable comme une fonction en escalier approche une courbe. Mais toutes les positions dans l'élément ne se valent pas. Comparons, sur le même problème, l'erreur de la pente au milieu des éléments et aux nœuds (où l'on prend la valeur constante de l'élément voisin):
| erreur max au milieu | rapport | erreur max aux nœuds | rapport | |
|---|---|---|---|---|
| 2 | — | — | ||
| 4 | 2,991 |
Au milieu de l'élément, l'erreur sur la dérivée est d'ordre 2; aux nœuds, elle n'est que d'ordre 1, et plus de 120 fois plus grande pour . La raison est simple: les valeurs nodales étant exactes, la pente de l'élément est le taux d'accroissement de la solution exacte, qui est égal à en un point intérieur (théorème des accroissements finis) et, pour une fonction régulière, à (c'est la différence centrée de l', chapitre 8). Ce phénomène — un point de l'élément où la dérivée est plus précise qu'ailleurs — s'appelle la ; le chapitre 11 s'en servira pour estimer l'erreur.
La fonction post_traitement calcule l'effort normal de chaque élément et les réactions aux nœuds d'appui. Sa réaction est fausse dès qu'une charge a été assemblée sur le nœud d'appui, ce qui arrive avec toute charge répartie: elle oublie le terme -F_p de (4.16). Corrigez-la. Le programme reprend la tige suspendue de l'exemple 1.1 (60 m, 500 mm², son poids propre et 10 kN en bas), maillée en trois éléments; la réaction à l'accrochage doit valoir -(P + qL) = -12 310,3 N.
Les discontinuités
La barre étagée a été résolue exactement par deux éléments parce qu'un nœud était placé sur la marche. Cette section montre ce qui se passe quand ce n'est pas le cas, pour les trois sortes de discontinuités qu'un modèle de barre rencontre: un saut de section, une interface entre deux matériaux et une force concentrée à l'intérieur du domaine. Chacune produit une solution exacte dont la dérivée saute, et chacune pose la même question: le maillage peut-il représenter ce saut?
Un saut de section
La règle de l'exemple 4.2 — un nœud sur la marche — n'est pas toujours respectée. Un mailleur automatique qui découpe la barre en éléments égaux ignore la géométrie. Que vaut alors l'élément qui chevauche la marche? On lui donne la rigidité moyenne du théorème 4.1, qui est ce que donne l'intégration exacte de , et l'on répartit la force de 20 kN de la marche sur ses deux nœuds par ses fonctions de forme. Le programme suivant résout la barre étagée sur éléments égaux:
def barre_etagee(n, xs=1.2, L=2.0):
"""Barre etagee sur n elements egaux; renvoie u_3 en mm."""
coord = [L * i / n for i in range(n + 1)]
connect = [[e, e + 1] for e in range(n)]
EA, F, place = [], [0.0] * (n +
1 u3 = 1.03175 mm erreur -4.8e-02 mm
2 u3 = 1.05820 mm erreur -2.1e-02 mm
3 u3 = 1.06996 mm erreur -9.4e-03 mm
4 u3 = 1.06576 mm erreur -1.4e-02 mm
5 u3 = 1.07937 mm erreur -6.7e-16 mm
6 u3 = 1.07143 mm erreur -7.9e-03 mm
7 u3 = 1.07332 mm erreur -6.0e-03 mm
8 u3 = 1.07584 mm erreur -3.5e-03 mm
9 u3 = 1.07332 mm erreur -6.0e-03 mm
10 u3 = 1.07937 mm erreur -1.6e-15 mm
La valeur exacte est mm. Les erreurs relatives valent % pour , % pour , % pour , puis % pour : . Pour et , l'erreur tombe à l'arrondi: ce sont les seuls maillages de la liste dont un nœud tombe sur la marche, puisque . Entre ces deux valeurs, l'erreur dépend de la qui la chevauche, plus que de la taille des éléments.
Trois observations précisent le mécanisme.
- Les efforts normaux sont justes hors de l'élément qui chevauche. La barre est isostatique: les efforts se déduisent de la statique seule. Avec trois éléments, le même calcul donne , et kN: les éléments 1 et 3, entièrement d'un côté de la marche, portent l'effort exact; l'élément 2 porte une valeur intermédiaire, parce que la force de 20 kN a été répartie sur ses deux nœuds (4 kN et 16 kN).
- L'élément qui chevauche est trop raide. Sa rigidité est la moyenne arithmétique de sur l'élément, alors qu'un tronçon fait de deux sections en série a la souplesse somme , c'est-à-dire la moyenne de la rigidité. La moyenne arithmétique est toujours plus grande que l'harmonique: la barre calculée est trop raide, et trop petit. C'est exactement la remarque de l'exercice 1.2 sur la barre conique, et c'est une conséquence du chapitre 3: la solution de Galerkin minimise l'énergie potentielle dans , et une approximation conforme est toujours au sens de l'énergie.
Une interface entre deux matériaux
Le saut de section de la barre étagée est un saut de , et le théorème 4.3 dit qu'un tel saut ne gêne rien s'il tombe sur un nœud. Le cas d'une barre faite de deux matériaux est le même, mais il permet de séparer la discontinuité du coefficient de celle de la charge: dans la barre étagée, la force de 20 kN et le changement de section tombaient au même point.
Une force concentrée à l'intérieur
Le troisième cas est le plus instructif, parce que le théorème 4.3 y prédit un résultat surprenant: une force concentrée entre deux nœuds n'altère pas l'exactitude nodale, à condition que sa charge nodale cohérente soit calculée exactement. L'erreur est ailleurs.
Dans l'exemple 4.4, on maille la barre en éléments égaux de m. Quel effort normal , en kN, le calcul affiche-t-il dans l'élément qui porte la force?
La règle, et l'explorateur
Les trois cas se résument en une règle de maillage, qui vaut bien au-delà de la dimension 1.
L'explorateur suivant réunit tout ce chapitre. Il assemble et résout le problème sur bloqué aux deux bouts, avec trois charges — uniforme, sinusoïdale, ponctuelle en — et deux matériaux possibles: homogène, ou sur et au-delà, image sans dimension de la barre acier–aluminium de l'exemple 4.3. À droite, la matrice complète est dessinée comme une bande colorée; les lignes et colonnes grisées sont celles des deux nœuds bloqués, que l'élimination met de côté.
Faites varier le nombre d'éléments , puis la charge et le matériau. La matrice reste tridiagonale quel que soit : trois termes non nuls par ligne. Avec un seul matériau, l'erreur nodale reste au niveau de l'arrondi pour les trois charges, même quand la force ponctuelle tombe entre deux nœuds. Avec deux matériaux, elle n'y reste que si un nœud tombe sur l'interface, c'est-à-dire pour , ou .
Quatre expériences. D'abord, faites varier avec la charge uniforme: la matrice garde ses trois diagonales, et seule sa taille change; l'erreur nodale reste au niveau de l'arrondi, et l'erreur maximale, entre les nœuds, vaut ( pour ). Ensuite, choisissez la charge ponctuelle: pour , , , ou , un nœud tombe sous la force, et l'erreur maximale s'effondre à l'arrondi; pour les autres valeurs, les nœuds restent exacts mais la ligne brisée coupe le sommet. passez aux deux matériaux: la diagonale de la matrice change de teinte à l'interface, plus foncée dans la partie raide, et l'erreur nodale n'est plus nulle, sauf pour , et . , combinez les deux matériaux et la charge ponctuelle: il faut que soit multiple de 3 de 5 pour que tout soit exact, ce qui, jusqu'à 16, n'arrive que pour .
Synthèse
- Un maillage 1D est décrit par une table des coordonnées et une table de connectivité; celle-ci est le seul lien entre la numérotation locale des fonctions de forme et la numérotation globale de , et la numérotation des nœuds détermine la largeur de bande de la matrice.
- L'élément linéaire a deux fonctions de forme, et sur l' , ramené à l'élément physique par de jacobien ; ce sont les restrictions des fonctions chapeau.
Un élément linéaire va de m à m. Quelle est l'abscisse du point de coordonnée de référence , en mètres?
Exercices
Vous pouvez afficher le corrigé directement sous chaque énoncé après avoir cherché la solution.
Une barre de longueur 3 m, de rigidité constante, est maillée en trois éléments de 1 m, de raideur chacun. Le mailleur a numéroté les nœuds ainsi: nœud 1 en , nœud 3 en m, nœud 4 en m, nœud 2 en m.
Une barre d'acier de longueur m et de rigidité kN est encastrée en et libre en . Elle porte une charge axiale répartie qui décroît linéairement de kN/m à l'encastrement à zéro à l'extrémité: .
On impose à la barre étagée de l'exemple 4.2, sans aucune force, un allongement mm à son extrémité, l'encastrement restant en place: c'est un essai à déplacement imposé.
- Partager les nœuds en libres et prescrits et écrire le système réduit (4.15).
- Calculer et l'effort normal dans la barre.
On considère une barre de rigidité constante, bloquée à ses deux extrémités, portant une force en un point intérieur à l'élément d'un maillage quelconque, et aucune autre charge. La charge nodale cohérente est utilisée.
On considère la barre encastrée en et libre en , , , , de rigidité constante sur chaque élément d'un maillage , et dont le vecteur des charges est intégré exactement.
Références
- Fish, J. et Belytschko, T., A First Course in Finite Elements, Wiley, chap. 2 à 5 (barres, fonctions de forme, assemblage, conditions aux limites et réactions en une dimension).
- Hughes, T. J. R., The Finite Element Method: Linear Static and Dynamic Finite Element Analysis, Dover, chap. 1 (le problème en une dimension, l'élément de référence, l'exactitude nodale et la superconvergence des dérivées).
- Cook, R. D., Malkus, D. S., Plesha, M. E. et Witt, R. J., Concepts and Applications of Finite Element Analysis, Wiley (assemblage, conditions de Dirichlet, pénalisation, numérotation et largeur de bande).
- Dhatt, G., Touzot, G. et Lefrançois, E., Méthode des éléments finis, Hermès / Lavoisier (en français; programmation de l'assemblage par tables de connectivité).
- Strang, G. et Fix, G., An Analysis of the Finite Element Method, Prentice Hall / Wellesley-Cambridge (l'exactitude nodale en dimension 1 et la fonction de Green).
- Zienkiewicz, O. C., Taylor, R. L. et Zhu, J. Z., The Finite Element Method: Its Basis and Fundamentals, Butterworth-Heinemann (fonctions de forme, élément de référence et charges nodales cohérentes).