2022-Inp-Info-Psi-Corrigé

Assemblage et déformation de pièces automobiles

Introduction

L’objectif de cette étude est la modélisation du comportement géométrique de pièces déformables par la méthode des éléments finis. La partie I déterminera le posage des pièces et sera illustrée par un modèle 2 dimensions.

La partie II permettra d’établir le maillage d’une surface plane pour pouvoir appliquer dans la partie III la méthode de résolution par éléments finis. La partie IV présentera la synthèse de l’étude et une application sur le cas d’une portière et de la caisse d’un véhicule.

Dans tout le sujet, il sera supposé que les modules python numpy, matplotlib.pyplot sont déjà importés dans le programme. Les fonctions et méthodes autorisées sont indiquées dans l’annexe à la fin du Document Réponse .

De manière générale, l’utilisation des méthodes autres que les méthodes append et remove ne sera pas acceptée. L’utilisation des fonctions min ou max sera proscrite.

Partie I – Posage des pièces

I.1 – Présentation

Le posage consiste à mettre deux pièces en contact. Les pièces n’étant jamais parfaites, le nombre des points en contact n’est jamais infini comme cela pourrait être le cas lorsque l’on pose 2 surfaces parfaitement planes l’une sur l’autre.

Ces points de contacts ponctuels doivent être identifiés car ils peuvent être à l’origine de variations géométriques et d’efforts particuliers.

Pour illustrer la question du posage, on travaille en 2 dimensions et on se limite au posage d’une pièce sur une autre.

Les géométries ont été discrétisées en un nombre fini de points. Un prétraitement, qui ne sera pas abordé ici, permet d’extraire les coordonnées des points concernés par le posage.

En 2 dimensions, pour chaque point d’abscisse {x_k=k\Delta x}, est associée son altitude {z_k} où {\Delta x} correspond à une longueur élémentaire associée à la discrétisation spatiale.

Ces données sont stockées dans une liste P de N nombres dont l’indice k des éléments permet de définir l’abscisse des points et les valeurs leur altitude. Ainsi, {x_k=k\Delta x} et {z_k=} P[k] où {0\le k\lt N} . La variable Delta_x sera associée à la longueur élémentaire {\Delta x} et sera supposée déjà définie.

Il sera inutile de la passer en argument des différentes fonctions.

Par exemple, sur la figure 2, on associe à la pièce de gauche la liste P1 =[2.,3.,2.,3.].

Pour répondre au problème du posage, il faut commencer par retourner la pièce en effectuant une rotation autour de l’axe {(O,\vec y)}, pour que les points soient bien en vis-à-vis.

On appellera posage optimal, le posage qui vérifie les deux contraintes suivantes : ⋄ le minimum des distances algébriques entre les points en vis-à-vis doit être positif pour valider le fait qu’il n’y ait pas de collision, ⋄ le maximum des distances entre les points en vis-à-vis doit être le plus petit possible.

Pour caractériser le posage, la méthode mise en œuvre consiste à repérer la position des points du maillage de la deuxième pièce par rapport à ceux de la première pièce. On se placera dans le cadre d’une pièce faiblement inclinée (sur la figure 4, l’inclinaison est amplifiée).

Pour chacune des deux pièces, on associe une droite passant par les 2 points de contacts. Les deux droites ainsi définies seront coïncidentes lors du contact entre les deux pièces comme l’illustre la figure 4.

On considère le posage sur les points n et p.

On note {z_{rn}} et {z_{rp}} les altitudes des points de contact sur la pièce référence, et {z_{pn}} et {z_{pp}} les altitudes des points de contact sur la pièce à poser.

L’équation d’une droite qui passe par deux points de contact ayant pour coordonnées {(n\Delta x,z_n)} et {(p\Delta x,z_p)} est : {z=ax+b\tag{1}}avec pour coefficient directeur {a} et pour ordonnée à l’origine {b} tels que :

{a=\frac{z_n-z_p}{(n-p)\Delta x},\quad b=z_n-\frac{n(z_n-z_p)}{n-p}.}

On applique cette fonction au bâti (pièce 1) et à la pièce 2 retournée. Pour un point donné d’une pièce donnée, on définit son altitude comme la distance algébrique entre ce point et le point sur la droite de même abscisse.

