8. Problème instationnaire (hyperbolique)#

8.1. Etude de la vibration d’une corde souple#

Considérons une corde tendue suivant l’axe Ox, initialement au repos avec une longueur L et soumise à une tension T0[1], que l’on déforme à l’instant t=0. C’est le problème de la corde pincée d’un instrument de musique (clavecin), dont on se propose d’étudier la vibration (figure Fig. 8.1).

../_images/corde.png

Fig. 8.1 corde tendue souple au repos et à un instant t#

On suppose que la corde est sans raideur (ou souple), c’est à dire que la résultante des contraintes est la tension T(x) qui reste tangente à la corde. On néglige en particulier le moment de flexion qui apparaît si la corde possède une raideur (cas d’une poutre). On tiens cependant compte de l’élasticité de la corde, qui peut s’allonger proportionnellement à une variation de tension.

Considérons un élément de corde MM de longueur ds (figure Fig. 8.1) situé en un point M de la corde. Au repos ce brin de corde MM(0) est situé au point M(x,0) d’abscisse x, d’ordonnée y=0 et a pour longueur dx : M=M(x+dx,0)

A l’instant t il a subit un déplacement longitudinal u(x,t) et un déplacement transversal y(x,t). Le point M se trouve donc en x+u et y: M=M(x+u,y) et le point M en x+u+du et y+dy: M=M(x+dx+u+du,y+dy). Le brin MM fait un angle θ(x,t) avec l’axe Ox et sa longueur ds vaut donc ds2=(dx+du)2+dy2, soit puisque les variables ne dépendent que de x et du temps t:

ds=dx(1+ux)2+(yx)2

L’angle θ vérifie les relations suivantes:

sin(θ)=dyds=yx(1+ux)2+(yx)2 et cos(θ)=dx+duds=1+ux(1+ux)2+(yx)2

Pour de petits déplacements (dudx,dydx), et à un instant t fixé on en déduit les relations suivantes:

sx1+ux+12(yx)2
sin(θ)yx et cos(θ)1

L’allongement relatif du brin MM s’écrit:

dsdxdx=dsdx1ux+12(yx)2

La tension dans la corde est donc la somme de la tension initiale T0 et d’une contrainte élastique, qui d’après la théorie de l’élasticité linéaire s’écrit:

T(x,t)=T0+ES(dsdx1)T0+ES(ux+12(yx)2)

E est le module d’Young de la corde et S sa section.

L’équation d’équilibre dynamique pour le brin MM de masse volumique ρ:

(8.1)#ρSdx2yt2=x(Tsin(θ))dxx(T0yx+ES(ux+12(yx)2)yx)dxρSdx2ut2=x(Tcos(θ))dxx(T0+ES(ux+12(yx)2))dx

En ne conservant que les termes du premier ordre (petites oscillations), on obtiens un système de 2 équations découplées :

(8.2)#ρS2yt2T02yx2ρS2ut2x(ESux)

qui sont 2 équations des ondes, correspondant respectivement à des ondes de flexion ( y(x,t)) de célérité cy=T0ρS et des ondes de compression (u(x,t)) de célérité cu=ESρS.

Dans la cas d’une corde de section S constante, ces 2 équations s’écrivent sous la forme d’une équation des ondes avec une célérité constante c0 (c0=cy ou c0=cu):

(8.3)#2vt2=c02v2x2

Dans le cas d’une corde de section S(x) non constante (de section moyenne S0), cette équation s’écrit:

(8.4)#α2vt2=c02x(βvx)

avec α=SS0, β=1 et c0=T0ρS0 pour l’onde de flexion y(x,t) et α=SS0, β=α et c0=ES0ρS0 pour l’onde de compression u(x,t).

A cette équation on ajoutte des conditions aux limites de déplacement nul aux 2 extrémitées:

(8.5)#v(0,t)=0,v(L,t)=0

et des conditions initiales (CI) : déformation v0(x) sans vitesse initiale

(8.6)#v(x,0)=v0(x),vt(x,0)=0

8.2. Formulation faible#

La formulation faible de l’équation des ondes (8.3) s’obtient classiquement en multipliant par une fonction test w(x)=δv (variation de v(x,t) à un instant t fixé), et en intégrant l’équation à t fixé sur le domaine de calcul:

0L2vt2wdx=0Lc02v2x2wdx

En intégrant par partie le membre de droite, et en utilisant les conditions aux limites de Dirichlet (CL) sur v (8.5), qui imposent des conditions homogènes sur la variation w

w(0)=0 et w(L)=0

on obtiens la formulation faible suivante:

(8.7)#{Trouvez v(x,t) vérifiant les CI et CL et à chaque instant t0L2vt2wdx+c020Lvxwxdx=0w(x) t.q. w(0)=w(L)=0

Cette formulation faible correspond au théorème des travaux virtuels, qui traduit la condition de minimisation de l’action A=0τL(v,v˙)dt, avec un Lagrangien:

