Banque commune École Polytechnique - ENS de Cachan PSI
Session 2008
Épreuve de MODÉLISATION
Durée : 5 heures
Aucun document n'est autorisé
L'usage de calculatrices de poche à alimentation autonome, non imprimante et sans document d'accompagnement, est autorisé, une seule à la fois étant admise sur la table ou le poste de travail (circulaire no99 018 du 1 er février 1999).
Si, au cours de l'épreuve, un candidat repère ce qui lui semble être un erreur d'énoncé, il le signale sur sa copie et poursuit sa composition en expliquant les raisons des initiatives qu'il est amené à prendre.
Modélisation de la propagation des ondes de choc dans la structure d'Ariane 5
Figure 1 : photographie de la fusée Ariane 5
Le support d'étude proposé est la fusée Ariane 5 (voir figure 1). Créée en 1996, elle est la dernière génération d'une série de lanceurs de satellites initiée en 1979. Elle a pour fonction de mettre en orbite deux satellites. Depuis son premier vol, 31 lancements ont été effectués, avec un taux de réussite supérieur à 90%.
Cette fusée est constituée des éléments suivants (voir figure 2) :
deux propulseurs latéraux : ils permettent de fournir 90% de la poussée au décollage et sont largués lorsqu'ils sont vides (au bout de deux minutes).
un étage principal : il contient l'oxygène et l'hydrogène liquides qui permettent de faire fonctionner le moteur principal pendant 10 minutes après le décollage. Une fois vide, il est largué.
un étage secondaire : il contient les satellites et un moyen de propulsion qui se met en route lorsque l'étage principal est largué. Il permet d'amener les satellites à leur destination finale.
Figure 2 : constitution de la fusée Ariane 5.
Le diagramme FAST partiel de la figure 3 résume les fonctions que doit assurer chaque élément de la fusée.
Figure 3 : diagramme FAST partiel d'une fonction de la fusée Ariane 5.
Pour séparer l'étage principal de l'étage secondaire, on utilise trois cordons explosifs (pyrotechniques) qui découpent la fusée (voir figure 4).
Q1 : Expliquer pourquoi il est plus judicieux d'utiliser trois cordons explosifs plutôt qu'un seul.
Figure 4 : séparation de l'étage principal et de l'étage secondaire, à l'aide de trois cordons explosifs.
Cette découpe explosive engendre des ondes de chocs qui peuvent se propager dans l'étage secondaire jusqu'aux satellites et éventuellement les endommager. Pour éviter de tels dommages, une solution est d'amortir ces ondes de choc, mais cela nécessite la compréhension de leur propagation au sein de l'étage secondaire. Les ingénieurs du CNES (Centre National d'Étude Spatiale, organisme chargé de maîtriser les technologies spatiales) ont donc lancé un travail de modélisation de cette propagation d'ondes de choc, en s'appuyant largement sur de la simulation numérique afin d'éviter trop de tests coûteux sur des structures réelles.
La démarche de modélisation adoptée est présentée sur la figure 5 .
Figure 5 : démarche de modélisation en Sciences de l'Ingénieur.
Cette démarche se compose symboliquement de deux domaines : le domaine physique où l'on fait des mesures sur un système réel, et le domaine virtuel où l'on réalise des simulations sur un système virtuel (analytique, numérique, ...) qui doit représenter le comportement du système réel. Le domaine physique est nécessaire car il est indispensable de s'appuyer sur l'observation des phénomènes réels pour proposer un modèle de comportement et une modélisation du milieu extérieur. Le solveur est un outil qui permet de résoudre les équations et qui donne le résultat de la simulation. La modélisation est pertinente si l'écart entre les réponses expérimentales et les résultats issus de la simulation est petit.
Le sujet est divisé en trois parties. La première partie amène le candidat à se positionner dans le domaine physique en étudiant les mesures effectuées lors d'un essai de découpe par cordon pyrotechnique. Les deuxième et troisième parties amènent le candidat à se positionner dans le domaine virtuel pour modéliser la propagation des ondes de choc hors plan de flexion (deuxième partie) et membranaires (troisième partie). S'il le souhaite, le candidat pourra tirer profit de l'indépendance des trois parties.
Partie 1 : essais expérimentaux et observation des phénomènes mis en jeux
L'objectif de cette partie est d'étudier les mesures effectuées lors d'un essai de découpe pyrotechnique afin de proposer certaines hypothèses de modélisation. Toutes les étapes présentes dans le domaine physique de mesure de la démarche de modélisation (voir figure 5 ), rappelées dans la figure 6 , seront traitées.
Figure 6 : étapes du domaine physique de mesure de la démarche de modélisation.
Système étudié et chargement imposé
Faire un essai sur une structure réelle étant très coûteux, les ingénieurs ont décidé de réaliser un essai de découpe par cordon pyrotechnique sur un système simplifié, constitué d'une simple plaque rectangulaire, découpée par un seul cordon pyrotechnique (voir figure 7). Ce système peut paraître simpliste face à la complexité d'Ariane 5 , mais il permet de bien comprendre les phénomènes mis en jeu lors de la découpe pyrotechnique et de poser les bonnes hypothèses pour modéliser la propagation des ondes de choc dans le système réel muni de ses trois cordons pyrotechniques. La découpe par le cordon pyrotechnique se fait à une vitesse de 7 km/s.
Figure 7 : moyens d'essai sur une plaque rectangulaire.
Capteurs utilisés
Pour mesurer les déformations de la plaque, des jauges de déformation sont réparties sur la surface de la plaque ( 80 jauges en tout, 40 de chaque côté de la plaque).
Q2 : Expliquer le principe de fonctionnement d'une jauge de déformation.
Il est classique de séparer les mouvements de la plaque en mouvements membranaires de traction compression cisaillement (mouvements selon les directions x⃗ et z⃗ ) et en mouvement hors plan de flexion (mouvements selon y⃗ ). Pour mesurer séparément ces deux types de mouvements, on utilise les mesures de deux jauges de déformation, l'une collée sur un côté de la plaque, l'autre collée sur l'autre côté de la plaque, et on fait la demi somme ou la demi différence des signaux mesurés par les jauges.
Q3 : Déterminer si le signal mesuré par la demi somme correspond au mouvement membranaire et le signal mesuré par la demi différence correspond au mouvement hors plan, ou inversement.
Réponses expérimentales obtenues
Les réponses expérimentales et quelques simulations numériques permettent de bien comprendre ce qui se passe à la suite de l'explosion dans cet essai de découpe pyrotechnique. La figure 8 regroupe les réponses obtenues pour la propagation des ondes de flexion dans la plaque (mouvements selon y⃗ ).
Figure 8 : propagation des ondes de flexion dans la plaque (amplitude selon y⃗ ).
Q4 : Expliquer en quoi le côté dispersif de la propagation des ondes de flexion apparaît sur la figure 8.
La figure 9 regroupe des résultats similaires mais concernant les ondes membranaires dans la plaque (mouvements selon les directions x⃗ et z⃗ ).
Figure 9 : propagation des ondes membranaires dans la plaque (amplitude selon z⃗ ).
Sur la figure 9 , on voit apparaître deux ondes de choc : l'onde P et l'onde S . La mécanique des milieux continus appliquée aux plaques permet de montrer que ces ondes, non dispersives, ont une vitesse de propagation C_P = √((λ + 2μ)/ρ) et C_S = √(μ/ρ), où λ, μ et ρ sont des paramètres matériaux donnés en annexes.
Q5 : À partir des données, retrouver les angles θ_P = 54^∘ et θ_S = 26^∘ que l'on peut mesurer sur la figure 9. Ne pas hésiter à s'appuyer sur un dessin pour les explications.
Q6 : Expliquer pourquoi il est possible de modéliser la propagation des ondes de flexion par un phénomène de propagation à 1 dimension, alors que la modélisation de la propagation des ondes membranaires nécessite l'établissement d'un modèle à 2 dimensions.
Partie 2 : modélisation de la propagation des ondes de flexion
L'objectif de cette partie est de proposer un modèle et de simuler la propagation des ondes de flexion dans la plaque découpée. Toutes les étapes présentes dans le domaine virtuel de simulation de la démarche de modélisation (voir figure 5), rappelées dans la figure 10, seront traitées.
Figure 10 : étapes du domaine virtuel de la démarche de modélisation.
Modélisation du comportement de la plaque
On considère la plaque rectangulaire de l'essai de la partie 1 (voir figure 11). Les dimensions de cette plaque sont L selon x⃗, e selon y⃗ et b selon z⃗. Le repère d'étude R(O; x⃗, y⃗, z⃗) est tel que tous les points M de la plaque sont repérés par OM^(→−) = xx⃗ + yy⃗ + zz⃗, avec x ∈ [0, L], y ∈ [(− e)/2, e/2] et z ∈ [(− b)/2, b/2].ρ est la masse volumique de la plaque. Les valeurs numériques des différents paramètres sont données en annexe du sujet. On appelle "plan médian" le plan à y = 0 et "section droite" le plan de la plaque à x constant.
Figure 11 : géométrie de la plaque.
Suite aux résultats d'essais de la partie 1, on s'intéresse à la propagation des ondes de flexion dans le plan (O; x⃗, y⃗) : les points du plan médian ne se déplacent que sur x⃗ et sur y⃗. On suppose que le déplacement des points est décrit par le vecteur u⃗(x, y, z) = − yθ(x, t)x⃗ + v(x, t)y⃗.
Q7 : Les ingénieurs ne voulant étudier que les phénomènes vibratoires de flexion dont les fréquences sont inférieures à 10 kHz (pour lesquels les ondes se propagent à une vitesse inférieure à 760 m.s^(− 1) ), calculer la longueur d'onde des vibrations à simuler, et justifier le fait que le déplacement selon x⃗ puisse être considéré comme linéaire dans l'épaisseur de la plaque (épaisseur selon y⃗ ). Proposer alors une interprétation physique de la variable θ(x, t).
Les efforts à l'intérieur de la plaque sont décrits par un torseur : {T(x, t)} = {F(x, t)y⃗; M(x, t)z⃗}, où G_x (voir figure 12) est le centre de gravité de la section droite à l'abscisse x, et F(x, t)y⃗ et M(x, t)z→ sont l'effort et le moment qu'exerce la partie de la plaque d'abscisses [ x; L ] sur la partie de la plaque d'abscisses [0; x]. Le volume élémentaire de plaque dV_x compris entre les abscisses x et x + dx est donc soumis aux deux torseurs − {T(x, t)} et {T(x + dx, t)}.
Figure 12 : géométrie de la plaque déformée.
Q8 : En appliquant le théorème de la résultante dynamique sur dV_x en projection sur y⃗, déterminer une équation entre v(x, t) et F(x, t).
Q9 : En appliquant le théorème du moment dynamique sur dV_x, en projection sur ( G_x; z⃗ ), déterminer une équation entre F(x, t) et M(x, t). On supposera que l'inertie de dV_x est négligeable.
Q10: En combinant les deux équations précédentes, montrer que l'équation du mouvement de dV_x s'écrit : ρ be (∂^2 v)/(∂t^2) = − (∂^2 M)/(∂x^2).
Le torseur {T(x, t)} est décrit par un champ surfacique d'effort df⃗(x, y, t) = (σ(x, y, t)x⃗ + τ(x, y, t)y⃗)dS. Le théorie de l'élasticité linéaire permet de montrer que σ = E/(1 − v^2)(∂u_x)/(∂x), où E et v sont des paramètres caractéristiques du matériau utilisé (voir annexe), et u_x est la composante sur x⃗ du déplacement u⃗.
Q11 : En intégrant les efforts surfaciques sur une section droite, déterminer une relation entre M(x, t) et θ(x, t).
On se place dans l'approximation d'Euler Bernoulli, dans laquelle la section droite est supposée rester perpendiculaire au plan médian (voir figure 13).
Figure 13 : déformation de la plaque.
Q12: En utilisant le fait que la section droite reste perpendiculaire au plan médian, déterminer une relation entre θ(x, t) et v(x, t).
Q13 : En combinant toutes les équations précédentes, montrer que le mouvement de la plaque est géré par l'équation ρ be (∂^2 v)/(∂t^2) = − E/(1 − v^2)(be^3)/(12)(∂^4 v)/(∂x^4).
On considère une onde de flexion qui se propage selon la direction x⃗ dans la plaque: v(x, t) = e^(i(ωt − kx)).
Q14 : Déterminer l'équation de dispersion reliant k à ω et montrer que cette relation est dispersive (de manière cohérente aux mesures faites dans la partie 1).
Modélisation du milieu extérieur
La force créée par l'explosion peut être modélisée par la fonction triangle de la figure 14.
Figure 14 : modélisation adoptée pour la force créée par l'explosion du cordon pyrotechnique.
Utilisation d'un solveur analytique
On se propose d'utiliser un solveur analytique basé sur les séries de Fourier pour résoudre les équations de propagation des ondes. Le temps d'observation des phénomènes sera noté T.
On rappelle que toute fonction g(t), périodique de période T , peut être décomposée en séries de Fourier, sur la base de fonctions e^(iω_n t), ω_n = 2πn/T, de la manière suivante :
Q15 : Déterminer les coefficients de Fourier F_n de la fonction F_d(t), fonction périodique de période T, dont la valeur sur [0; T[ vaut f_d(t).
On recherche le déplacement de la plaque sous la forme v(x, t) = ∑_(n = − ∞)^(n = + ∞)v_n(x)e^(iω_n t).
Q16 : Déterminer l'équation différentielle satisfaite par v_n(x).
Q17: Montrer que v(x, t) peut s'écrire sous la forme
où a est une constante à déterminer.
Q18: On ne considère que les ondes qui se propagent suivant les x positifs dans un milieu semi infini x ∈ [0; ∞[. Simplifier l'écriture de v(x), en expliquant les simplifications retenues.
Q19 : Expliquer l'origine de chaque condition aux limites.
Q20 : Exprimer le déplacement v(x) de la plaque en déterminant C_(1n) et C_(3n) pour n ≠ 0 (le cas n = 0 revient à ajouter une constante au déplacement et est sans intérêt pour cette étude).
Affichage du résultat
Le solveur utilisé a permis de trouver le résultat du déplacement en flexion de la plaque soumise à l'explosion. La superposition de ce résultat avec les mesures expérimentales est représentée sur la figure 15.
Figure 15 : comparaison du résultat issu du solveur et des mesures expérimentales à x = 0, 3 m.
Q21: Commenter l'allure générale des courbes (apparition d'oscillations rapides puis d'oscillations lentes). Commenter également l'écart entre les deux courbes.
Partie 3 : modélisation de la propagation des ondes membranaires
L'objectif de cette partie est de proposer un modèle et de simuler la propagation des ondes membranaires dans la plaque. Toutes les étapes présentes dans le domaine virtuel de simulation de la démarche de modélisation (voir figure 5), rappelées dans la figure 16, seront traitées.
Figure 16 : étapes du domaine virtuel de la démarche de modélisation.
Pour représenter la propagation des ondes membranaires, un modèle unidimensionnel n'est pas suffisant car les deux ondes P et S ne se propagent pas dans la même direction (voir partie 1 , les ondes P se propagent à la vitesse C_P = √((λ + 2μ)/ρ), et les ondes S à la vitesse C_S = √(μ/ρ), λ, μ et ρ étant des paramètres matériaux donnés en annexes). Le modèle qui va être construit et le solveur qui sera utilisé seront donc 2 dimensions (2D).
Pour cette partie 3 , la plaque est décrite en 2 D par les variables x et y : x ∈ [0, b], y ∈ [0, L] (voir figure 17). On note Ω = [0, b] × [0, L]. L'épaisseur de la plaque est notée e.
Les bords de la plaque sont notés respectivement Γ_0, Γ_1, Γ_2 et Γ_3.M(x, y) est le point de coordonnées x et y. u⃗ = u(x, y, t)x⃗ + v(x, y, t)y⃗ est le déplacement de M(x, y).
La masse volumique de la plaque est notée ρ.
Les matrices seront indiquées entre crochets : [A] = [A_(xx), A_(xy); A_(yx), A_(yy)]_((x⃗, y⃗)), avec A_(xy) = A_(yx) si la matrice est symétrique.
Figure 17 : géométrie de la plaque.
Modélisation du comportement membranaire de la plaque en 2D
Pour modéliser la plaque, on adopte le point de vue de la mécanique des milieux continus. Pour trouver le déplacement u⃗ = u(x, y, t)x⃗ + v(x, y, t)y⃗ de la plaque soumise à une densité surfacique d'effort F⃗_d(x, t) = F_(dx)(x, t)x⃗ + F_(dy)(x, t)y⃗ sur Γ_0, les équations à résoudre sont les suivantes :
équations d'équilibre intérieur
{(∂σ_(xx))/(∂x) + (∂σ_(xy))/(∂y) = ρ(∂^2 u)/(∂t^2) sur Ω; (∂σ_(xy))/(∂x) + (∂σ_(yy))/(∂y) = ρ(∂^2 v)/(∂t^2) sur Ω
conditions aux limites
{− σ_(xy)(x, y = 0, t) = F_(dx)(x, t); − σ_(yy)(x, y = 0, t) = F_(dy)(x, t); σ_(xx)(x = b, y, t) = σ_(xy)(x = b, y, t) = 0; σ_(xy)(x, y = L, t) = σ_(yy)(x, y = L, t) = 0; σ_(xx)(x = 0, y, t) = σ_(xy)(x = 0, y, t) = 0
Dans l'équation (3), F_(dx)(x, t) et F_(dy)(x, t) sont les composantes de l'effort surfacique imposé sur Γ_0. Dans l'équation (4), Tr est la trace d'une matrice, et [I] = [1, 0; 0, 1] est la matrice identité. [σ] et [ε] correspondent respectivement aux contraintes et déformations engendrées par les efforts surfaciques dans la plaque.
Déterminer une solution exacte de ce groupe d'équations aux dérivées partielles est, dans la plupart des cas, impossible, vu sa complexité. Il faut donc trouver une solution approchée. Pour cela, il est nécessaire d'écrire ce système d'équations sous une autre forme.
où on utilise la notation simplifiée ∫_Ω()dS = ∫_(x = 0)^(x = b)∫_(y = 0)^(y = L)()dxdy.
Q23: En utilisant des intégrations par parties, montrer que ∀u^∗, ∀v^∗
Les questions précédentes ont permis de montrer que {(1), (2), (3), (4)} ⇒ (5). On admettra la réciproque.
Utilisation d'un solveur numérique
1 - Construction d'une approximation de la solution exacte
Pour trouver une solution approchée ( u~, v~ ) de la solution exacte ( u, v ), il suffit de vérifier (5) dans un espace mathématique d'approximation U~ : ∀u~^∗ ∈ U~, ∀v~^∗ ∈ U~
Dans la suite, on prendra comme espace d'approximation l'espace Ũ de dimension n, généré par la base de fonctions d'espace {φ_i}_(i = 1..n). Ainsi, u~(x, y, t) = ∑_(i = 1)^n φ_i(x, y)U_i(t) et v~(x, y, t) = ∑_(i = 1)^n φ_i(x, y)V_i(t). On peut donc écrire
Pour simplifier les écritures, on pose φ = (φ_1; φ_2; ⋮; φ_n), U_– = (U_1; U_2; ⋮; U_n) et V_– = (V_1; V_2; ⋮; V_n). On a donc
u→ = ((^t φ ⋅ U_–)/(^t φ ⋅ V_–))_((x⃗, y⃗))
où ^t φ correspond à la transposée du vecteur φ.
De la même manière
u⃗^∗ = ((^t φ ⋅ U^∗)/(^t φ ⋅ V^∗))_((x⃗, y⃗))
Q26 : Montrer que ^t φ ⋅ ∪ ^∗ = ^t(∪ ^∗) ⋅ φ
Q27 : En introduisant (7) et (8) dans (6), montrer que que le problème s'écrit sous la forme suivante : trouver U_– et V_– tels que ∀U_–^∗, ∀V_–^∗
Donner l'expression détaillée de la matrice [M].
Q28: On pose finalement X_– = (U_–_–). Montrer qu'on aboutit à la résolution du système suivant : [A]X_–_(, tt) + [B]X_– = F_– où l'on précisera le vecteur F, les matrices [A] et [B].
Q29 : Déterminer la dimension des matrices et vecteurs intervenant dans ce système. Montrer que les matrices [A] et [B] sont symétriques.
On s'intéresse maintenant au cas où la base de fonctions {φ_i}_(i = 1..n) est affine par morceaux, prenant appui sur un maillage. Ce cas particulier est appelé méthode des éléments finis en Sciences de l'Ingénieur. Il nécessite la construction d'un maillage (voir figure 18), dans lequel on trouve des éléments et des noeuds. Il y a autant de fonctions φ_i que de noeuds (représentant des points géométriques de coordonnées (xi,yi)). Les fonctions {φ_i}_(i = 1..n) sont affines sur chaque élément. Elles sont appelées fonctions de forme.
Figure 18 : exemples de maillages de la plaque.
De plus, on impose la propriété φ_i(x_i, y_i) = 1 et φ_i(x_j, y_j) = 0 pour i ≠ j. Ainsi u~→_i = U_i(t)x⃗ + V_i(t)y⃗ = ((U_i(t))/(V_i(t)))_((x⃗, y⃗)) représente la valeur du déplacement approché u~→ au noeud i de coordonnées (x_i, y_i).
Ces propriétés et le fait que les fonctions soient affines montrent que les fonctions φ_i(x, y) ne sont non nulles que sur un nombre restreint d'éléments entourant chaque noeud i .
Le calcul des différents termes de chaque matrice ou du vecteur second membre intervenant dans le système s'effectue d'une manière identique. Pour ne pas alourdir les calculs, on s'intéresse uniquement à la détermination de la matrice [M] = ∫_Ω(ρφ^t φ)dS. La composante M_(ij) de cette matrice vaut M_(ij) = ∫_Ω(ρφ_i φ_j)dS et fait intervenir les fonctions de forme des noeuds iet j.
Q30 : Déterminer la valeur des M_(ij) si les noeuds i et j ne sont pas sur un même élément.
Q31 : En utilisant le maillage défini sur la figure 19, préciser la forme de la matrice [M] en indiquant par des croix les termes non nuls et en mettant des 0 pour les termes nuls.
Figure 19 : maillage à 2 éléments et 6 noeuds.
Le calcul de chaque terme non nul M_(ij) = ∫_Ω(ρφ_i φ_j)dS peut s'écrire de la façon suivante : M_(ij) = ∫_Ω(ρφ_i φ_j)dS = ∑_(Eléments en lien avec i et j)(∫_(Ω_E)(ρφ_i φ_j)dS) où Ω_E est l'élément E (en lien avec les noeuds i et j).
On remarque donc qu'il est nécessaire d'intégrer un produit de fonctions sur chaque élément. Le calcul de ces intégrales est le même pour tous les éléments à un changement de variable prêt. Pour ne pas avoir à calculer toutes les intégrales, on construit un élément de référence pour lequel les intégrales sont entièrement déterminées. Par un changement de variables approprié, on peut ensuite obtenir les valeurs des intégrales sur chaque élément réel, quelle que soit sa position, à partir des valeurs obtenues sur l'élément de référence.
On s'intéresse donc à l'élément défini sur la figure 20 (élément de référence) dont les noeuds sont repérés par les lettres a, b, c et d. Cet élément a pour dimensions l_x et l_y. On appelle φ_a(ξ, η), φ_b(ξ, η), φ_c(ξ, η), φ_d(ξ, η) les fonctions de forme associées à chaque noeud : ce sont des fonctions affines des variables ξ (coordonnée selon x⃗ ) et η (coordonnée selon y⃗ ), donc des polynômes de degré 1 en ξ et η. Par exemple, on a φ_a(ξ, η) = a_3 ξη + a_2 ξ + a_1 η + a_0, (a_(0,)a_1, a_(2,)a_3) étant les coefficients constants du polynôme φ_a(ξ, η).
Figure 20 : élément de référence.
Q32 : Déterminer φ_a(ξ, η) (donc les 4 constantes a_0, a_1, a_2 et a_3 ) pour que φ_a soit égale à 1 au noeud a et égale à 0 aux noeuds b, c et d. Tracer en perspective l'allure de cette fonction sur l'élément de référence en justifiant succinctement le tracé.
La matrice de masse sur l'élément de référence s'écrit de la façon suivante :
avec M_(aa) = M_(bb) = M_(cc) = M_(dd)
Un changement de variable permet de passer de l'élément de référence à un élément du maillage de la plaque. On suppose que les longueurs des côtés de chaque élément réel sont I_x et I_y = I_x. Il est donc simplement nécessaire de déterminer la correspondance entre les noeuds de chaque élément du maillage de la figure 19 et ceux de l'élément de référence (le respect de l'ordre des noeuds permet d'utiliser directement la forme de la matrice de masse élémentaire [M^e] ).
Q34 : Recopier sur votre copie et compléter le tableau de correspondance entre les noeuds.
Élément E 1
Élément E 2
Numérotation des
noeuds des éléments
réels
Noeuds de l'élément
de référence
a − b − c − d
a − b − c − d
Nous disposons à présent de toutes les données nécessaires pour calculer les termes non nuls de la matrice [ M ]. Pour cela, on utilise le tableau de correspondance pour positionner la matrice [M^e] dans la matrice [M]. Par exemple, pour l'expression de M_(33), on peut écrire M_(33) = ∫_Ω φ_3^2 dS = ∫_(E_1)φ_3^2 dS + ∫_(E_2)φ_3^2 dS = M_(dd) + M_(aa). On constate donc que M_(33) est obtenu par assemblage d'un terme provenant de l'élément E 1 (terme d -d) et d'un terme provenant de l'élément E2 (terme a-a). Étant donné que la matrice [ M^e ] est la même pour tous les éléments (ils ont tous la même forme et les mêmes dimensions), M_(33) est simplement obtenu en sommant deux termes de la matrice [ M^e ] choisis en fonction de la correspondance entre les noeuds.
Les autres termes non nuls sont déterminés de la même manière.
Q35 : Donner l'expression de la matrice [ M ] en utilisant cette technique d'assemblage et les valeurs de la matrice [M^e] (relation 9).
2 - Étude de l'influence des paramètres sur la solution à un problème donné
Les différentes matrices et le second membre étant déterminés, il est possible de résoudre le système différentiel matriciel [A]X_–_(, tt) + [B]X_– = F_–. On souhaite étudier l'influence des paramètres du solveur de manière à le configurer correctement pour la simulation de la découpe de la plaque.
On s'intéresse alors à la plaque rectangulaire de côté 1 m × 2 m soumise à un chargement sinusoïdal selon la direction y⃗ uniquement sur une petite zone du côté Γ_0. On étudie la propagation des ondes avant réflexion sur les autres bords dans la partie inférieure de cette plaque. La fréquence du chargement est f = 30kHz et son amplitude est unitaire (figure 20).
Figure 20 : cas test : plaque rectangulaire soumise à un chargement sinusoïdal.
Comme indiqué dans la partie 1 , deux ondes se propagent dans ce milieu 2 D : I'onde P à la vitesse C_P et l'onde S à la vitesse C_S, avec C_P > C_S. On note λ_P = (C_P)/f la longueur d'onde des ondes P et λ_s = (C_s)/f la longueur d'onde des onde S.
On étudie l'influence de la taille des éléments choisis ( l_x ) sur la résolution (qualité de la solution approchée, temps de résolution,...). La figure 21 montre un exemple de maillage lorsque I_x = (λ_s)/8. Les résultats représentés sur les figures 22 et 23 montrent le déplacement selon les directions y⃗ et x⃗ pour un même instant mais en modifiant la taille des éléments.
Le tableau ci-dessous indique la taille des éléments choisis pour chaque simulation et le nombre de noeuds n total dans la zone étudiée.
Nom du maillage
a
b
c
d
e
Taille des
éléments
I_x = (λ_s)/2
I_x = (λ_s)/4
I_x = (λ_s)/8
I_x = (λ_s)/(16)
I_x = (λ_s)/(20)
Nombre de
noeuds n
231
840
3160
12324
19110
Q36 : Donner la propriété que semble vérifier la solution éléments finis lorsque la taille des éléments tend vers zéro.
Figure 21 : exemple de maillage utilisé pour la simulation, ici pour I_x = (λ_s)/8.
Figure 22 : déplacement selon y⃗ pour différentes tailles d'éléments.
Figure 23 : déplacement selon x⃗ pour différentes tailles d'éléments.
Le choix de la taille des éléments doit être réalisé de manière astucieuse. En effet, on peut montrer que le temps de calcul est proportionnel au nombre d'inconnues ( 2n ) élevé au carré. Le temps de calcul du cas a) est de 0.01s par simulation.
Q37 : Calculer le temps nécessaire pour les 5 cas de calculs a), b), c), d) et e), et choisir le maillage qui semble le plus pertinent pour obtenir rapidement un résultat juste.
Q38 : Compte tenu du rapport (λ_s)/(I_x) pour les 5 cas de calculs, de la forme sinusoïdale des ondes se propageant dans la plaque et du caractère affine des fonctions de forme, expliquer pourquoi les calculs a) et b) ne semblent pas précis, alors que les calculs c), d) et e) semblent l'être. Ne pas hésiter à s'aider de schémas pour répondre à cette question.
Q39 : Expliquer pourquoi la taille des éléments est choisie en fonction de la longueur d'onde des ondes S et non pas en fonction de celle des ondes P.
Pour résoudre le système différentiel [A]X_–_(, t) + [B]X_– = E, on utilise un schéma temporel d'intégration. On introduit pour cela un pas de temps Δt, et on calcule la solution à différents pas de temps: 0, Δt, 2Δt, 3Δt… Ce pas de temps doit être choisi assez petit pour pouvoir observer les phénomènes engendrés par la sollicitation. Ainsi il doit être plus petit que les temps caractéristiques de la sollicitation (période, durée de présence de la sollicitation...).
On utilise de plus une règle supplémentaire pour le schéma temporel d'intégration qui est la suivante: Δt < (I_x)/c où c est la célérité des ondes dans la plaque et I_x la dimension minimale d'un élément.
Q40: Déterminer le pas de temps Δt maximal à utiliser pour la taille des éléments retenue précédemment.
Après ces différentes analyses effectuées, les paramètres du solveur utilisé (solveur élément fini) sont bien réglés pour obtenir un résultat juste. On s'intéresse maintenant à la résolution numérique du problème de découpe de la plaque. La bande de fréquence intéressante est [0, 30kHz].
Q41 : Déduire des analyses précédentes la taille des éléments à utiliser pour obtenir une représentation correcte de la solution sur la plaque entière. Donner alors le nombre d'inconnues du problème (on arrondira le nombre d'éléments selon x et y à la valeur supérieure).
Q42 : Montrer enfin en quoi un pas de temps de 1, 5.10^(− 6) s semble raisonnable pour étudier le phénomène de découpe (la découpe se fait à 7 km/s et le signal triangulaire de la force de découpe, décrite sur la figure 24 , est défini pour T_0 = 1, 5.10^(− 5) s ).
Figure 24 : modélisation adoptée pour la force créée par l'explosion du cordon pyrotechnique.
Conclusion : Le solveur ainsi paramétré permet d'obtenir des résultats numériques pertinents qui peuvent être comparés aux essais expérimentaux.
Annexes : notations utilisées et paramètres géométriques et matériaux de la plaque étudiée
i
Nombre purement imaginaire ( i^2 = − 1 ) ω e = 6 mm b = 1 m L = 2 m E = 70GPa v = 0, 28 λ = (vE)/((1 + v)(1 − 2v)) μ = E/(2(1 + v)) ρ = 2800 kg.m^(− 3)
k
Pulsation
Epaisseur de la plaque
Largeur de la plaque
Longueur de la plaque
Paramètre matériau de la plaque, appelé module d'Young Paramètre matériau de la plaque, appelé coefficient de Poisson
Paramètre matériau de la plaque
Paramètre matériau de la plaque
Masse volumique de la plaque
Vecteur d'onde
Pas de description pour le moment
Commentaires• X ENS Modélisation PSI 2008
Connectez-vous pour participer aux discussions
Partagez vos avis, posez des questions et échangez avec la communauté
Sujet proche avec corrigéCCINP Modélisation PSI 2024CCINP · 2024Corrigé
.jpg)
.jpg)
.jpg)
.jpg)
.jpg)
.jpg)
.jpg)