Cette altitude sera comptée positivement si le point est au-dessus de la droite et négativement s’il est situé en-dessous.

Pour une abscisse donnée, la distance algébrique entre deux points en vis-à-vis (chacun appartenant à une des deux pièces) est la différence entre leur altitude respective par rapport à la droite.

Partie II – Maillage de la surface

On souhaite appliquer la méthode des éléments finis pour déterminer la déformation de la surface suite aux contraintes qui découlent du posage. Cette méthode nécessite, au préalable, un découpage de la surface en surfaces élémentaires. Cette étape s’appelle le maillage.

Le maillage proposé lors de cette étude est un maillage triangulaire. Chaque surface élémentaire est alors un triangle dont les sommets sont des nœuds de la surface à étudier.

Le maillage est une étape déterminante permettant d’obtenir des résultats cohérents ou non par la méthode des éléments finis. En effet, un mauvais maillage comme illustré à gauche figure 5 peut être à l’origine d’instabilités numériques ou encore d’erreurs d’arrondis.

Au contraire, une triangulation de Delaunay illustrée à droite figure 5 présente de bons résultats lors d’un traitement numérique.

Dans ce maillage, les triangles ont des arêtes courtes, des angles ni trop petits, ni trop grands.

Une triangulation est dite de Delaunay lorsque le cercle circonscrit à chaque triangle ne contient aucun autre point du maillage.

Une triangulation peut devenir une triangulation de Delaunay par une suite de basculements (flips) successifs : quand 2 triangles ne respectent pas la condition de Delaunay, il faut remplacer l’arête commune {P_2P_4} par l’arête {P_1P_3} comme l’illustre la figure 6.

L’algorithme mis en œuvre pour réaliser une triangulation de Delaunay est un algorithme itératif. La construction s’appuie sur un maillage initial de Delaunay.

Ce maillage constitué d’un faible nombre de nœuds est construit de telle sorte que la surface ainsi délimitée contienne l’ensemble des nœuds qui devront être ajoutés à chaque itération.

Ce maillage initial peut être composé d’un ou de plusieurs triangles vérifiant le critère de Delaunay.

Chaque nœud {P_k} restant du maillage est alors ajouté en suivant l’algorithme : ⋄ repérer à quel triangle {T_{\mathrm{in}}} appartient le point {P_k} (par exemple sur la figure 7, {P_4} appartient au triangle {P_0} {P_1} {P_2} ), ⋄ construire 3 triangles à partir du point {P_k} et des sommets du triangle {T_{\mathrm{in}}}.

Ajouter ces triangles à la liste des triangles T. Supprimer le triangle {T_{\mathrm{in}}} de la liste T. La triangulation ainsi obtenue n’est pas nécessairement de Delaunay. Pour chaque nouveau triangle : → vérifier que les triangles voisins respectent la triangulation de Delaunay.

Si ce n’est pas le cas, basculer l’arête commune ( figure 6). → répéter cette opération pour les triangles modifiés après le basculement.

Structure des données

On dispose de N nœuds {P_k} d’une surface plane repérés par leurs coordonnées cartésiennes {(x_k,y_k)}. Ces nœuds sont stockés dans le tableau noeud qui est un tableau Numpy de flottants de taille {N\times2}. Chaque nœud est repéré par son indice dans le tableau noeud.

L’élément k du tableau noeud contient les coordonnées cartésiennes du nœud d’indice k. Ainsi noeud[k] = array([xk,yk]). On souhaite construire un tableau T qui contient des listes à 3 éléments. Ces listes représentent l’indice des nœuds d’un même triangle.

Dans le cas de la figure 7, la liste T s’écrit : T = [[0,1,2],[2,3,0]]. L’ordre des indices des nœuds est sans importance.

L’algorithme de Delaunay doit permettre de construire la liste T des triangles du maillage. Cette liste est initialisée par 2 triangles contenant tous les nœuds du maillage. Ensuite, à chaque nouveau nœud ajouté, on applique l’algorithme de triangulation, pour mettre à jour la liste T.

L’implémentation de cet algorithme nécessite l’utilisation de plusieurs fonctions : ⋆ addT(indN, indT) : fonction qui prend en arguments l’indice d’un nœud indN à l’intérieur du triangle d’indice indT.