L(v,v˙)=120L(vt)2dxT=énergie cinétique12c020L(vx)2dxU=énergie potentielle

La formulation variationnelle associée s’écrit alors:

{Trouvez v(x,t) vérifiant les CI et CLminimisant l'action A=0τL(v,v˙)dt

8.3. Formulation éléments finis#

Pour construire une approximation vh par éléments finis de la solution de (8.6), on utilise un maillage Mh du domaine de calcul Ω=[0,L] en ne éléments {ek=[xk1,xk]}k=1,ne associé à ne+1 noeuds {xi}i=0,ne. On choisit ensuite une interpolation en espace de la solution. Pour notre étude nous utiliserons des éléments finis P1. La solution approchée s’écrit alors comme une combinaison linéaire des fonctions de base {Φi(x)}associées à ce maillage Mh et à l’interpolation P1:

vh(x,t)=i=0nevi(t)Φi(x)

En tenant compte des conditions aux limites (8.5), et des propriétés des fonctions de base Φi(xj)=δij:

vh(0,t)=v0(t)=0 et vh(L,t)=vne(t)=0

l’approximation possède ne1 degrés de liberté et s’écrit:

vh(x,t)=i=1ne1vi(t)Φi(x)

On note que, contrairement aux problèmes statiques étudiés précédemment, les coefficients vi de la combinaison linéaire (valeurs nodales) dépendent maintenant du temps. Les fonctions tests wh(x) étant des variations (à t fixé) de la solution vh, c.a.d. une combinaison des fonctions de base, on choisit comme fonctions tests les fonctions de base Φi. La formulation faible discrète s’écrit alors:

j=1ned2dt2(vj)0LΦj(x)Φi(x)dx+c02j=1ne1vj(t)0LΦjxΦixdx=0i=1,ne1

C’est un système de ne1 équations différentielles linéaires du second ordre qui s’écrit sous la forme matricielle:

Md2Vdt2+c02KV=0

M est la matrice de masse, K la matrice de rigidité et V le vecteur inconnu des valeurs nodales:

Mij=0LΦj(x)Φi(x)dx,Kij=0LΦjxΦixdx,Vi=vi(t)

Ce système s’écrit sous la forme suivante (noté le signe - devant la matrice A):

dV2dt2=AV avec A=c02M1K

auquel on ajoute les conditions initiales (8.6). La matrice A est une matrice symétrique et définie positive[2], puisque M et K sont des matrices symétriques définies positives. Apossède donc ne1 valeurs propres réelles positives, que l’on note ωi2.

La solution générale de ce système d’équations différentielles linéaires est donc la somme de ne1 solutions élémentaires où modes propres Λi :

V(t)=i=1ne1αiΛieIωit

les Λi étant les vecteurs propres associés aux valeurs propres ωi2 de la matrice A:

AΛi=ωi2Λi

D’un point de vue mécanique, la vibration du système est une combinaison de modes propres de vibration Λi associés à des fréquences propres 2πωi. On obtiens ainsi les premiers modes de vibrations de la corde.

8.4. Assemblage et calcul des modes propres#

Le calcul des matrices de masse et de raideur utilise la technique d’assemblage décrite dans les chapitres précédents.

On calcule tout d’abord une matrice de masse et de raideur élémentaire sur l’élément de référence, ce qui permet d’obtenir les matrices élémentaires de masse Mk et de raideur Kk sur un élément ek de longueur h:

Kk=c02h[1111],Mk=h[13161613]

On effectue ensuite un assemblage élément par élément pour calculer les matrices globales de raideur K et de masse M, puis on applique les conditions aux limites.

La prise en compte des conditions aux limites est effectuée en éliminant les lignes et les colonnes des 2 noeuds frontières (noeud 1 et noeud N).

Pour calculer les valeurs propres et vecteurs propres de la matrice A=M1K, on résout le problème équivalent de valeurs propres généralisés:

Mx=λKx

en utilisant une fonction de calcul de valeurs propres en spécifiant que les matrices sont symétriques définie positives.

8.5. Résultats#

Le calcul a été fait avec Ne=10 éléments, soit 9 degrés de liberté. Le nombre de modes propres calculés est donc 9.

L’amplitude des modes propres de vibration calculés est tracé sur les figures ci-dessous pour les modes k=1 k=6 et k=8 et comparé à la solution analytique sinkπx.

../_images/mode1.png

Fig. 8.2 mode 1#

../_images/mode6.png

Fig. 8.3 mode 6#

../_images/mode8.png

Fig. 8.4 mode 8#

On constate que la précision diminue lorsque k augmente.

De même, on a tracé la courbe des fréquences propres calculés comparées aux fréquences analytiques k/2

../_images/frequence.png

Fig. 8.5 fréquences propres#

On constate la encore une diminution importante de la précision lorsque k augmente.

Conclusion: on vérifie le théorème de Shanon sur l’échantillonnage !