Cette fonction modifie la liste T en insérant les nouveaux triangles et en supprimant le triangle d’indice indT ; ⋆ circle(T1) : fonction qui prend en argument une liste à 3 éléments (les indices des sommets d’un triangle) et qui renvoie les coordonnées du centre du cercle circonscrit ainsi que le rayon [[xc,yc],R] ; ⋆ insideT(indN, T1) : fonction qui prend en arguments l’indice d’un nœud indN et une liste correspondant aux indices des sommets d’un triangle T1.

Cette fonction renvoie un booléen suivant si le nœud appartient ou non au triangle considéré ; ⋆ insideC([xc,yc], R, [x,y]) : fonction qui prend en argument les coordonnées du centre, le rayon et les coordonnées d’un point.

Cette fonction renvoie un booléen suivant si le point considéré appartient ou non au cercle de centre (xc,yc) de rayon R ; ⋆ edge(indN, indT) : fonction qui prend en arguments l’indice d’un nœud indN et l’indice d’un triangle indT.

Cette fonction renvoie l’indice du triangle ayant pour arête commune l’arête à l’opposé du nœud considéré ; ⋆ flip(indT1, indT2) : fonction qui prend en arguments 2 indices (associés à des triangles de la liste T).

Cette fonction permet de basculer l’arête commune en modifiant en conséquence la liste T. Les tableaux noeud et T sont des variables globales qu’il sera inutile de passer en arguments des fonctions et qui pourront être utilisées et modifiées à l’intérieur de celles-ci.

Dans la suite du sujet, on considérera que les différents nœuds sont toujours strictement inclus dans un triangle de la liste T et qu’ils ne sont jamais situés sur l’une des arêtes.

Pour vérifier si un point M du plan cartésien appartient à un triangle ABC, on propose de faire le produit vectoriel {\overrightarrow{AB}\wedge\overrightarrow{AM}} et de vérifier que le résultat est du même signe que le produit vectoriel {\overrightarrow{AB}\wedge\overrightarrow{AC}}.

On procédera de la même manière pour chaque sommet en comparant les signes de {\overrightarrow{BC}\wedge\overrightarrow{BM}} à {\overrightarrow{BC}\wedge\overrightarrow{BA}} et enfin de {\overrightarrow{CA}\wedge\overrightarrow{CM}} à {\overrightarrow{CA}\wedge\overrightarrow{CB}}. L’annexe du sujet rappelle la syntaxe et la définition du produit vectoriel sous numpy.

Une fois le nœud d’indice indN inséré, pour se ramener à une configuration de Delaunay, on applique la fonction récursive checkedge(indN, indT) où indN est l’indice du nœud inséré, et indT un triangle qui a pour sommet le nouveau nœud.

Cette fonction doit d’abord vérifier si le nœud inséré est inclus dans les cercles circonscrits des triangles voisins, c’est-à-dire, les triangles qui possèdent une arête commune aux nouveaux triangles.

Par exemple, sur la figure 7, {P_0P_2} est une arête commune à un nouveau triangle formé par l’insertion du nœud {P_4} et au triangle défini par les nœuds {P_0}, {P_2} et {P_3}.

Si le nœud inséré appartient effectivement à l’un des cercles circonscrits, il faut alors procéder à un basculement grâce à la fonction flip.

Une fois ce basculement réalisé, il faut à nouveau s’assurer que le nœud d’indice indN et les triangles ainsi modifiés vérifient le critère de Delaunay en appliquant la fonction checkedge.

Partie III – Détermination de la matrice de rigidité d’une pièce

À partir de ce maillage, l’objectif est maintenant de calculer numériquement la déformation de la pièce étudiée. Dans la suite du sujet, on considérera que toutes les constantes N, E, S, L, Fx, p0 et n ont été préalablement définies de manière globale.

Pour illustrer cette partie, nous travaillons sur un modèle 1D où 2 nœuds du maillage sont reliés par une barre.

III.1 – Modélisation mécanique

Chaque barre se déforme proportionnellement à l’effort subit, comme un ressort.

La raideur dépend de la section de la barre S de sa longueur L et de l’élasticité du matériau E (module d’Young) : {\vec F=k\Delta L\vec v=\frac{ES}{L}\Delta L\vec v.\tag{2}} où {\vec v} est un vecteur unitaire de même direction que la barre ( figure 8).

On considère que les efforts extérieurs ne s’exercent que sur les nœuds.

Le problème consiste à déterminer les déplacements des nœuds. Sur des structures complexes, le système d’équations peut atteindre plusieurs centaines de milliers voire des millions d’inconnues. Une mise en œuvre numérique est alors indispensable.

III.2 – Calcul de la déformation d’une poutre en traction

Pour introduire la méthode de résolution par éléments finis, nous allons travailler sur un problème simple à une dimension : la poutre en traction. La poutre représentée figure 9 est encastrée en O.

Cette poutre de longueur L de section S est soumise à un effort {F_x=100\,\mathrm N} en A et à un effort réparti, également suivant {\vec x}, représenté par une densité linéique de force {p(x)=p_0=500\,\mathrm{N/m}} (non réprésenté).

La théorie de la résistance des matériaux nous fournit les équations suivantes : {\left\{\begin{array}{l}\dfrac{dN}{dx}+p(x)=0\\[8pt]N=ES\dfrac{du}{dx}\end{array}\right.}avec comme conditions aux limites

{\left\{\begin{array}{l}u(0)=0\\[8pt]N(L)=F_x\end{array}\right.} où N correspond à l’effort normal dans la poutre et u le déplacement longitudinal dû à la traction. En combinant ces deux premières équations, on obtient : {ES\frac{d^2u}{dx^2}=-p(x).\tag{3}} Après intégration, on obtient simplement l’expression analytique suivante : {u_{\mathrm{th}}(x)=-\frac{p_0}{2ES}x^2+\frac{F_x+Lp_0}{ES}x.}

III.3 – Tracé de la solution analytique

L’équation différentielle peut être résolue analytiquement dans ce cas simple, mais dans le cadre général des problèmes d’élasticité 2D ou 3D, ce n’est pas possible. On envisage alors une résolution numérique en appliquant la méthode d’approximation que sont les éléments finis.

III.4 – Recherche de la solution par éléments finis

La méthode de résolution par éléments finis consiste à transformer l’équation différentielle en un système linéaire d’équations.

En repartant de l’équation ( 3), la multiplication de l’équation différentielle par une fonction test Φ, puis son intégration permettent d’aboutir à la relation suivante : {\int_0^L ES\frac{d^2u}{dx^2}\Phi(x)\,dx=-\int_0^L p(x)\Phi(x)\,dx} où la fonction test Φ (dont un exemple est donné figure 10) doit vérifier les mêmes conditions aux limites que la fonction u, soit {\Phi(0)=0}.

Une intégration par parties nous amène au résultat suivant : {\int_0^L ES\frac{du}{dx}\frac{d\Phi}{dx}\,dx=\int_0^L p(x)\Phi(x)\,dx+N(L)\Phi(L).\tag{4}} Nous allons alors chercher une solution linéaire par morceaux à l’équation ( 4) ce qui va nous amener à la discrétisation du problème et à son écriture matricielle.

Nous cherchons alors u(x) sous la forme : {u(x)=\sum_{k=0}^n u_k\Phi_k(x)\tag{5}} où les fonctions {\Phi_k} sont appelées fonctions de forme illustrées sur la figure 10.

Ces fonctions ont les propriétés suivantes : ⋄ {x_k=k\,dx=k\dfrac{L}{n}} et {0\le k\le n} ; ⋄ {\Phi_k} est continue ; ⋄ {\Phi_k} est affine sur {[x_{k-1},x_k]} et sur {[x_k,x_{k+1}]}; ⋄ {\Phi_k(x_k)=1} ; ⋄ {\forall x\in[0,x_{k-1}]\cup[x_{k+1},L],\ \Phi_k(x)=0} ; ⋄ même si les fonctions {\Phi_k} ne sont pas dérivables en tout point de {[0,L]}, on considèrera les dérivées {\dfrac{d\Phi_k}{dx}} comme des fonctions continues par morceaux, en les prolongeant de manière quelconque aux points de ”cassure”.

Le choix du prolongement n’aura de toute façon aucun impact sur les calculs d’intégrales demandés.

Les scalaires {u_k} constituent les inconnues du problème et représentent le déplacement associé à chaque nœud k.

Une fois les déplacements {u_k} déterminés, il est possible de calculer la fonction déplacement u en n’importe quel point x de la poutre grâce à l’équation ( 5).

Nous pouvons alors réécrire le problème sous forme matricielle. En notant U et Phi les vecteurs respectivement associés aux déplacements {u_k} et aux fonctions {\varphi_k(x)}, on peut écrire la fonction u de la manière suivante : {U=\left[\begin{array}{c}u_0\\[8pt]u_1\\[8pt]\vdots\\[8pt]u_n\end{array}\right],\quad \Phi(x)=\left[\begin{array}{c}\Phi_0(x)\\[8pt]\Phi_1(x)\\[8pt]\vdots\\[8pt]\Phi_n(x)\end{array}\right]} {u(x)={}^tU\Phi(x).\tag{6}} où {{}^tU} correspond à la transposée du vecteur U .

La dérivée de u par rapport à x s’écrit alors simplement : {\frac{du}{dx}={}^tU\frac{d\Phi}{dx}.\tag{7}} Pour résoudre l’équation ( 4) il faut alors déterminer la matrice {M=\displaystyle\int_0^L\frac{d\Phi}{dx}{}^t\!\frac{d\Phi}{dx}\,dx}.

Les fonctions de forme ainsi choisies permettent d’aboutir à une matrice M tridiagonale : {M=C\left[\begin{array}{ccccc}\alpha/2&\beta&&&\\[8pt]\beta&\alpha&\beta&&\\[8pt]&\ddots&\ddots&\ddots&\\[8pt]&&\beta&\alpha&\beta\\[8pt]&&&\beta&\alpha/2\end{array}\right].}

Le problème peut alors se réécrire, à partir de l’équation ( 4) de la manière suivante : {ESMU=P+F_x\left[\begin{array}{c}0\\[8pt]\vdots\\[8pt]0\\[8pt]1\end{array}\right].\tag{8}} où P est le vecteur associé à l’effort appliqué aux nœuds de la poutre. Par la suite, on considère que le vecteur P et la matrice M ont été définis comme des tableaux numpy.

Cette équation peut se mettre sous la forme suivante et correspond à un système linéaire de n équations à n inconnues : {AU=B.\tag{9}}

La résolution d’un tel système linéaire peut se faire à l’aide de l’algorithme du pivot partiel de Gauss. La fonction resolution prend en arguments un tableau R de dimensions {(n\times n)} de type flottant et un vecteur à n éléments de type flottant.

Cette fonction renvoie un vecteur à n éléments solution de l’équation {RX=Y}. Les fonctions combinaison et echanger_ligne permettent respectivement de faire une combinaison linéaire entre deux lignes du système et d’échanger deux lignes du système.

Une fois le problème mis sous forme triangulaire, il reste à créer le vecteur des solutions X. C’est la phase de remontée.

Partie IV – Analyse et comparaison des résultats

IV.1 – Illustration de la question du posage

On modélise des défauts de position des points de la caisse afin de visualiser leur impact sur la position de la portière. Ainsi, on peut déterminer pour chaque type de défaut dans quel ordre il est préférable d’assembler la portière sur la caisse.

On va proposer deux scénarios d’assemblage : − scénario 1 : commencer par poser la caisse au niveau des points 1, 2 et 5 ( figure 12 en haut) avant de la déformer pour assembler la charnière basse ; − scénario 2 : commencer par poser la caisse au niveau des points 3, 4 et 5 ( figure 12 en bas) avant de la déformer pour assembler la charnière haute.

IV.2 – Illustration de l’influence de la finesse du maillage

Pour visualiser l’impact de la taille des triangles du maillage, on réalise plusieurs simulations numériques.

Un déplacement est imposé aux 5 points de liaison avec la caisse et un chargement est imposé sur tous les nœuds du pourtour de la portière comme l’illustre la figure 11(b).

Tableau 1 – Tableau de comparaison des différentes simulations

Annexe – Rappels des syntaxes en Python

Remarque : sous Python, l’import du module numpy permet de réaliser des opérations pratiques sur les tableaux : from numpy import *. Les indices de ces tableaux commencent à 0.