Ecole Nationale d’Ingénieurs
Abderhamane Baba Touré
ENI-ABT
Cours de Dynamique des Structures
M1_S2
2023-2024
Dynamique des Structures
Introduction
1. Concepts de base de la Vibration :
Le sujet de la vibration traite du mouvement oscillatoire des systèmes
dynamiques. Un système dynamique est une combinaison de matière qui
possède une masse et dont les parties sont capables de mouvement relatif.
Tous les corps possédant une masse et une élasticité sont capables de
vibrations. La masse est inhérente au corps, et l'élasticité est due au
mouvement relatif des parties du corps. Le système considéré peut-être très
simple ou complexe. Il peut se présenter sous la forme d'une structure, d'une
machine ou de ses composants, ou d'un groupe de machines. Le mouvement
oscillatoire du système peut être répréhensible, insignifiant ou nécessaire à
l'exécution d'une tâche.
L'objectif du concepteur est de contrôler la vibration lorsqu'elle est
répréhensible et pour améliorer la vibration quand c'est utile, bien que les
vibrations en général sont indésirables. Les vibrations répréhensibles dans
une machine peuvent provoquer le desserrage de pièces, son
dysfonctionnement ou sa défaillance éventuelle. D'autre part, les fonderies
d'agitateurs et les vibrateurs dans les machines d'essai nécessitent des
vibrations. La performance de nombreux instruments dépend du bon contrôle
des caractéristiques vibratoires du dispositif.
L'objectif principal de notre étude est d'analyser le mouvement oscillatoire
des systèmes dynamiques et des forces associées au mouvement.
Le but ultime de l'étude des vibrations est de déterminer son effet sur le la
performance et la sécurité du système sous considération. L'analyse du
mouvement oscillatoire est une étape importante vers cet objectif
Tout mouvement qui se répète après un intervalle de temps est appelé
vibration ou oscillation. Le balancement d’un pendule et le mouvement d’une
corde pincée sont des exemples typiques de vibration. La théorie des
vibrations traite de l’étude des mouvements oscillatoires des corps et forces
associées à eux.
Un système vibratoire, en général, inclut moyen de stockage d’énergie
potentielle (ressort ou élasticité), un moyen de stockage d’énergie cinétique
(masse ou inertie) et un moyen par lequel l’énergie est graduellement perdue
(amortisseur).
La vibration d’un système implique le transfert de son énergie potentielle en
énergie cinétique de l’énergie cinétique en énergie potentielle,
alternativement. Si le système est amorti, une partie de l’énergie est dissipée
à chaque cycle de vibration et doit être remplacée par une source externe si
un état de vibration constante doit être maintenue.
Comme exemple, considérons d’un pendule simple, comme montré sur la
figure. Considérons la masse relâché après lui avoir donné un déplacement
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 2
Dynamique des Structures
angulaire 𝜃. A la position 1 la vitesse de la masse, et donc son énergie
cinétique, est nulle. Mais elle possède une énergie potentielle d’amplitude
𝑚𝑔𝑙(1 − cos 𝜃) par rapport à la position 2. Dès lors que la force gravitationnelle
𝑚𝑔 induit un couple 𝑚𝑔𝑙 sin 𝜃 par rapport au point 𝑂, la masse commence à
balancer vers la gauche à partir de la position 1. Cela donne a la masse une
certaine accélération angulaire dans le sens horaire, au moment où il atteint
la position 2, toute son energie potentielle sera convertie en énegie cinétique.
Mais la masse ne s’arrêtera à la position 2, mais continuera à balancer vers la
position 3. Cependant, comme il passe par la position médiane 2, un couple,
dans le sens antihoraire dû à la gravitation, commence à agir sur la masse et
crée une décélération de cette dernière.
Figure 1 : Pendule Simple
Ce processus continue à se répéter, et le pendule aura un mouvement
oscillatoire. Cependant, dans la pratique, l’amplitude de l’oscillation dimunie
graduellement et le pendule finit par s’arrêter due à la résistance
(amortissement) exercé par le milieu environnant (air). Cela signifie qu’une
partie de son énergie est dissipée à chaque cycle de vibration due à
l’amortissement de l’air.
Le nombre minimal de coordonnées indépendantes nécessaires pour
déterminer complètement les positions de toutes composantes du système à
tout instant définie le nombre degré de liberté.
Le pendule et tout le système qui le compose représente un système à un
degré de liberté.
Par exemple, le mouvement du pendule simple (figure 1) peut être exprimé en
termes d’angle 𝜃 ou en termes de coordonnées cartésiennes 𝑥 et 𝑦. Si les
coordonnées cartésiennes sont utilisées pour décrire le mouvement, il faut
reconnaitre que ces coordonnées ne sont pas indépendantes. Elles sont reliées
par la relation 𝑥 2 + 𝑦 2 = 𝑙 2 , ou 𝑙 est la longueur constante du pendule. Chaque
type de coordonnées peut décrire le mouvement du pendule. Dans l’exemple
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 3
Dynamique des Structures
du pendule simple le choix de 𝜃 comme coordonnée indépentante est plus
convenable que le choix de de 𝑥 ou 𝑦.
Pour le coulisseau de la figure 2 (a) aussi bien la coordonnée angulaire 𝜃 que
la coordonnée 𝑥 peuvent être utiliser pour décrire le mouvement. Sur la figure
2 (b) la coordonnée linéaire 𝑥 peutêtre utiliser pour spécifier le mouvement.
Pour le système torsionnelle (longue barre avec un disque lourd à l’extrémité)
figure 2 (c), la coordonnée angulaire 𝜃 peut-être utiliser pour décrire
le mouvement.
Figure 2 : Systèmes à un degré de liberté
Des exemples de systèmes à deux degrés de liberté sont montrés
respectivement la figure 3.
Figure 3 : Systèmes à deux degrés de liberté
La figure 3 (a) montre un système de deux masses et deux ressorts qui est
décrit par deux coordonnées linéaires 𝑥 1 et 𝑥2 . La figure 3 (b) dénote un
système de deux rotors dont le mouvement peut être spécifier en termes de
𝜃1 et 𝜃2 . Le mouvement du système montré sur la figure 3 (c) peut être
complètement décrit soit par 𝑋 et 𝜃 soit par 𝑥, 𝑦 et 𝑋. Dans le dernier cas
𝑥 𝑒𝑡 𝑦 sont liés par la relation 𝑥 2 + 𝑦 2 = 𝑙 2 où 𝑙 est constante.
Des exemples de systèmes à trois degrés de liberté sont montrés
respectivement la figure 4.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 4
Dynamique des Structures
Figure 4 : Systèmes à trois degrés de liberté
Pour les systèmes montrés sur la figure 4, les coordonnées 𝑥 𝑖 (𝑖 = 1; 2; 3) et
𝜃𝑖 (𝑖 = 1; 2; 3) peuvent être utilisées respectivement pour décrire le
mouvement.
Les coordonnées nécessaires pour décrire le mouvement d’un système
constituent un groupe de coordonnées généralisées, ces coordonnées sont
généralement notées 𝑞1 , 𝑞2 , … et peuvent représenter les coordonnées
cartésiennes et/ou non cartésiennes.
Un grand nombre de systèmes pratiques peuvent être décrites en utilisant un
nombre finis de coordonnées. Certains systèmes, spécialement ceux
impliquant des membres élastiques, ont un nombre infini de degrés de
libertés et ont besoin d’un nombre infini de coordonnées pour être décrites.
C’est le cas des poutres, cordes, membranes…
Les systèmes ayant un nombre fini de degrés de liberté sont dit discrets et les
systèmes ayant un nombre infini de degrés de libertés sont dits continus.
2. Classification des Vibrations
Vibration libre
Si un système, après une perturbation initiale, est laissé en oscillation a lui-
même, la vibration qui s’ensuit est appelée vibration libre. Aucune force
extérieure n’agit sur le système. L’oscillation d’un pendule simple est un
exemple de vibration libre.
Vibration Forcée
Si un système est soumis à une force extérieure (souvent répétitive) la
vibration résultante est dite forcée. L’oscillation qui survient des machines
comme le moteur à combustion interne est un exemple de vibration forcée.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 5
Dynamique des Structures
Si la fréquence de la force extérieure coïncide avec une des fréquences
naturelles du système, le système passe dans un état de vibration qu’on
appelle résonance. Le système subit alors de larges et dangereuse
oscillations. Des défaillances dans les structures comme les buildings, ponts,
turbines, les ailes d’avions ont été associées à l’occurrence de résonance.
Vibration amortie et non amortie
S’il n’y a pas d’énergie dissipée durant l’oscillation, alors la vibration est
appelée vibration non amortie. Tandis que s’il ya une quelconque perte
d’énergie de quelque manière que ce soit, on parle de vibration amortie. Dans
de nombreux systèmes physiques, l’amortissement est si petite qu’elle peut
etre negligée. Cependant la considération de l’amortissement est
extremement importante au voisinage de la résonance.
3. Procédure d’Analyse de Vibration
Un système vibratoire est un système dynamique pour lequel les variables
telles que les excitations (entrées) et réponses (sorties) dépendent du temps.
La réponse d’un système vibrant dépend généralement aussi bien des
conditions initiales que des excitations extérieures. La plupart des systèmes
vibrants sont très complexes, et il est impossible de considérer tous les
détails pour une analyse mathématique. Seules les caractéristiques
importantes sont considérées dans l’analyse pour prédire le comportement
du système dans des conditions d’entrées spécifiques.
Etape 1: Modélisation Mathématique
Le but de la modélisation mathématique est de représenter les
caractéristiques importantes du système en vue de déterminer les équations
mathématiques qui gouvernent le comportement du système.
Etape 2 : Détermination des équations qui gouvernent le mouvement
Une fois qu’on a la modélisation mathématique, on applique le principe
fondamental de la dynamique pour déduire les équations qui décrive la
vibration du système. Ces équations peuvent être obtenues en dessinant le
diagramme du corps libre de toutes les masses concernées.
Le diagramme du corps libre d’une masse peut être obtenu en isolant la masse
et en indiquant toutes les forces extérieures qui s’exercent sur cette masse,
les forces de réactions et les forces d’inertie. Les équations du mouvement
d’un system vibrant sont généralement sous la forme d’un groupe d’équations
différentielles pour un système discret et sous forme d’équations aux
dérivées partielles pour un système continu. Plusieurs approches sont
habituellement utilisées pour déterminer les équations du mouvement, parmi
lesquelles Le seconde loi de Newton du mouvement, le principe de
D’Alembert, et le principe de conservation de l’énergie.
Etape 3 : Solution des équations du mouvement
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 6
Dynamique des Structures
Les équations qui régissent le mouvement doivent être résolues pour trouver
la réponse d’un système vibratoire. En fonction de la nature du problème, on
peut utiliser les techniques suivantes pour trouver la solution : méthodes
standards de résolution d’équations différentielles, les méthodes de la
transformation de Laplace, les méthodes matricielles et les méthodes
numériques utilisant des ordinateurs.
Etape 4 : Interprétation des résultats
La solution des équations qui régissent le mouvement donnent les
déplacements, vitesses et accélérations des différentes masses du système.
Ces résultats doivent être interpréter avec une vision claire du but de l’analyse
vibratoire et des possibles conséquences en conception.
4. Exemple de modélisation masse-ressort-amortisseur
La figure 5(a) montre un motocycliste sur une moto. Développons une
séquence de trois modèles mathématiques du système pour analyser les
vibrations dans la direction verticale. Considérer l’élasticité des pneus,
l’élasticité et l’amortissement des entretoises (dans la direction verticale) les
masses des roues, l’élasticité, l’amortissement et la masse du motocycliste.
Nous commençons avec le modèle le plus simple que nous affinerons
graduellement. Lorsque nous utilisons les valeurs équivalentes de la masse,
l’élasticité et l’amortissement du système, on obtient un modèle à un degré
de liberté figure 5(b). Dans ce modèle, la valeur équivalente de la raideur
(𝑘𝑒𝑞 ) inclus les raideurs des pneus, entretoise et motocycliste.
L’amortissement équivalente constante (𝑐𝑒𝑞 ) inclus les amortissements de
l’entretoise et motocycliste. La masse équivalente (𝑚 𝑒𝑞 ) inclus les masses des
roues, corps du véhicule et motocycliste.
Ce modèle peut être raffiner en représentant les masses des roues, élasticité
des roues, et l’élasticité et l’amortissement des entretoises séparément figure
(c). Dans ce modèle, la masse du véhicule 𝑚 𝑣 et la masse du motocycliste
𝑚 𝑚 sont montrées en une seule masse 𝑚 𝑣 + 𝑚 𝑚 . Lorsque l’élasticité et
l’amortissement du motocycliste sont considéré on obtient la modèle raffiné
de la figure (d).
A noter que ces modèles ne sont pas uniques et peuvent être raffiner en
fonction des hypothèses de notre étude.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 7
Dynamique des Structures
(a)
(b) (c)
(d) (e)
m: motocycliste; v: vehicule; e: entretoise; r: roue; p: pneu; eq: équivalent
Figure 5 : Un motocycliste sur une moto ; un système physique et un
modèle mathématique
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 8
Dynamique des Structures
Chapitre 1 : Vibrations des Systèmes
Discrets à plusieurs degrés de liberté
1. Mise en équation
Nous avons un système typique à 𝑛 degré de liberté sur la figure suivante :
Figure 6 : Système à 𝑛 degrés de liberté typique
Pour chaque masse 𝑚 𝑖 nous avons les caractéristiques suivantes :
xi ∶ le déplacement de la masse par rapport à sa position d’équilibre statique
k i ∶ la raideur d’élasticité
ci ∶ le coefficient d’amortissement
fi ∶ la force d’excitation extérieure.
Les équations du mouvement de ce système à 𝑛 degrés de liberté peuvent être
déterminées à partir du diagramme du corps libre, en appliquant la relation
fondamentale de la dynamique.
Figure 7 : Diagramme du corps libre des masses 𝑚𝑖
Pour la masse 𝑚 1 , on a :
𝑘1 𝑥1 + 𝑐1 𝑥̇ 1 − 𝑘2 (𝑥 2 − 𝑥 1 ) − 𝑐2 (𝑥̇ 2 − 𝑥̇ 1 ) = 𝑚 1 𝑥̈ 1 (1.1)
Pour la masse 𝑚 2 , on a :
𝑘2 (𝑥 2 − 𝑥 1 ) + 𝑐2 (𝑥̇ 2 − 𝑥̇ 1 ) − 𝑘3 (𝑥 3 − 𝑥 2 ) − 𝑐3 (𝑥̇ 3 − 𝑥̇ 2 ) = 𝑚 2 𝑥̈ 2 (1.2)
Pour la masse 𝑚 3 , on a :
𝑘3 (𝑥 3 − 𝑥 2 ) + 𝑐3 (𝑥̇ 3 − 𝑥̇ 2 ) − 𝑘4 (𝑥 4 − 𝑥 3 ) − 𝑐4 (𝑥̇ 4 − 𝑥̇ 3 ) = 𝑚 3 𝑥̈ 3 (1.3)
Ainsi de suite, pour la masse 𝑚 𝑛 on a :
𝑘𝑛 (𝑥 𝑛 − 𝑥 𝑛−1 ) + 𝑐𝑛 (𝑥̇ 𝑛 − 𝑥̇ 𝑛−1 ) − 𝑘𝑛+1 (𝑥 𝑛+1 − 𝑥 𝑛 ) − 𝑐𝑛+1 (𝑥̇ 𝑛+1 − 𝑥̇ 𝑛 ) = 𝑚 𝑛 𝑥̈ 𝑛 (1.4)
On obtient un système de 𝑛 équations différentielles que l’on peut
réorganiser comme suit :
𝑚 1 𝑥̈ 1 + (𝑐1 + 𝑐2 )𝑥̇ 1 − 𝑐2 𝑥̇ 2 + (𝑘1 + 𝑘2 )𝑥 1 − 𝑘2 𝑥 2 = 𝑓1
𝑚 2 𝑥̈ 2 + (𝑐2 + 𝑐3 )𝑥̇ 2 − 𝑐2 𝑥̇ 1 + (𝑘2 + 𝑘3 )𝑥 2 − 𝑘2 𝑥 1 = 𝑓2
𝑚 3 𝑥̈ 3 + (𝑐3 + 𝑐4 )𝑥̇ 3 − 𝑐3 𝑥̇ 2 + (𝑘3 + 𝑘4 )𝑥 3 − 𝑘3 𝑥 2 = 𝑓3
……………………………………………………………….
𝑚 𝑛 𝑥̈ 𝑛 + (𝑐𝑛 + 𝑐𝑛+1 )𝑥̇ 𝑛 − 𝑐𝑛 𝑥̇ 𝑛−1 + (𝑘𝑛 + 𝑘𝑛+1 )𝑥 𝑛 − 𝑘𝑛 𝑥 𝑛−1 = 𝑓𝑛
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 9
Dynamique des Structures
Définissons les matrices [𝑚], [𝑐 ] et [𝑘] suivantes, respectivement matrice de
masse, matrice d’amortissement et matrice de raideur :
Alors le système de 𝑛 équations différentielles peut s’écrire en utilisant
l’écriture matricielle :
[𝑚]𝑥⃗̈ + [𝑐 ]𝑥⃗̇ + [𝑘]𝑥⃗ = 𝑓⃗ (1.5)
𝑥̈1 𝑥̇1 𝑥1 𝑓1
𝑥̈2 𝑥̇2 𝑥2 𝑓2
On a : 𝑥⃗̈ = ( 𝑥̈3 ) , 𝑥⃗̇ = ( 𝑥̇3 ), 𝑥⃗ = (𝑥⋮3 ) et 𝑓⃗ = (𝑓3 ).
⋮ ⋮ ⋮
𝑥̈𝑛 𝑥̇𝑛 𝑥𝑛 𝑓𝑛
Noter que le système masse-ressort-amortisseur de la figure est un cas
particulier d’un système général à 𝑛 degrés de liberté.
Dans la forme la plus générale, les matrices de masse, d’amortissement et de
raideur sont entièrement remplies de coefficient non nuls.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 10
Dynamique des Structures
L’équation (1.5) dénote d’un système de 𝑛 équations différentielles ordinaires
du second ordre couplées. Ces équations peuvent être découplées en utilisant
une procédure appelée analyse modale. Elle requiert les fréquences naturelles
et les modes normaux ou modes naturels du système. Pour déterminer les
fréquences naturelles et les modes normaux, le problème aux valeurs propres
correspondant à la vibration du système non-amorti doit être résolu.
2. Problème aux valeurs propres
La vibration libre du système non-amorti est gouvernée par l’équation (1.6) :
[𝑚]𝑥⃗̈ + [𝑘]𝑥⃗ = 0 (1.6)
En supposant la solution de l’équation (1.6) harmonique, posons là :
𝑥⃗ = 𝑋⃗ sin(𝜔𝑡 + 𝜑) (1.7)
Alors on a :
𝑥⃗̈ = −𝜔2 𝑋⃗ sin (𝜔𝑡 + 𝜑) (1.8)
où 𝑋⃗ est l’amplitude de 𝑥⃗ (𝑡), 𝜑 est l’angle de phase et 𝜔 la fréquence de
vibration.
En substituant les équations (1.7) et (1.8) dans l’équation (1.6), on obtient :
[[𝑘] − 𝜔2 [𝑚]]𝑋⃗ = 0 ⃗⃗. (1.9)
L’équation (1.9) représente un système de 𝑛 équations homogènes en
coefficients inconnus 𝑋1 , 𝑋2 , … 𝑋𝑛 (amplitude de 𝑥 1 , 𝑥 2 , … , 𝑥 𝑛) avec 𝜔2 jouant le
rôle de paramètre.
Pour avoir une solution non triviale du vecteur de coefficients 𝑋⃗, le
déterminant de la matrice de coefficient doit être nulle.
[[𝑘] − 𝜔2 [𝑚]] = 0 (1.10)
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 11
Dynamique des Structures
L’équation (1.10) est une équation polynomiale de degré 𝑛 en 𝜔2 (𝜔2 est
appelé valeur propre) et est appelée équation caractéristique ou équation de
fréquence.
Les racines du polynôme donnent les 𝑛 valeurs propres 𝜔1 2 , 𝜔2 2 , … , 𝜔𝑛 2 . Les
racines carrées positives des valeurs propres donnent les fréquences
naturelles du système 𝜔1 , 𝜔2 , … , 𝜔𝑛 .
Les fréquences naturelles sont habituellement arrangées par ordre croissant,
tel que 𝜔1 ≤ 𝜔2 ≤ ⋯ ≤ 𝜔𝑛 . La plus petite fréquence est référée comme étant la
fréquence fondamentale.
Pour chaque fréquence naturelle 𝜔𝑖 un vecteur non trivial correspondant 𝑋⃗ (𝑖)
peut être obtenu à partir de l’équation (1.9) :
[[𝑘] − 𝜔2 [𝑚]]𝑋⃗ (𝑖) = ⃗0⃗ (1.11)
Le vecteur 𝑋⃗ (𝑖) est appelé vecteur propre, vecteur caractéristique, ou vecteur
modal ou vecteur normal associé à la valeur propre 𝜔𝑖 .
Des lors que l’équation (1.11) représente un système d’équations homogènes,
si 𝑋⃗ (𝑖) est solution de l’équation (1.11) alors 𝑐𝑖 𝑋⃗ (𝑖) aussi est une solution, ou
𝑐𝑖 est une constante arbitraire. Cela indique que la forme d’un mode naturel
est unique mais pas son amplitude. Habituellement une amplitude est
assignée au vecteur propre 𝑋⃗ (𝑖) pour le rendre unique en utilisant une
procédure appelée normalisation. On utilise communément la normalisation
dite normalisation par rapport à la matrice des masses, consistant à poser :
𝑇
𝑋⃗ (𝑖) [𝑚]𝑋⃗ (𝑖) = 1, 𝑖 = 1; 2; … ; 𝑛. (1.12)
( ) 𝑇 ( )
où 𝑋⃗ représente la transposée de 𝑋⃗ .
𝑖 𝑖
3. Orthogonalité des vecteurs modaux
Les vecteurs modaux possèdent une propriété importante appelée
orthogonalité par rapport à la matrice d’inertie (respectivement par rapport à
la matrice de raideur).
En effet, considérons deux valeurs propres distinctes 𝜔𝑖 2 et 𝜔𝑗 2 et les vecteurs
propres correspondants respectivement 𝑋⃗ (𝑖) et 𝑋⃗ (𝑗) . Ces solution vérifient
l’équation (1.9), donc :
[𝑘]𝑋⃗ (𝑖) = 𝜔𝑖 2 [𝑚]𝑋⃗ (𝑖) (1.13)
[𝑘]𝑋⃗ (𝑗) = 𝜔𝑗 2 [𝑚]𝑋⃗ (𝑗) (1.14)
𝑇
Multiplions les deux membres de l’équation (1.13) par 𝑋⃗ (𝑗) et les deux
𝑇
membres de l’équation (1.14) par 𝑋⃗ (𝑖) , on obtient :
𝑇 𝑇
𝑋⃗ (𝑗) [𝑘]𝑋⃗ (𝑖) = 𝜔𝑖 2 𝑋⃗ (𝑗) [𝑚]𝑋⃗ (𝑖) (1.15)
𝑇 𝑇
𝑋 ⃗ ( 𝑖)
[𝑘]𝑋⃗ (𝑗) = 𝜔𝑗 𝑋 2 ⃗ ( 𝑖)
[𝑚]𝑋⃗ (𝑗) (1.16)
En notant que les matrices [𝑚] et [𝑘] sont symétriques, effectuons la
transposée de l’équation (1.16) et retranchons le résultat de l’équation (1.15),
membre à membre. On obtient :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 12
Dynamique des Structures
𝑇
(𝜔𝑖 2 − 𝜔𝑗 2 )𝑋⃗ (𝑗) [𝑚]𝑋⃗ (𝑖) = 0 (1.17)
Du fait que les valeurs propres sont distinctes, 𝜔𝑖 2 ≠ 𝜔𝑗 2 alors l’équation (1.17)
donne :
𝑇
𝑋⃗ (𝑗) [𝑚]𝑋⃗ (𝑖) = 0, 𝑖≠𝑗 (1.18)
En substituant l’équation (1.18) dans l’équation (15), on obtient :
𝑇
𝑋⃗ (𝑗) [𝑘]𝑋⃗ (𝑖) = 0, 𝑖≠𝑗 (1.19)
Les équations (1.18) et (1.19) illustrent la propriété d’orthogonalité des
vecteurs propres par rapport aux matrices d’inertie et de raideur
respectivement. Lorsque 𝑖 = 𝑗, les équations (1.15) et (1.16) deviennent :
𝑇 𝑇
𝑋⃗ (𝑖) [𝑘]𝑋⃗ (𝑖) = 𝜔𝑖 2 𝑋⃗ (𝑖) [𝑚]𝑋⃗ (𝑖) (1.20)
Si les vecteurs propres sont normalisés suivant l’équation (1.12), alors
l’équation (1.20) donne :
𝑇
𝑋⃗ (𝑖) [𝑘]𝑋⃗ (𝑖) = 𝜔𝑖 2 (1.21)
En considérant tous les vecteurs propres, les équations (1.12) et (1.21)
s’écrivent :
1 0
𝑇
[𝑋] [𝑚][𝑋] = [𝐼 ] = [ 1 ] (1.22)
⋱
0 1
2
𝜔1 0
2
[𝑋]𝑇 [𝑘][𝑋] = [𝐼 ] = [ 𝜔2 ] (1.23)
⋱
0 𝜔𝑛 2
où la matrice carrée d’ordre 𝑛 [𝑋] appelée matrice modale, contient les
vecteurs propres 𝑋⃗ (1) ; 𝑋⃗ (2) ; … ; 𝑋⃗ (𝑛) comme colonnes :
[𝑋] = [𝑋⃗ (1) 𝑋⃗ (2) … 𝑋⃗ (𝑛) ] (1.24)
4. Analyse des vibrations pour un système en vibrations libres et non
amorties
La vibration libre d’un système libre non amorti à 𝑛 degrés de liberté est régie
par les équations :
[𝑚]𝑥⃗̈ + [𝑘]𝑥⃗ = 0⃗⃗ (1.25)
Les 𝑛 équations différentielles homogènes couplées du second ordre peuvent
être découpler en utilisant l’analyse modale. Dans l’analyse la solution, 𝑥⃗, est
exprimée comme une superposition des modes normaux :
𝑛
𝑥⃗ = ∑ 𝑞𝑖 (𝑡)𝑋⃗ (𝑖) = [𝑋]𝑞⃗(𝑡) (1.26)
𝑖 =1
où [𝑋] est la matrice modale, 𝑞𝑖 (𝑡) sont des fonctions inconnues du temps
appelées coordonnées modales ou coordonnées généralisées, le vecteur 𝑞⃗(𝑡) est
le vecteur des coordonnées modales :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 13
Dynamique des Structures
𝑞1 (𝑡)
𝑞 (𝑡)
𝑞⃗(𝑡) = { 2 } (1.27)
⋮
𝑞𝑛 (𝑡)
L’équation (1.26) représente le Théorème d’expansion et est basée sur le fait
que les vecteurs propres sont orthogonaux et forment une base de l’espace
de dimension 𝑛. Cela implique que tout vecteur, tel que 𝑥⃗ (𝑡), dans l’espace de
dimension 𝑛 peut être générer par une combinaison linéaire de vecteurs
linéairement indépendants, tels que les vecteurs propres 𝑋⃗ (𝑖) , 𝑖 = 1; 2; … ; 𝑛.
En substituant l’équation (1.26) dans l’équation (1.25), on obtient :
[𝑚][𝑋]𝑞⃗̈ + [𝑘][𝑋]𝑞⃗ = ⃗0⃗ (1.28)
En multipliant, à gauche, l’équation (1.28) par 𝑋 [ ] 𝑇
[𝑋]𝑇 [𝑚][𝑋]𝑞⃗̈ + [𝑋]𝑇 [𝑘][𝑋]𝑞⃗ = 0 ⃗⃗ (1.29)
En vertu des équations (1.22) et (1.23), l’équation (1.29) se simplifie en :
𝑞⃗̈ + 𝜔𝑖 2 𝑞⃗ = 0
⃗⃗ (1.30)
L’équation (1.30) représente un système de 𝑛 équations différentielles
homogènes découplées :
𝑑 2 𝑞𝑖 (𝑡)
+ 𝜔𝑖 2 𝑞𝑖 (𝑡) = 0, 𝑖 = 1, 2, … , 𝑛 (1.31)
𝑑𝑡 2
Si les conditions initiales du système sont données par :
𝑥 1,0
𝑥 2,0
𝑥⃗ (𝑡 = 0) = 𝑥⃗ 0 = { ⋮ } (1.32)
𝑥 𝑛,0
𝑥̇ 1,0
𝑥̇
𝑥⃗̇ (𝑡 = 0) = 𝑥⃗ 0̇ = { 2,0 } (1.33)
⋮
𝑥̇ 𝑛,0
Les conditions initiales sur 𝑞⃗(𝑡) peuvent être déterminées comme suit :
Multiplier, à gauche, l’équation (1.26) par [𝑋]𝑇 [𝑚] et utiliser l’équation (1.22)
pour obtenir :
𝑞⃗(𝑡) = [𝑋]𝑇 [𝑚]𝑥⃗ (𝑡) (1.34)
Alors :
𝑞 1 ( 0)
𝑞 ( 0)
{ 2 } = 𝑞⃗(0) = [𝑋]𝑇 [𝑚]𝑥⃗ 0 (1.35)
⋮
𝑞 𝑛 ( 0)
𝑞̇ 1 (0)
𝑞̇ (0)
{ 2 } = 𝑞⃗̇(0) = [𝑋]𝑇 [𝑚]𝑥⃗̇ 0 (1.36)
⋮
𝑞̇ 𝑛 (0)
La solution de l’équation (1.31) peut être exprimé par :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 14
Dynamique des Structures
𝑞̇ 𝑖 (0)
𝑞𝑖 (𝑡) = 𝑞𝑖 (0) cos 𝜔𝑖 𝑡 + sin 𝜔𝑖 𝑡 , 𝑖 = 1, 2, … , 𝑛 (1.37)
𝜔𝑖
où 𝑞𝑖 (0) et 𝑞̇ 𝑖 (0) sont donnés par les équations (1.35) et (1.36) tel que :
𝑇
𝑞⃗(0) = 𝑋 (𝑖) [𝑚]𝑥⃗0 (1.38)
̇𝑞⃗(0) = 𝑋 (𝑖) 𝑇 [𝑚]𝑥⃗̇ (1.39)
0
Une fois les 𝑞𝑖 (𝑡) déterminés, la solution 𝑥⃗ (𝑡) de la vibration libre peut être
trouver en utilisant l’équation (1.26).
5. Analyse des vibrations forcées d’un système non amorti
L’équation du mouvement peut être exprimer par l’équation (40) :
[𝑚]𝑥⃗̈ + [𝑘]𝑥⃗ = 𝑓⃗ (1.40)
Les valeurs propres 𝜔𝑖 et les vecteurs propres correspondants 𝑋⃗ , 𝑖 =
2 (𝑖)
1, 2, … , 𝑛 sont supposés connus. On suppose également que la solution de
l’équation (1.40) est donnée par une combinaison linéaire des vecteurs
propres telle que :
𝑛
𝑥⃗ = ∑ 𝑞𝑖 (𝑡)𝑋⃗ (𝑖) = [𝑋]𝑞⃗(𝑡) (1.41)
𝑖 =1
où 𝑞𝑖 (𝑡) représente les coordonnées modales et [𝑋] la matrice modale. En
substituant l’équation (1.41) dans l’équation (1.40) et en multipliant à gauche
par [𝑋]𝑇 , on obtient :
[𝑋]𝑇 [𝑚][𝑋]𝑞⃗̈ + [𝑋]𝑇 [𝑘][𝑋]𝑞⃗ = 𝑓⃗ (1.42)
En utilisant les équations (1.22) et (1.23), l’équation (1.42) peut être exprimée
par :
𝑞⃗̈ + [𝜔𝑖 2 ]𝑞⃗ = 𝑄⃗⃗ (1.43)
Ou 𝑄⃗⃗ est appelé vecteur des forces modales (ou forces généralisées), donné
par :
𝑄⃗⃗ = [𝑋]𝑇 𝑓⃗(𝑡) (1.44)
Les 𝑛 équations différentielles découplées indiquées par Eq. (1.43) peuvent
être exprimées sous forme scalaire par Eq. (1.45) :
𝑑 2 𝑞𝑖 (𝑡)
+ 𝜔𝑖 2 𝑞𝑖 (𝑡) = 𝑄𝑖 (𝑡), 𝑖 = 1, 2, … , 𝑛 (1.45)
𝑑𝑡 2
où
𝑇
𝑄𝑖 (𝑡) = 𝑋 (𝑖) 𝑓⃗(𝑡), 𝑖 = 1, 2, … , 𝑛 (1.46)
Chacune des équations dans Eq. (1.45) peut être considérée comme
l’équation du mouvement d’un système en vibration forcée non amorti à un
seul degré de liberté. Des lors, on peut écrire :
𝑡
𝑞𝑖 (𝑡) = ∫ 𝑄𝑖 (𝜏)ℎ(𝑡 − 𝜏)𝑑𝜏 + 𝑔 (𝑡)𝑞𝑖,0 + ℎ(𝑡)𝑞̇ 𝑖,0 (1.47)
0
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 15
Dynamique des Structures
où
1
ℎ(𝑡) = sin 𝜔𝑖 (𝑡) (1.48)
𝜔𝑖
𝑔(𝑡) = cos 𝜔𝑖 (𝑡) (1.49)
Les valeurs initiales 𝑞𝑖,0 et 𝑞̇ 𝑖,0 peuvent être déterminées à partir des conditions
𝑥 et 𝑥⃗̇ en utilisant les équations (1.35) et (1.36).
initiales ⃗⃗⃗⃗⃗
0 0
6. Analyse des vibrations forcées d’un système avec un amortissement
proportionnel
Dans l’amortissement proportionnel, la matrice d’amortissement [𝑐 ] dans
l’équation (1.5) peut être exprimée par une combinaison linéaire des matrices
d’inertie et de raideur telle que :
[𝑐 ] = 𝛼[𝑚] + 𝛽[𝑘] (1.50)
où 𝛼 et 𝛽 sont des constantes connues. En substituant l’équation (1.50) dans
l’équation (1.5), on obtient :
[𝑚]𝑥⃗̈ + (𝛼[𝑚] + 𝛽[𝑘])𝑥⃗̇ + [𝑘]𝑥⃗ = 𝑓⃗ (1.51)
Comme indiqué plutôt, dans l’analyse modale, la solution est supposée être
de la forme :
𝑥⃗ = [𝑋]𝑞⃗(𝑡) (1.52)
En substituant l’équation (1.52) dans l’équation (1.51) et en multipliant à
gauche par [𝑋]𝑇 , on obtient :
[𝑋]𝑇 [𝑚][𝑋]𝑞⃗̈ + 𝛼[𝑋]𝑇 [𝑚][𝑋]𝑞⃗̇ + 𝛽[𝑋]𝑇 [𝑘][𝑋]𝑞⃗̇ + [𝑋]𝑇 [𝑘][𝑋]𝑞⃗ = [𝑋]𝑇 𝑓⃗ (1.53)
Lorsqu’on utilise les équations (1.22) et (1.23), on obtient :
𝑞⃗̈ + (𝛼[𝐼 ] + 𝛽[𝜔𝑖 2 ])𝑞⃗̇ + [𝜔𝑖 2 ]𝑞⃗ = 𝑄⃗⃗ (1.54)
où
𝑄⃗⃗ = [𝑋]𝑇 𝑓⃗(𝑡) (1.55)
Posons :
𝛼 + 𝛽𝜔𝑖 2 = 2𝛾𝑖 𝜔𝑖 , 𝑖 = 1, 2, … , 𝑛 (1.56)
où 𝛾𝑖 est appelé le facteur d’amortissement visqueux dans le i-ème mode, Eq.
(1.54) peut être réécrite sous forme scalaire telle que :
𝑑 2 𝑞𝑖 (𝑡) 𝑑 𝑞𝑖 (𝑡)
+ 2𝛾 𝜔
𝑖 𝑖 + 𝜔𝑖 2 𝑞𝑖 (𝑡) = 𝑄𝑖 (𝑡), 𝑖 = 1, 2, … , 𝑛 (1.57)
𝑑𝑡 2 𝑑𝑡
Chacune des équations dans Eq. (1.57) peut être considérée comme
l’équation du mouvement d’un système à un degré de liberté à
amortissement visqueux, dont la solution est donnée par :
𝑡
𝑞𝑖 (𝑡) = ∫ 𝑄𝑖 (𝜏)ℎ(𝑡 − 𝜏)𝑑𝜏 + 𝑔 (𝑡)𝑞𝑖,0 + ℎ(𝑡)𝑞̇ 𝑖,0 (1.58)
0
1 −𝛾 𝜔
ℎ (𝑡) = 𝑒 𝑖 𝑖 sin 𝜔𝑑𝑖 (𝑡) (1.59)
𝜔𝑑𝑖
𝛾𝜔
𝑔(𝑡) = 𝑒 −𝛾𝑖 𝜔𝑖 (cos 𝜔𝑑𝑖 (𝑡) + 𝑖 𝑖 sin 𝜔𝑑𝑖 (𝑡)) (1.60)
𝜔𝑑𝑖
Et 𝜔𝑑𝑖 est la i-ème fréquence de vibration amortie :
𝜔𝑑𝑖 = √1 − 𝛾𝑖 2 𝜔𝑖 (1.61)
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 16
Dynamique des Structures
7. Analyse des vibrations forcées d’un système avec un amortissement
visqueux général
L’équation du mouvement d’un système à 𝑛 degrés de liberté avec un
amortissement quelconque peut être exprimer par l’équation (1.5), reprise :
[𝑚]𝑥⃗̈ + [𝑐 ]𝑥⃗̇ + [𝑘]𝑥⃗ = 𝑓⃗ (1.65)
Dans ce cas, la matrice modale ne peut pas diagonaliser la matrice
d’amortissement, et une solution analytique n’est pas possible.
Posons :
𝑥⃗̇ = 𝑥⃗̇ (1.66)
𝑥⃗̈ = −[𝑚]−1 [𝑐 ]𝑥⃗̇ − [𝑚]−1 [𝑘]𝑥⃗ + [𝑚]−1 𝑓⃗ (1.67)
Définissons le vecteur 𝑦⃗(𝑡), de dimension 2𝑛, par :
𝑥⃗
𝑦⃗ (𝑡) = { ̇ } (1.68)
𝑥⃗
Les équations (1.66) et (1.67) peuvent être exprimées par :
𝑦⃗̇(𝑡) = [𝐴]𝑦⃗(𝑡) + [𝐵 ]𝑓⃗(𝑡) (1.69)
où les matrices [𝐴] et [𝐵 ], d’ordre 2𝑛 × 2𝑛, sont définies par :
[ 0] [𝐼 ]
[𝐴] = [ −1 [ ] ] (1.70)
−𝑚[ ] 𝑘 − 𝑚]−1 [𝑐 ]
[
[ 0]
[𝐵] = [ ] (1.71)
[𝑚]−1
Analyse modale
Pour l’analyse modale, nous considérons d’abord le problème de la vibration
libre avec 𝑓⃗ = ⃗0⃗. L’équation (1.69) se réduit en :
𝑦⃗̇(𝑡) = [𝐴]𝑦⃗(𝑡) (1.72)
Eq. (1.72) représente un système de 2𝑛 équations différentielles ordinaires du
premier ordre dont on peut supposer la solution sous la forme :
𝑦⃗(𝑡) = 𝑌⃗⃗𝑒 𝜆𝑡 (1.72)
où 𝑌⃗⃗ est un vecteur constant et 𝜆 est un scalaire constant. En substituant Eq.
(1.72) dans Eq. (1.71) et en simplifiant le terme 𝑒 𝜆𝑡 , on obtient :
[𝐴]𝑌⃗⃗ = 𝜆𝑌⃗⃗ (1.73)
L’équation (1.73) peut être un problème algébrique de valeur propre standard
avec une matrice réelle non symétrique [𝐴] dont la solution donne les valeurs
propres 𝜆 𝑖 et les vecteurs propres 𝑌⃗⃗ (𝑖) correspondants 𝑖 = 1, 2, … , 𝑛. Ces valeurs
propres et vecteurs propres peuvent être réels ou complexes. Si 𝜆 𝑖 est une
valeur propre complexe, on peut démontrer que son conjugué 𝜆̅𝑖 sera aussi
une valeur propre. Aussi les vecteurs propres 𝑌⃗⃗ (𝑖) et ̅̅̅̅̅𝑌⃗⃗ (𝑖) correspondants à
𝜆 𝑖 et 𝜆̅𝑖 seront complexes, conjugués l’un de l’autre.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 17
Dynamique des Structures
Les vecteurs propres 𝑌⃗⃗ (𝑖) correspondants au problème aux valeurs propres (Eq.
(1.73)) sont appelés les vecteurs propres à droite de la matrice [𝐴]. Les vecteurs
propres correspondants à la transposée de la matrice [𝐴] sont appelés les
vecteurs propres à gauche de la matrice [𝐴]. Donc les vecteurs propres à
gauche de la matrice [𝐴] sont obtenus en résolvant le problème aux valeurs
propres :
[𝐴]𝑇 𝑍⃗ = 𝜆𝑍⃗ (1.74)
Dès lors que les déterminants des matrices [𝐴] et [𝐴] sont égaux, les 𝑇
équations caractéristiques correspondants aux équations (1.73) et (1.74)
seront identiques :
|[𝐴] − 𝜆[𝐼 ]| ≡ |[𝐴]𝑇 − 𝜆[𝐼 ]| = 0 (1.75)
Alors les valeurs propres des équations (1.73) et (1.74) seront identiques.
Cependant les vecteurs propres de [𝐴] et de [𝐴]𝑇 seront différents. Pour trouver
la relation entre 𝑌⃗⃗ (𝑖), 𝑖 = 1, 2, … , 2𝑛 et 𝑍⃗ (𝑗) , 𝑗 = 1, 2, … , 2𝑛, écrivons les problèmes
aux valeurs propres correspondant à 𝑌⃗⃗ (𝑖) et 𝑍⃗ (𝑗) comme suit :
[𝐴]𝑌⃗⃗ (𝑖) = 𝜆 𝑖 𝑌⃗⃗ (𝑖) et [𝐴]𝑇 𝑍⃗ (𝑗) = 𝜆𝑗 𝑍⃗ (𝑗) (1.76)
où
𝑇 𝑇
𝑍⃗ (𝑗) [𝐴] = 𝜆𝑗 𝑍⃗ (𝑗) (1.77)
𝑇
Multiplions à gauche l’équation (1.76) par 𝑍⃗ (𝑗) et multiplions à droite
l’équation (1.77) par 𝑌⃗⃗ (𝑖) :
𝑇 𝑇
𝑍⃗ (𝑗) [𝐴]𝑌⃗⃗ (𝑖) = 𝜆 𝑖 𝑍⃗ (𝑗) 𝑌⃗⃗ (𝑖) (1.78)
𝑇 𝑇
𝑍⃗ (𝑗) [𝐴]𝑌⃗⃗ (𝑖) = 𝜆𝑗 𝑍⃗ (𝑗) 𝑌⃗⃗ (𝑖) (1.79)
Effectuons la soustraction de l’équation (1.79) de l’équation (1.78), on
obtient :
𝑇
(𝜆 𝑖 − 𝜆𝑗 )𝑍⃗ (𝑗) 𝑌⃗⃗ (𝑖) = 0 (1.80)
En supposant que 𝜆 𝑖 ≠ 𝜆𝑗 , il vient que :
𝑇
𝑍⃗ (𝑗) 𝑌
⃗⃗ (𝑖) = 0, 𝑖, 𝑗 = 1, 2, … , 2𝑛 (1.81)
Ce qui montre que le i-ème vecteur propre à droite de la matrice [𝐴] est
orthogonal au j-ème vecteur à gauche de la matrice [𝐴], avec les valeurs
propres 𝜆 𝑖 et 𝜆𝑗 différentes. En substituant Eq. (1.81) dans Eq. (1.78) ou Eq.
(1.79), on trouve que :
𝑇
𝑍⃗ (𝑗) [𝐴]𝑌
⃗⃗ (𝑖) = 0, 𝑖, 𝑗 = 1, 2, … , 2𝑛 (1.82)
Lorsque les vecteurs propres à droite et à gauche sont normalisés tels que :
𝑇
𝑍⃗ (𝑖) 𝑌
⃗⃗ (𝑖) = 1 (1.83)
L’équation (1.82) donne :
𝑇
𝑍⃗ (𝑖) [𝐴]𝑌
⃗⃗ (𝑖) = 𝜆 𝑖 , 𝑖 = 1, 2, … , 2𝑛 (84)
Les équations (1.83) et (1.84) peuvent être exprimées sous forme matricielle
par :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 18
Dynamique des Structures
[𝑍]𝑇 [𝑌] = [𝐼 ] (1.85)
[𝑍]𝑇 [𝐴][𝑌] = [𝜆𝑖 ] (1.86)
où les matrices [𝑍] et [𝑌] sont définies par :
[𝑌] ≡ [𝑌⃗⃗ (1) 𝑌⃗⃗ (2) ⋯ 𝑌⃗⃗ (2𝑛) ] (1.87)
[𝑍] ≡ [𝑍⃗ (1) 𝑍⃗ (2) ⋯ 𝑍⃗ (2𝑛) ] (1.88)
Et la matrice diagonale des valeurs propres est définie par :
𝜆1 0
𝜆 2
[𝜆 𝑖 ] = [ ] (1.89)
⋱
0 𝜆𝑛
Dans l’analyse modale la solution de l’équation d’état, Eq. (1.69), est supposée
être une combinaison linéaire des vecteurs propres :
2𝑛
𝑦⃗(𝑡) = ∑ 𝑞𝑖 (𝑡)𝑌⃗⃗ (𝑖) = [𝑌]𝑞⃗(𝑡) (1.90)
𝑖=1
où 𝑞𝑖 (𝑡), 𝑖 = 1, 2, … , 2𝑛, sont les coordonnées modales et 𝑞⃗(𝑡) est le vecteur des
coordonnées modales :
𝑞1 (𝑡)
𝑞 (𝑡)
𝑞⃗(𝑡) = { 2 } (1.91)
⋮
𝑞2𝑛 (𝑡)
En substituant Eq. (1.90) dans Eq. (1.69) et multipliant à gauche par [𝑍]𝑇 , on
obtient :
[𝑍]𝑇 [𝑌]𝑞⃗̇(𝑡) = [𝑍]𝑇 [𝐴][𝑌]𝑞⃗(𝑡) + [𝑍]𝑇 [𝐵 ]𝑓⃗(𝑡) (1.92)
Au regard des équations (1.85) et (1.86), l’équation (1.92) se réduit à :
𝑞⃗̇(𝑡) = [𝜆 𝑖 ] 𝑞⃗(𝑡) + 𝑄⃗⃗(𝑡) (1.93)
L’équation (1.93) s’écrit sous forme scalaire par :
𝑑𝑞𝑖 (𝑡)
= 𝜆 𝑖 𝑞𝑖 (𝑡) + 𝑄𝑖 (𝑡), 𝑖 = 1, 2, … , 2𝑛 (1.94)
𝑑𝑡
où le vecteur des forces modales est donné par :
𝑄⃗⃗ (𝑡) = [𝑍]𝑇 [𝐵 ]𝑓⃗(𝑡) (1.95)
Et la i-ème force modale est donnée par :
𝑇
𝑄𝑖 (𝑡) = 𝑍⃗ (𝑖) [𝐵]𝑓⃗(𝑡), 𝑖 = 1, 2, … , 2𝑛 (1.96)
Les solutions des équations différentielles ordinaires du premier ordre, Eq.
(1.94), peuvent être exprimées par :
𝑡
𝑞𝑖 (𝑡) = ∫ 𝑒 𝜆𝑖 (𝑡−𝜏) 𝑄𝑖 (𝜏)𝑑𝜏 + 𝑒 𝜆𝑖 𝑡 𝑞𝑖 (0), 𝑖 = 1, 2, … , 2𝑛 (1.97)
0
Sous forme matricielle, on a :
𝑡
⃗⃗⃗(𝑡) = ∫ 𝑒 [𝜆𝑖 ](𝑡−𝜏) 𝑄⃗⃗𝑖 (𝜏)𝑑𝜏 + 𝑒 [𝜆𝑖 ]𝑡 𝑞⃗(0),
𝑞 𝑖 = 1, 2, … , 2𝑛 (1.98)
0
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 19
Dynamique des Structures
où 𝑞⃗(0) désigne la valeur initiale de 𝑞
⃗⃗⃗(𝑡). Pour déterminer 𝑞
⃗⃗⃗(0), multiplions à
𝑇
gauche l’équation (1.90) par 𝑍⃗ (𝑖) pour obtenir :
𝑇 𝑇
𝑍⃗ (𝑖) 𝑦⃗(𝑡) = 𝑍⃗ (𝑖) [𝑌]𝑞⃗(𝑡) (1.99)
Au regard des relations d’orthogonalité, Eq. (1.81), l’équation (1.99) donne :
𝑇
𝑞𝑖 (𝑡) = 𝑍⃗ (𝑖) 𝑦⃗(𝑡), 𝑖 = 1, 2, … , 2𝑛 (1.100)
En posant 𝑡 = 0 dans Eq. (1.100), la valeur initiale de 𝑞𝑖 (0) peut être trouvée
telle que :
𝑇
𝑞𝑖 (0) = 𝑍⃗ (𝑖) 𝑦⃗(0), 𝑖 = 1, 2, … , 2𝑛 (1.101)
Finalement, la solution de l’équation (1.69) peut être exprimée, en utilisant
Eq. 1.90 et Eq. 1.98, par :
𝑡
⃗⃗⃗(𝑡) = ∫ [𝑌]𝑒 [𝜆𝑖 ](𝑡−𝜏) 𝑄⃗⃗𝑖 (𝜏)𝑑𝜏 + [𝑌]𝑒 [𝜆𝑖]𝑡 𝑞⃗(0),
𝑦 𝑖 = 1, 2, … , 2𝑛 (1.102)
0
Détermination des équations du mouvement en utilisant les équations de
Lagrange
Equations de Lagrange :
𝑑 𝜕𝐸𝑐 𝜕𝐸𝑐 𝜕𝐸𝑝
[ ]− + = 𝑄𝑖
𝑑𝑡 𝜕𝑥̇ 𝑖 𝜕𝑥 𝑖 𝜕𝑥 𝑖
1 1 1
Energie cinétique : 𝐸𝑐 = 2 𝑚𝑥̇ 1 2 + 2 𝑚𝑥̇ 2 2 + 2 𝑚𝑥̇ 3 2
1 1 1 1
Energie potentielle : 𝐸𝑝 = 𝑘𝑥 1 2 + 𝑘(𝑥 2 − 𝑥 1 )2 + 𝑘(𝑥 2 − 𝑥 3 )2 + 𝑘𝑥 3 2
2 2 2 2
Forces non conservatives : 𝑄𝑖 = 𝑓𝑖
Pour 𝑖 = 1 on a :
𝑑 𝜕𝐸𝑐 𝜕𝐸𝑐 𝜕𝐸𝑝
[ ]− + = 𝑓1 (𝑡)
𝑑𝑡 𝜕𝑥̇ 1 𝜕𝑥 1 𝜕𝑥 1
1
Ce qui nous donne : 𝑚𝑥̈ 1 + 𝑘𝑥 1 − 2 𝑘(𝑥 2 − 𝑥 1 ) = 𝑓1 (𝑡)
De même pour 𝑖 = 2 et 3, on obtient :
𝑚𝑥̈ 2 + 𝑘 (𝑥 2 − 𝑥 1 ) + 𝑘(𝑥 2 − 𝑥 3 ) = 𝑓2 (𝑡)
𝑚𝑥̈ 3 + 𝑘𝑥 3 + 𝑘(𝑥 3 − 𝑥 2 ) = 𝑓3 (𝑡)
On obtient donc les mêmes équations.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 20
Dynamique des Structures
Chapitre 2 : Vibrations des Systèmes Continus
I. Vibrations des cordes
Il est bien connu que la plupart des instruments de musique importants, y
compris le violon et la guitare, impliquent des cordes dont les fréquences
naturelles et les formes de mode jouent un rôle important dans leur
performance. Les caractéristiques de nombreux systèmes d'ingénierie, tels
que les haubans, les lignes de transmission électrique, les câbles et les
courroies utilisés dans les machines et la fabrication de fils, peuvent être
dérivées d'une étude de la dynamique des cordes tendues. La vibration
transversale libre et forcée des cordes est étudiée dans ce chapitre. Comme
on le verra dans les chapitres suivants, l'équation régissant la vibration
transversale des cordes sera de la même forme que les équations de
mouvement des vibrations longitudinales des barres et des vibrations de
torsion des arbres.
1. Equation du mouvement :
a. Approche Equilibre :
Considérons une corde de longueur 𝑙 fortement tendue soumise à une force
transversale répartie 𝑓(𝑥, 𝑡) par unité de longueur (voir Figure 1).
La corde est supposée être supportée aux extrémités par deux ressorts de
raideur 𝑘1 et 𝑘2 .
Figure 8 : Corde Vibrante
En faisant l’hypothèse de petits déplacements transversaux 𝑤(𝑥, 𝑡),
appliquons la seconde loi de Newton dans le mouvement d’un élément de
longueur 𝑑𝑥 (voir Figure 2) dans la direction 𝑧, on a, si 𝑇 est la tension sur la
corde, 𝜌 la masse par unité de longueur et 𝜃 l’angle formé par la corde fléchie
avec l’axe 𝑥:
𝜕2 𝑤(𝑥, 𝑡) (2.1)
(𝑇 + 𝑑𝑇 ) 𝑠𝑖𝑛 (𝜃 + 𝑑𝜃) + 𝑓(𝑥, 𝑡)𝑑𝑥 − 𝑇 𝑠𝑖𝑛 𝜃 = 𝜌 𝑑𝑥
𝜕𝑡 2
On peut écrire l’équation (1) plus simplement :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 21
Dynamique des Structures
Figure 9 : Elément différentiel de corde
𝜕2 𝑤 (2.2)
(𝑇 + 𝑑𝑇 ) sin (𝜃 + 𝑑𝜃) + 𝑓 𝑑𝑥 − 𝑇 sin 𝜃 = 𝜌 𝑑𝑥
𝜕𝑡 2
Notons que :
𝜕𝑇 𝜕𝑤(𝑥, 𝑡)
𝑑𝑇 = 𝑑𝑥 ; sin 𝜃 ≈ tan 𝜃 =
𝜕𝑥 𝜕𝑥
𝜕𝑤 𝜕2 𝑤
sin(𝜃 + 𝑑𝜃) ≈ tan(𝜃 + 𝑑𝜃) = + 𝑑𝑥
𝜕𝑥 𝜕𝑥 2
L’équation (2.2) devient alors :
𝜕 𝜕𝑤(𝑥, 𝑡) 𝜕2 𝑤 (𝑥, 𝑡) (2.3)
[𝑇 ( ) (
] + 𝑓 𝑥, 𝑡 = 𝜌 𝑥 )
𝜕𝑥 𝜕𝑥 𝜕𝑡 2
Si la corde est uniforme et la tension constante, l’équation (2.3) prend la
forme :
𝜕2 𝑤(𝑥, 𝑡) 𝜕2 𝑤(𝑥, 𝑡)
𝑇 ( )
+ 𝑓 𝑥, 𝑡 = 𝜌(𝑥)
𝜕𝑥 2 𝜕𝑡 2
Pour une fibration libre, 𝑓(𝑥, 𝑡) = 0, l’équation devient :
𝜕2 𝑤(𝑥, 𝑡) 𝜕2 𝑤(𝑥, 𝑡)
𝑇 = 𝜌(𝑥)
𝜕𝑥 2 𝜕𝑡 2
Posons 𝑐 = √𝑇/𝜌 , on peut alors écrire :
𝝏𝟐 𝒘(𝒙, 𝒕) 𝝏𝟐 𝒘(𝒙, 𝒕) 𝝏𝟐 𝒘 𝝏𝟐 𝒘 (2.4)
𝒄𝟐 𝟐
= 𝟐
𝐬𝐨𝐢𝐭 𝒄𝟐
𝟐
= 𝟐
𝝏𝒙 𝝏𝒕 𝝏𝒙 𝝏𝒕
Equation du mouvement d’une corde tendue soumise à une vibration
transversale
2. Conditions Initiales (C.I) et Conditions aux Limites (C.L)
L’équation du mouvement de la corde est une équation différentielle d’ordre
2 en 𝑥 et en 𝑡. Dès lors deux conditions initiales et deux conditions aux limites
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 22
Dynamique des Structures
sont nécessaires pour trouver la solution 𝑤(𝑥, 𝑡). Si la corde est donnée avec
une déflexion initiale 𝑤0 (𝑥 ) et une vitesse initiale 𝑤̇0 (𝑥 ) alors les conditions
initiales peuvent être énoncées comme suit :
𝜕𝑤
𝑤(𝑥, 𝑡 = 0) = 𝑤0 (𝑥 ) et (𝑥, 𝑡 = 0) = 𝑤̇0 (𝑥 ).
𝜕𝑡
Le tableau 1 donne des indications supplémentaires sur les conditions aux
limites.
Tableau 1 : Conditions aux limites en fonction du support de la corde
Support de fixation de la corde Conditions aux limites à
satisfaire
1. Les deux extrémités de la corde sont fixées
𝑤 (0, 𝑡) = 0
𝑤(𝑙, 𝑡) = 0
2. Les deux extrémités de la corde sont libres
𝜕𝑤
(0, 𝑡) = 0
𝜕𝑡
𝜕𝑤
(𝑙, 𝑡) = 0
𝜕𝑡
3. Les deux extrémités de la corde sont fixées à
des masses
𝜕2 𝑤 𝜕𝑤
𝑚 1 2 (0, 𝑡) = 𝑇 (0, 𝑡)
𝜕𝑡 𝜕𝑥
𝜕2 𝑤 𝜕𝑤
−𝑚 2 2 (𝑙, 𝑡) = 𝑇 (𝑙, 𝑡)
𝜕𝑡 𝜕𝑥
4. Les deux extrémités attachées à des ressorts
𝜕𝑤
𝑘1 𝑤(0, 𝑡) = 𝑇 (0, 𝑡)
𝜕𝑥
𝜕𝑤
−𝑘1 𝑤(𝑙, 𝑡) = 𝑇 (𝑙, 𝑡)
𝜕𝑥
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 23
Dynamique des Structures
Si la corde est fixée à 𝑥 = 0, alors le déplacement est nul, les conditions aux
limites seront :
𝑤 (𝑥 = 0, 𝑡) = 0, 𝑡 ≥ 0 .
Si la corde est connectée à une broche dont le mouvement transversal est libre
alors l’autre bout de la corde ne pourra supporter aucune force transversale,
la condition à la limite sera :
𝜕𝑤
𝑇 (𝑥 ) (𝑥, 𝑡) = 0.
𝜕𝑥
Si la force axiale est constante et l’extrémité à 𝑥 = 0 est libre, on a la
condition :
𝜕𝑤
(0, 𝑡) = 0; 𝑡 ≥ 0.
𝜕𝑥
3. Vibration Libre d’une corde Infinie :
Considérons une corde de longueur infinie, l’équation de la vibration libre de
la corde peut être résolue de trois manières différentes.
a. Solution ondulatoire :
Nous pouvons chercher la solution sous forme d’onde progressive :
𝑤 (𝑥, 𝑡) = 𝐹1 (𝑥 − 𝑐𝑡) + 𝐹2 (𝑥 + 𝑐𝑡)
𝐹1 et 𝐹2 Sont des fonctions arbitraires de 𝑥 − 𝑐𝑡 et 𝑥 + 𝑐𝑡 respectivement. La
solution obtenue est appelée solution de D’Alembert. On vérifie la validité de
cette solution en effectuant les dérivées d’ordre 2 de la solution qui sont :
𝜕2 𝑤
= 𝐹1 ′′ (𝑥 − 𝑐𝑡) + 𝐹2 ′′ (𝑥 + 𝑐𝑡)
𝜕𝑥 2
𝜕2 𝑤
= 𝑐 2 𝐹1 ′′ (𝑥 − 𝑐𝑡) + 𝑐 2 𝐹2 ′′ (𝑥 + 𝑐𝑡)
𝜕𝑡 2
En substituant ces expressions dans l’équation, on trouve que cette dernière
est satisfaite.
Les fonctions 𝐹1 (𝑥 − 𝑐𝑡) et 𝐹2 (𝑥 + 𝑐𝑡) dénotent les ondes qui se propagent dans
les directions positive et négative de l’axe 𝑥, respectivement, avec une vitesse
𝑐.
Les fonctions 𝐹1 et 𝐹2 peuvent être déterminées à partir des conditions
initiales de la corde :
𝐹1 (𝑥 ) + 𝐹2 (𝑥 ) = 𝑤0 (𝑥 )
−𝑐𝐹1 ′ (𝑥 ) + 𝑐𝐹2 ′ (𝑥 ) = 𝑤̇0 (𝑥 )
En intégrant cette dernière expression on a :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 24
Dynamique des Structures
1 𝑥
−𝐹1 (𝑥 ) + 𝐹2 (𝑥 ) = ∫ 𝑤̇ (𝑢)𝑑𝑢
𝑐 𝑥0 0
Où 𝑥 0 est une constante. On peut alors former le système suivant :
𝐹1 (𝑥 ) + 𝐹2 (𝑥 ) = 𝑤0 (𝑥 )
{ 1 𝑥
−𝐹1 𝑥 + 𝐹2 𝑥 = ∫ 𝑤̇0 (𝑢)𝑑𝑢
( ) ( )
𝑐 𝑥0
Ce qui nous donne finalement :
1 1 𝑥 1 1 𝑥
𝐹1 (𝑥 ) = [𝑤0 (𝑥 ) − ∫ 𝑤̇0 (𝑢)𝑑𝑢 ] et 𝐹2 (𝑥 ) = [𝑤0 (𝑥 ) − ∫ 𝑤̇0 (𝑢)𝑑𝑢]
2 𝑐 𝑥0 2 𝑐 𝑥0
En remplaçant respectivement 𝑥 par 𝑥 − 𝑐𝑡 et par 𝑥 + 𝑐𝑡, nous pouvons
exprimer la solution ondulatoire de la corde :
𝑤 (𝑥, 𝑡) = 𝐹1 (𝑥 − 𝑐𝑡) + 𝐹2 (𝑥 + 𝑐𝑡)
1 1 𝑥+𝑐𝑡
𝑤 (𝑥, 𝑡) = [𝑤0 (𝑥 − 𝑐𝑡) + 𝑤0 (𝑥 + 𝑐𝑡) − ∫ 𝑤̇ (𝑢)𝑑𝑢 ]
2 2𝑐 𝑥−𝑐𝑡 0
Posons :
𝑤(𝑥, 𝑡) = 𝑤𝑑 (𝑥, 𝑡) + 𝑤𝑣 (𝑥, 𝑡)
Où 𝑤𝑑 (𝑥, 𝑡) représente une onde progressive due au déplacement initial
𝑤0 (𝑥 ) avec une vitesse initiale nulle et 𝑤𝑣 (𝑥, 𝑡) représente une onde progressive
due à la vitesse initiale 𝑤̇0 (𝑥 ) avec un déplacement initial nul.
b. Solution utilisant une Transformation de Fourier :
Pour trouver la réponse libre de la vibration d’une corde infinie
(−∞ < 𝑥 < +∞) sous les conditions initiales données, nous utilisons une
transformation de Fourier de l’équation (2.4). Pour cela multiplions cette
équation par 𝑒 𝑖𝑎𝑥 et intégrons de 𝑥 = −∞ à 𝑥 = +∞:
+∞ 2 (
𝜕 𝑤 𝑥, 𝑡) 𝑖𝑎𝑥 1 𝜕2 +∞ (2.5)
∫ 2
𝑒 𝑑𝑥 = 2 2
∫ 𝑤(𝑥, 𝑡)𝑒 𝑖𝑎𝑥 𝑑𝑥
−∞ 𝜕𝑥 𝑐 𝜕𝑡 −∞
+∞ 2 (
𝜕 𝑤 𝑥, 𝑡) 𝑖𝑎𝑥 𝜕𝑤 𝑖𝑎𝑥 +∞ +∞
𝜕𝑤 𝑖𝑎𝑥
∫ 2
𝑒 𝑑𝑥 = [ 𝑒 ] − 𝑖𝑎 ∫ 𝑒 𝑑𝑥
−∞ 𝜕𝑥 𝜕𝑥 −∞ −∞ 𝜕𝑥
+∞ +∞
𝜕𝑤 +∞
= [ 𝑒 𝑖𝑎𝑥 ] −= [𝑖𝑎𝑤𝑒 𝑖𝑎𝑥 ] − (𝑖𝑎)2 ∫ 𝑤𝑒 𝑖𝑎𝑥 𝑑𝑥
𝜕𝑥 −∞
−∞
−∞
+∞ +∞
𝜕𝑤
= [( − 𝑖𝑎𝑤) 𝑒 𝑖𝑎𝑥 ] − 𝑎2 ∫ 𝑤𝑒 𝑖𝑎𝑥 𝑑𝑥
𝜕𝑥 −∞ −∞
𝜕𝑤
En supposant que 𝑤 et tendent vers 0 quand |𝑥 | ⟶ ∞ le premier terme du
𝜕𝑥
second membre de l’équation, la transformation de Fourier de 𝑤(𝑥, 𝑡) est
définie par :
+∞
1
𝑊 (𝑎, 𝑡) = ∫ 𝑤(𝑥, 𝑡)𝑒 𝑖𝑎𝑥 𝑑𝑥
√ 2𝜋 −∞
Alors l’équation (2.5) devient :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 25
Dynamique des Structures
1 𝜕2 𝑊 (2.6)
+ 𝑎2 𝑊 = 0
𝑐 2 𝜕𝑡 2
Notons que l’utilisation de la transformation de Fourier réduit l’équation aux
dérivées partielles (2.4) en une équation différentielle ordinaire (2.6).
Les solutions de l’équation (6) peuvent s’exprimer par :
𝑊 (𝑎, 𝑡) = 𝐶1 𝑒 𝑖𝑎𝑐𝑡 + 𝐶2 𝑒 −𝑖𝑎𝑐𝑡 (2.7)
Où 𝐶1 et 𝐶2 sont des constantes que l’on peut déterminer en utilisant les
conditions initiales. En prenant les transformations de Fourier du
𝜕𝑤
déplacement [𝑤 = 𝑤0 (𝑥 )] initial et de la vitesse initiale [ = 𝑤̇0 (𝑥 )], on obtient : 𝜕𝑡
1 +∞
(2.8)
𝑊 (𝑎, 𝑡 = 0) = ∫ 𝑤0 (𝑥 )𝑒 𝑖𝑎𝑥 𝑑𝑥 = 𝑊0 (𝑎)
√2𝜋 −∞
𝑑𝑊 1 +∞
(2.9)
(𝑎, 𝑡 = 0) = ∫ 𝑤̇0 (𝑥 )𝑒 𝑖𝑎𝑥 𝑑𝑥 = 𝑊̇0 (𝑎)
𝑑𝑡 √ 2𝜋 −∞
En utilisant les expressions (2.8) et (2.9) dans l’expression (2.6) on a :
𝑊0 (𝑎) = 𝐶1 + 𝐶2 𝑒𝑡 𝑊̇0 (𝑎) = 𝑖𝑎𝑐 (𝐶1 − 𝐶2 ).
Ce qui donne :
𝑊 1 𝑊 1
𝐶1 = 0 + 𝑊̇0 et 𝐶2 = 0 − 𝑊̇
2 2𝑖𝑎𝑐 2 2𝑖𝑎𝑐 0
Finalement l’expression (2.7) nous donne :
𝑊 1 𝑊 1
𝑊 (𝑎, 𝑡) = ( 0 + 𝑊̇0 ) 𝑒 𝑖𝑎𝑐𝑡 + ( 0 − 𝑊̇ ) 𝑒 −𝑖𝑎𝑐𝑡
2 2𝑖𝑎𝑐 2 2𝑖𝑎𝑐 0
1 1 (2.10)
𝑊(𝑎, 𝑡) = 𝑊0 (𝑎)(𝑒 𝑖𝑎𝑐𝑡 + 𝑒 −𝑖𝑎𝑐𝑡 ) + 𝑊̇0 (𝑎)(𝑒 𝑖𝑎𝑐𝑡 − 𝑒 −𝑖𝑎𝑐𝑡 )
2 2𝑖𝑎𝑐
En utilisant la transformation inverse de Fourier, on obtient :
+∞
1
𝑤(𝑥, 𝑡) = ∫ 𝑊(𝑎, 𝑡) 𝑒 −𝑖𝑎𝑥 𝑑𝑎
√ 2𝜋 −∞
Ce qui va donner, en utilisant l’expression (2.10):
1 1 +∞
(𝑥−𝑐𝑡) (𝑥+𝑐𝑡) 1 1 +∞
𝑊̇0 (𝑎)
𝑤 (𝑥, 𝑡) = { ∫ 𝑊0 (𝑎 ) [𝑒 −𝑖𝑎 + 𝑒 −𝑖𝑎 ] 𝑑𝑎} + { ∫ [𝑒 −𝑖𝑎(𝑥−𝑐𝑡) − 𝑒 −𝑖𝑎 (𝑥+𝑐𝑡) ]𝑑𝑎}
2 √2𝜋 −∞ 2𝑐 √2𝜋 −∞ 𝑖𝑎
4. Vibration libre d’une corde de longueur finie :
L’équation (4) peut être résolue en utilisant la méthode de séparation des
variables. Dans cette méthode, la solution se pose sous la forme :
𝑤(𝑥, 𝑡) = 𝑊(𝑥 )𝑇 (𝑡)
Où 𝑊 (𝑥 ) est une fonction de 𝑥 seulement et 𝑇 (𝑡) une fonction de 𝑡
seulement. En substituant cette solution dans l’équation (2.4), on obtient :
𝑐 2 𝜕2 𝑊 1 𝜕2 𝑇
=
𝑊 2 𝜕𝑥 2 𝑇 𝜕𝑡 2
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 26
Dynamique des Structures
En remarquant que le premier membre ne dépend que de 𝑥 et le second
membre ne dépend que de 𝑡, alors les valeurs commune doit être égale à une
constante 𝑎. On a donc :
𝑐 2 𝜕2 𝑊 1 𝜕2 𝑇 (2.11)
2
= 2
=𝑎
𝑊 𝜕𝑥 𝑇 𝜕𝑡
L’équation (2.11) nous permet d’avoir deux équations séparées :
𝑐 2 𝜕2 𝑊 1 𝜕2 𝑇
= 𝑎 et =𝑎
𝑊 𝜕𝑥 2 𝑇 𝜕𝑡 2
Soit :
𝜕2 𝑊 𝑎 (2.12)
− 𝑊 =0
𝜕𝑥 2 𝑐 2
𝜕2 𝑇 (2.13)
− 𝑎𝑇 = 0
𝜕𝑡 2
Choisissons la constante 𝑎 négative et posons : 𝑎 = −𝜔2 . Alors les
équations (2.12) et (2.13) deviennent :
𝜕2 𝑊 𝜔 2 (2.14)
− 2 𝑊 =0
𝜕𝑥 2 𝑐
𝜕2 𝑇 (2.15)
2
− 𝜔2 𝑇 = 0
𝜕𝑡
Les solutions des équations (2.14) et (2.15) peuvent être exprimées par :
𝜔𝑥 𝜔𝑥 (2.16)
𝑊 (𝑥 ) = 𝐴 cos + 𝐵 sin
𝑐 𝑐
𝑇(𝑡) = 𝐶 cos 𝜔𝑡 + 𝐷 sin 𝜔𝑡 (2.17)
où 𝜔 est la fréquence de vibration, les constantes 𝐴 et 𝐵 peuvent être
déterminées à partir des conditions aux limites et les constantes 𝐶 et 𝐷 à
partir des conditions initiales de la corde.
Etude de cas : Vibration libre d’une avec les deux extrémités fixées
Lorsque les deux extrémités de la corde sont fixées, les conditions aux limites
sont données par :
𝑤 (0, 𝑡) = 𝑤 (𝑙, 𝑡) = 0 , 𝑡≥0 (2.18)
Les équations (2.16) et (2.18) nous donnent alors : 𝑊 (0) = 0 et 𝑊(𝑙 ) = 0. Ce qui
implique :
𝜔𝑙
𝐴 = 0 et 𝐵 sin = 0.
𝑐
Pour une solution non triviale, la constante 𝐵 ne peut être nulle, alors :
𝜔𝑙 (2.19)
sin =0
𝑐
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 27
Dynamique des Structures
L’équation (2.19) est appelée équation de fréquence ou équation
caractéristique.
Les solutions de cette équations sont appelées les fréquences naturelles de la
corde. L’équation (2.19) nous donne :
𝜔𝑙 𝑛𝑐𝜋 (2.20)
= 𝑛𝜋 soit 𝜔𝑛 = 𝑛 = 1, 2, …
𝑐 𝑙
Le déplacement transversale de la corde correspondant à 𝜔𝑛 , appelé n-ieme
mode de vibration ou n-ie’me mode normal est donné par :
𝑛𝜋𝑥 𝑛𝑐𝜋𝑡 𝑛𝑐𝜋𝑡 (2.21)
𝑤𝑛 (𝑥, 𝑡) = 𝑊𝑛 (𝑥 )𝑇𝑛 (𝑡) = sin (𝐶𝑛 cos + 𝐷𝑛 sin )
𝑙 𝑙 𝑙
Dans le n-ieme mode, chaque point de la corde vibre avec une amplitude
proportionnelle à la valeur d 𝑊𝑛 à ce point avec une fréquence circulaire 𝜔𝑛 .
Le mode correspondant à 𝑛 = 1 est appelé mode fondamental, 𝜔1 est appelé
fréquence fondamentale et on définit la période fondamentale 𝜏1 par :
2𝜋 2𝑙
𝜏1 = = .
𝜔1 𝑐
La vibration libre de la corde qui satisfait aux conditions aux limites peut être
trouvée par superposition de tous les modes naturels 𝑤𝑛 (𝑥, 𝑡) , on a :
∞ ∞
𝑛𝜋𝑥 𝑛𝑐𝜋𝑡 𝑛𝑐𝜋𝑡 (2.22)
𝑤 (𝑥, 𝑡) = ∑ 𝑤𝑛 (𝑥, 𝑡) = ∑ 𝑠𝑖𝑛 (𝐶𝑛 𝑐𝑜𝑠 + 𝐷𝑛 𝑠𝑖𝑛 )
𝑙 𝑙 𝑙
𝑛=1 𝑛=1
Cette solution représente la solution générale de l’équation (2.4) et inclus
toutes les vibrations possibles de la corde.
5. Vibration forcée d’une corde :
L’équation du mouvement d’une corde soumise à une vibration forcée par une
distribution de charge 𝑓(𝑥, 𝑡) par unité de longueur est donnée par :
𝜕2 𝑤(𝑥, 𝑡) 𝜕2 𝑤 (𝑥, 𝑡) (2.23)
(
𝜌 𝑥 ) − 𝑇 = 𝑓 ( 𝑥, 𝑡)
𝜕𝑡 2 𝜕𝑥 2
Soit une corde fixée aux deux extrémités telle que les conditions aux limites
sont :
𝑤(0, 𝑡) = 0 et 𝑤(𝑙, 𝑡) = 0.
La solution de l’équation homogène (sans second membre, 𝑓(𝑥, 𝑡) = 0) peut
être exprimée par (voir l’équation (2.22) ) :
∞
𝑛𝜋𝑥 𝑛𝑐𝜋𝑡 𝑛𝑐𝜋𝑡 (2.24)
𝑤 (𝑥, 𝑡) = ∑ 𝑠𝑖𝑛 (𝐶𝑛 𝑐𝑜𝑠 + 𝐷𝑛 𝑠𝑖𝑛 )
𝑙 𝑙 𝑙
𝑛=1
La solution de l’équation non homogène peut être supposée de la forme :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 28
Dynamique des Structures
∞
𝑛𝜋𝑥 (2.25)
𝑤(𝑥, 𝑡) = ∑ 𝑠𝑖𝑛 𝜂 (𝑡)
𝑙 𝑛
𝑛=1
Où 𝜂𝑛 (𝑡) représente les coordonnées généralisées. En remplaçant la solution
(2.25) dans l’équation (2.23), on obtient :
∞ ∞
𝑛𝜋𝑥 𝑑 2 𝜂𝑛 (𝑡) 𝑛𝜋 2 𝑛𝜋𝑥 (2.26)
𝜌 ∑ sin + 𝑇 ∑ ( ) sin 𝜂 ( 𝑡) = 𝑓 ( 𝑥, 𝑡)
𝑙 𝑑𝑡 2 𝑙 𝑙 𝑛
𝑛=1 𝑛=1
𝑛𝜋𝑥
Multiplions l’équation (2.26) par sin et intégrons de 0 à 𝑙, de même que
𝑙
l’utilisation de l’orthogonalité des fonctions sin(𝑖𝜋𝑥/𝑙), on obtient :
𝑑 2 𝜂𝑛 (𝑡) 𝑛2 𝑐 2 𝜋 2 2 (2.27)
2
+ 2
𝜂𝑛 (𝑡) = 𝑄𝑛 (𝑡)
𝑑𝑡 𝑙 𝜌𝑙
Où :
𝑙
𝑛𝜋𝑥
𝑄𝑛 (𝑡) = ∫ 𝑓(𝑥, 𝑡) sin 𝑑𝑥
0 𝑙
La solution de l’équation (2.27) incluant la solution homogène et la solution
particulière peut être exprimé par :
𝑡
𝑛𝑐𝜋𝑡 𝑛𝑐𝜋𝑡 2 𝑛𝑐𝜋 (𝑡 − 𝜏)
𝜂𝑛 (𝑡) = 𝐶𝑛 cos + 𝐷𝑛 sin + ∫ 𝑄𝑛 (𝑡) sin 𝑑𝑡
𝑙 𝑙 𝑛𝑐𝜋𝜌 0 𝑙
Alors au vu de l’équation (2.26), la réponse forcée de la vibration de la corde
est donnée par :
∞ ∞
𝑛𝑐𝜋𝑡 𝑛𝑐𝜋𝑡 𝑛𝜋𝑥 2 1 𝑛𝜋𝑥 𝑡 𝑛𝑐𝜋(𝑡 − 𝜏)
𝑤(𝑥, 𝑡) = ∑ (𝐶𝑛 𝑐𝑜𝑠 + 𝐷𝑛 𝑠𝑖𝑛 ) 𝑠𝑖𝑛 + ∑ sin ∫ 𝑄𝑛 (𝑡) 𝑠𝑖𝑛 𝑑𝑡
𝑙 𝑙 𝑙 𝑐𝜋𝜌 𝑛 𝑙 0 𝑙
𝑛=1 𝑛=1
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 29
Dynamique des Structures
II. Vibration transversale des poutres
2.1 Equation du mouvement :
Figure 10 : Vibration transversale d’une poutre mince
Soit une poutre de faible épaisseur sujette à une force transversale (voir Fig.
3(a)). Considérons le diagramme du corps libre d’un élément de poutre de
longueur 𝑑𝑥 (voir Fig. 3(b)) ou 𝑀(𝑥, 𝑡) est le moment fléchissant, 𝑉(𝑥, 𝑡) l’effort
tranchant et 𝑓(𝑥, 𝑡) la force transversale extérieure par unité de longueur de la
poutre.
𝜕2 𝑤
La force d’inertie est égale à : 𝜌 ∙ 𝐴 (𝑥 ) ∙ 𝑑𝑥 ∙ (𝑥, 𝑡) .
𝜕𝑡2
L’application de la deuxième loi de Newton dans la direction 𝑧 donne :
𝜕2 𝑤 (2.28)
−(𝑉 + 𝑑𝑉 ) + 𝑓(𝑥, 𝑡) 𝑑𝑥 + 𝑉 = 𝜌 ∙ 𝐴(𝑥 ) ∙ 𝑑𝑥 ∙ 2 (𝑥, 𝑡)
𝜕𝑡
Où 𝜌 est la densité de masse et 𝐴 𝑥 la section de la poutre. L’équation
( )
d’équilibre des moments au point P donne :
𝑑𝑥 (2.29)
(𝑀 + 𝑑𝑀) − (𝑉 + 𝑑𝑉 )𝑑𝑥 + 𝑓 (𝑥, 𝑡) 𝑑𝑥 −𝑀=0
2
𝜕𝑀 𝜕𝑉
Posons : 𝑑𝑀 = 𝑑𝑥 et 𝑑𝑉 = 𝑑𝑥
𝜕𝑥 𝜕𝑥
Les équations (28) et (29) nous donnent alors :
𝜕𝑉 𝜕2 𝑤 (2.30)
− (𝑥, 𝑡) + 𝑓(𝑥, 𝑡) = 𝜌𝐴(𝑥 ) 2 (𝑥, 𝑡)
𝜕𝑥 𝜕𝑡
𝜕𝑀 (2.31)
(𝑥, 𝑡) − 𝑉 (𝑥, 𝑡) = 0
𝜕𝑥
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 30
Dynamique des Structures
En utilisant la relation : 𝑉 = 𝜕𝑀⁄𝜕𝑥 les équations (2.30) et (2.31) donnent alors :
𝜕2 𝑀 𝜕2 𝑤 (2.32)
− 2
( 𝑥, 𝑡) + 𝑓 ( 𝑥, 𝑡) = 𝜌𝐴 ( 𝑥 ) 2
(𝑥, 𝑡)
𝜕𝑥 𝜕𝑡
De la théorie élémentaire de la flexion des poutres, la relation entre le moment
fléchissant et la flèche (déplacement) peut être exprimée par :
𝜕2 𝑤
( )
𝑀 𝑥, 𝑡 = 𝐸𝐼 𝑥 ( ) (𝑥, 𝑡)
𝜕𝑥 2
où 𝐸 est le module de Young et 𝐼(𝑥 ) le moment d’inertie de la section droite
de la poutre par rapport à l’axe 𝑦.
L’équation (2.32) devient alors :
𝜕2 𝜕2 𝑤 𝜕2 𝑤 (2.33)
[𝐸𝐼( 𝑥 ) ( 𝑥, 𝑡) ] + 𝜌𝐴 ( 𝑥 ) ( 𝑥, 𝑡 ) = 𝑓 ( 𝑥, 𝑡)
𝜕𝑥 2 𝜕𝑥 2 𝜕𝑡 2
Pour une poutre uniforme, l’équation (33) devient :
𝜕4 𝑤 𝜕2 𝑤 (2.34)
𝐸𝐼 4 (𝑥, 𝑡) + 𝜌𝐴 2 (𝑥, 𝑡) = 𝑓 (𝑥, 𝑡)
𝜕𝑥 𝜕𝑡
L’équation (2.34) peut être résolue en utilisant les conditions initiales et aux
limites appropriées. Si la poutre est donnée avec un déplacement in itial
𝑤0 (𝑥 ) et une vitesse initiale 𝑤̇0 (𝑥 ), les conditions initiales peuvent être
exprimées par :
𝑤(𝑥, 𝑡 = 0) = 𝑤0 (𝑥 ), 0 ≤ 𝑥 ≤ 𝑙
𝜕
𝑤(𝑥, 𝑡 = 0) = 𝑤̇0 (𝑥 ), 0 ≤ 𝑥 ≤ 𝑙
𝜕𝑡
Si la poutre est fixe à 𝑥 = 0 et fixée à 𝑥 = 𝑙, la deflection et la pente sont nuls
à 𝑥 = 0 et la deflection et le moment fléchissant sont nuls à 𝑥 = 𝑙. Alors les
conditions aux limites sont :
𝑤(𝑥 = 0, 𝑡) = 0, 𝑡 ≥ 0
𝜕𝑤
(𝑥 = 0, 𝑡) = 0, 𝑡≥ 0
𝜕𝑥
𝑤 (𝑥 = 𝑙, 𝑡) = 0, 𝑡 ≥ 0
2
𝜕 𝑤
(𝑥 = 𝑙, 𝑡) = 0, 𝑡≥0
𝜕𝑥 2
2.2 Vibrations libres :
Dans ce cas la force d’excitation est nulle, posons 𝑐 = √𝐸𝐼/𝜌𝐴 l’équation
(2.34) devient alors :
𝜕4 𝑤 𝜕2 𝑤 (2.35)
𝑐2 ( 𝑥, 𝑡) + (𝑥, 𝑡) = 0
𝜕𝑥 4 𝜕𝑡 2
Cherchons des solutions pour l’équation (2.35) en utilisant la méthode de
séparation des variables.
Soit 𝑤(𝑥; 𝑡) = 𝑊 (𝑥 )𝑇 (𝑡)
En utilisant cette solution dans l’équation (2.35), on a :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 31
Dynamique des Structures
𝑐 2 𝜕4 𝑊 ( 𝑥 ) 1 𝜕2 𝑇 (𝑡)
= − = 𝑎 = 𝜔2
𝑊 (𝑥 ) 𝜕𝑥 4 𝑇 (𝑡) 𝜕𝑡 2
On obtient alors les deux équations suivantes :
𝑐 2 𝜕4 𝑊 ( 𝑥 )
= 𝜔2
𝑊(𝑥 ) 𝜕𝑥 4
1 𝜕2 𝑇 (𝑡)
= 𝜔2
𝑇(𝑡) 𝜕𝑡 2
Posons 𝛽4 = 𝜔2 /𝑐 2
On a finalement :
𝜕4 𝑊 ( 𝑥 ) (2.36)
− 𝛽4 𝑊 (𝑥 ) = 0
𝜕𝑥 4
𝜕2 𝑇 (𝑡) (2.37)
− 𝜔2 𝑇 (𝑡) = 0
𝜕𝑡 2
La solution de l’équation (2.37) est donnée sous la forme :
𝑇 (𝑡) = 𝐴 cos 𝜔𝑡 + 𝐵 sin 𝜔𝑡
où 𝐴 et 𝐵 sont des constantes déterminées à l’aide des conditions initiales.
Supposons que la solution de l’équation (2.36) soit de la forme exponentielle
suivante :
𝑊 (𝑥 ) = 𝐶𝑒 𝑠𝑥
où 𝐶 𝑒𝑡 𝑠 sont des constantes. Remplaçons cette solution dans l’équation
(36), on obtient l’équation auxiliaire suivante dont l’inconnue est 𝑠:
𝑠 4 − 𝛽4 = 0
Cette équation admet quatre solutions
𝑠1,2 = ±𝛽 et 𝑠3,4 = ±𝑖𝛽
Alors la solution de l’équation (36) est exprimée par :
𝑊 (𝑥 ) = 𝐶1 𝑒 𝛽𝑥 + 𝐶2 𝑒 −𝛽𝑥 + 𝐶3 𝑒 𝑖𝛽𝑥 + 𝐶4 𝑒 −𝑖𝛽𝑥
où 𝐶1 , 𝐶2 , 𝐶3 , et 𝐶4 sont des constantes.
La solution peut être exprimée sous forme :
𝑊 (𝑥 ) = 𝐶1 cos 𝛽𝑥 + 𝐶2 sin 𝛽𝑥 + 𝐶3 cosh 𝛽𝑥 + 𝐶4 sinh 𝛽𝑥
où
𝑊(𝑥) = 𝐶1 (cos𝛽𝑥 + cosh 𝛽𝑥) + 𝐶2(cos𝛽𝑥 − cosh 𝛽𝑥) + 𝐶3 (sin 𝛽𝑥 + sinh 𝛽𝑥 ) (2.38)
+ 𝐶4 (sin 𝛽𝑥 − sinh 𝛽𝑥)
où 𝐶1 , 𝐶2 , 𝐶3 , et 𝐶4 sont des constantes différentes à chaque fois.
Les fréquences naturelles de la poutre peuvent être exprimées à l’aide
l’équation :
4
𝜔2 𝐸𝐼 𝐸𝐼
𝛽 = 2 𝑠𝑜𝑖𝑡 𝜔 = 𝛽2 √ = (𝛽𝑙 )2 √
𝑐 𝜌𝐴 𝜌𝐴𝑙 4
La fonction 𝑊(𝑥 ) appelée mode normal ou fonction caractéristique de la
poutre et 𝜔 est appelé fréquence naturelle de la vibration.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 32
Dynamique des Structures
Pour n’importe quelle poutre, il y a un nombre infini de mode normal avec
une fréquence naturelle associée pour chaque mode normal. Les constantes
𝐶1 , 𝐶2 , 𝐶3 , et 𝐶4 et la valeur de 𝛽 peuvent être déterminées avec les conditions
aux limites de la poutre.
Si on note 𝜔𝑖 la i-eme fréquence naturelle et 𝑊𝑖 (𝑥 ) le mode normal
correspondant, la réponse totale de la vibration libre de la poutre peut être
trouvée par superposition des modes normaux :
∞
𝑤(𝑥, 𝑡) = ∑ 𝑊𝑖 (𝑥 )(𝐴𝑖 cos 𝜔𝑖 𝑡 + 𝐵𝑖 sin 𝜔𝑖 𝑡)
𝑖 =1
où les constantes 𝐴𝑖 et 𝐵𝑖 peuvent être déterminées avec les conditions
initiales de la poutre.
2.3 Fréquences et formes des modes des poutres uniformes :
a. Poutre sur Appuis Simples aux deux extrémités :
Le déplacement transversal (flèche) et le moment fléchissant sont nuls à une
extrémité sur appui simple. Ainsi les conditions aux limites peuvent être
décrites par :
𝑑2 𝑊
𝐴 𝑥 = 0, on a: 𝑊(0) = 0 et 𝐸𝐼 ( 0) = 0
𝑑𝑥 2
2
𝑑 𝑊
𝑊 ( 0) = 0 et ( 0) = 0
𝑑𝑥 2
𝑑2 𝑊
𝐴 𝑥 = 𝑙, on a: 𝑊(𝑙 ) = 0 et 𝐸𝐼 2 (𝑙 ) = 0
𝑑𝑥
2
𝑑 𝑊
𝑊 (𝑙 ) = 0 et (𝑙 ) = 0
𝑑𝑥 2
En utilisant les conditions à 𝑥 = 0 dans la solution (équation 2.38), on
obtient :
𝐶1 = 𝐶2 = 0
En utilisant les conditions à 𝑥 = 𝑙 dans la solution (équation 2.38), on obtient :
𝐶3 (sin 𝛽𝑙 + sinh 𝛽𝑙 ) + 𝐶4 (sin 𝛽𝑙 − sinh 𝛽𝑙 ) = 0 (2.39)
( ) (
−𝐶3 sin 𝛽𝑙 − sinh 𝛽𝑙 − 𝐶4 sin 𝛽𝑙 + sinh 𝛽𝑙 = 0 ) (2.40)
Nous pouvons donc former le système d’équation suivant :
𝐶 (sin 𝛽𝑙 + sinh 𝛽𝑙) + 𝐶4 (sin 𝛽𝑙 − sinh 𝛽𝑙) = 0
{ 3
−𝐶3 (sin 𝛽𝑙 − sinh 𝛽𝑙) − 𝐶4 (sin 𝛽𝑙 + sinh 𝛽𝑙) = 0
Pour avoir une solution non triviale, le déterminant du système doit être nul :
−(sin 𝛽𝑙 + sinh 𝛽𝑙 )2 + (sin 𝛽𝑙 − sinh 𝛽𝑙 )2 = 0
Soit : sin 𝛽𝑙 sinh 𝛽𝑙 = 0
Or sinh 𝛽𝑙 n’est nul que pour 𝛽 = 0, soit 𝜔 = 0.
Donc il en résulte l’équation de fréquence
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 33
Dynamique des Structures
sin 𝛽𝑙 = 0 (2.41)
Les racines 𝛽𝑛 𝑙 de cette équation sont données par :
𝛽𝑛 𝑙 = 𝑛𝜋, 𝑖 = 1,2, …
Les fréquences naturelles de la vibration deviennent alors :
𝐸𝐼 1/2 𝐸𝐼 1/2
𝜔𝑛 = (𝛽𝑛 𝑙 )2 ( ) = 𝑛2 2
𝜋 ( ) , 𝑛 = 1,2, ….
𝜌𝐴𝑙 4 𝜌𝐴𝑙 4
En substituant l’équation (41) dans l’équation (39) on obtient : 𝐶3 = 𝐶4 .
Alors la forme des modes est donnée par :
𝑛𝜋𝑥
𝑊𝑛 (𝑥 ) = 𝐶𝑛 sin 𝛽𝑛 𝑥 = 𝐶𝑛 sin , 𝑛 = 1,2, …
𝑙
Les modes normaux de la vibration sont données par :
𝑤𝑛 (𝑥, 𝑡) = 𝑊𝑛 (𝑥 )(𝐴𝑛 cos 𝜔𝑛 𝑡 + 𝐵𝑛 sin 𝜔𝑛 𝑡 ) 𝑛 = 1,2, …
La solution totale de la vibration libre peut être exprimée par :
∞ ∞
𝑛𝜋𝑥 (2.42)
𝑤(𝑥, 𝑡) = ∑ 𝑤𝑛 (𝑥, 𝑡) = ∑ sin (𝐴𝑛 cos 𝜔𝑛 𝑡 + 𝐵𝑛 sin 𝜔𝑛 𝑡)
𝑙
𝑛=1 𝑛=1
Si les conditions initiales sont données par :
𝜕𝑤
𝑤(𝑥, 0) = 𝑤0 (𝑥 ) et (𝑥, 0) = 𝑤̇0 (𝑥 )
𝜕𝑡
Alors on obtient en utilisant l’équation (42) :
∞
𝑛𝜋𝑥
∑ 𝐴𝑛 sin = 𝑤0 (𝑥 )
𝑙
𝑛=1
∞
𝑛𝜋𝑥
∑ 𝜔𝑛 𝐵𝑛 sin = 𝑤̇0 (𝑥 )
𝑙
𝑛=1
𝑛𝜋𝑥
Multiplions ces deux expressions par sin et intégrons de 0 à 𝑙, on obtient ;
𝑙
𝑙
2 𝑛𝜋𝑥
𝐴𝑛 = ∫ 𝑤0 (𝑥 ) sin 𝑑𝑥
𝑙 0 𝑙
𝑙
2 𝑛𝜋𝑥
𝐵𝑛 = ∫ 𝑤̇0 (𝑥 ) sin 𝑑𝑥
𝑙𝜔𝑛 0 𝑙
2.4 Orthogonalité des modes normaux :
Le problème de valeur propre correspondant à une poutre non uniforme peut
être obtenu en supposant une solution harmonique de fréquence 𝜔 de
l’équation (2.33)
𝑑2 𝑑 2 𝑊 (𝑥 ) (2.43)
2
[𝐸𝐼 ( 𝑥 ) 2
] = 𝜔2 𝜌𝐴 (𝑥 )𝑊 (𝑥 )
𝑑𝑥 𝑑𝑥
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 34
Dynamique des Structures
Pour déterminer les relations d’orthogonalité pour les poutres, considérons
deux valeurs propres 𝜔2𝑖 et 𝜔𝑗2 et les fonctions normales ou propres
correspondantes 𝑊𝑖 (𝑥 ) et 𝑊𝑗 (𝑥) respectivement, on a donc :
𝑑2 𝑑 2 𝑊𝑖 (𝑥 ) (2.44)
2
[𝐸𝐼( 𝑥 ) 2
] = 𝜔2 𝜌𝐴 (𝑥 )𝑊𝑖 (𝑥 )
𝑑𝑥 𝑑𝑥
Et
𝑑2 𝑑 2 𝑊𝑗 (𝑥 ) (2.45)
2
(
[𝐸𝐼 𝑥 ) 2
] = 𝜔𝑗2 𝜌𝐴 (𝑥 )𝑊𝑗 (𝑥)
𝑑𝑥 𝑑𝑥
Multiplions l’équation (2.44) par 𝑊𝑗 (𝑥) et intégrons sur la longueur de la
poutre, on obtient :
𝑙
𝑑2 𝑑 2 𝑊𝑖 (𝑥 ) 2
𝑙
(2.46)
∫ 𝑊𝑗 (𝑥 ) 2 [𝐸𝐼 (𝑥 ) 2
] 𝑑𝑥 = 𝜔𝑖 ∫ 𝜌𝐴(𝑥 )𝑊𝑗 (𝑥 )𝑊𝑖 (𝑥)𝑑𝑥
0 𝑑𝑥 𝑑𝑥 0
Par intégrant le premier membre de l’équation (2.46) deux fois par parties et
en utilisant une quelconque combinaison des conditions au x limites de la
poutre nous obtenons :
𝑙
𝑑2 𝑑 2 𝑊𝑖 (𝑥 )
∫ 𝑊𝑗 (𝑥 ) 2 [𝐸𝐼 (𝑥 ) ] 𝑑𝑥
0 𝑑𝑥 𝑑𝑥2
𝑙 𝑙
𝑑 𝑑 2 𝑊𝑖 (𝑥 ) 𝑑 𝑊𝑗 (𝑥 ) 𝑑 2 𝑊𝑖 (𝑥 )
= |𝑊𝑗 (𝑥 ) [𝐸𝐼 (𝑥 ) ]| − | 𝐸𝐼(𝑥 ) |
𝑑𝑥 𝑑𝑥 2 𝑑𝑥 𝑑 𝑥 2
0 0
𝑙
𝑑 2 𝑊𝑗 (𝑥 ) 𝑑 2 𝑊𝑖 (𝑥 )
+ ∫ 𝐸𝐼 (𝑥 ) 𝑑𝑥
0 𝑑𝑥 2 𝑑𝑥 2
𝑙
𝑑 2 𝑊𝑗 (𝑥 ) 𝑑 2 𝑊𝑖 (𝑥 ) (2.47)
= ∫ 𝐸𝐼(𝑥 ) 𝑑𝑥
0 𝑑𝑥2 𝑑𝑥2
Alors l’équation (2.46) peut être récrite :
𝑙
𝑑 2 𝑊𝑗 (𝑥 ) 𝑑 2 𝑊𝑖 (𝑥 ) 2
𝑙 (2.48)
∫ 𝐸𝐼 (𝑥 ) 𝑑𝑥 = 𝜔 𝑖 ∫ 𝜌𝐴(𝑥 )𝑊𝑗 (𝑥 )𝑊𝑖 (𝑥)𝑑𝑥
0 𝑑𝑥 2 𝑑𝑥 2 0
De façon similaire en multipliant par 𝑊𝑖 (𝑥) l’équation (2.45) et en intégrant
sur la longueur de la poutre, nous pouvons obtenir :
𝑙
𝑑 2 𝑊𝑖 (𝑥 ) 𝑑 2 𝑊𝑗 (𝑥 ) 2
𝑙 (2.49)
∫ 𝐸𝐼 (𝑥 ) 𝑑𝑥 = 𝜔𝑗 ∫ 𝜌𝐴( 𝑥 ) 𝑊𝑖 ( 𝑥 ) 𝑊𝑗 ( 𝑥 ) 𝑑𝑥
0 𝑑𝑥 2 𝑑𝑥 2 0
A noter que l’ordre des indices 𝑖 et 𝑗 sous les intégrales n’est pas déterminant.
En effectuant la soustraction membre à membre des équations (2.48) et (2.49),
on obtient :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 35
Dynamique des Structures
𝑙 (2.50)
(𝜔2𝑖 − 2
𝜔𝑗 ) ∫ 𝜌𝐴(𝑥 )𝑊𝑖 (𝑥 )𝑊𝑗 (𝑥)𝑑𝑥 =0
0
Dès lors que les valeurs propres sont différentes, on obtient :
𝑙 (2.51)
∫ 𝜌𝐴 (𝑥 )𝑊𝑖 (𝑥)𝑊𝑗 (𝑥)𝑑𝑥 = 0 𝑖, 𝑗 = 1,2, … . 𝜔2𝑖 ≠ 𝜔𝑗2
0
Il vient (voir équations (51), (49), (48)) :
𝑙
𝑑 2 𝑊𝑖 (𝑥 ) 𝑑 2 𝑊𝑗 (𝑥 ) (2.52)
∫ 𝐸𝐼 (𝑥 ) 𝑑𝑥 = 0 𝑖, 𝑗 = 1,2, …. 𝜔2𝑖 ≠ 𝜔𝑗2
0 𝑑𝑥 2 𝑑𝑥 2
L’équation (2.51) est appelée relation d’orthogonalité pour les fonctions
normales, l’équation (2.52) est une autre forme de la condition
d’orthogonalité pour les modes normaux.
En normalisant le i-ème mode normal tel que :
𝑙 (2.53)
∫ 𝜌𝐴(𝑥 )𝑊𝑖 2 (𝑥 )𝑑𝑥 = 1 𝑖 = 1,2, ….
0
Les équations (51) et (53) peuvent être exprimées dans la forme suivante :
𝑙 (54)
∫ 𝜌𝐴(𝑥 )𝑊𝑖 (𝑥 )𝑊𝑗 (𝑥 )𝑑𝑥 = 𝛿𝑖𝑗
0
Où 𝛿𝑖𝑗 est le symbole de Kronecker :
1 𝑠𝑖 𝑖 = 𝑗
𝛿𝑖𝑗 = {
0 𝑠𝑖 𝑖 ≠ 𝑗
En utilisant les équations (2.54) et (2.46), on obtient alors une autre forme de
la relation d’orthogonalité
𝑙
𝑑2 𝑑 2 𝑊𝑖 (𝑥 ) (2.55)
∫ 𝑊𝑗 (𝑥 ) 2 [𝐸𝐼 (𝑥 ) 2
] 𝑑𝑥 = 𝜔2𝑖 𝛿𝑖𝑗
0 𝑑𝑥 𝑑𝑥
Suivant le théorème de l’expansion, toute fonction 𝑊(𝑥 ) qui satisfait les
conditions aux limites de la poutre dénote un possible déplacement
transversal de la poutre et peut être exprimées comme étant une somme des
fonctions propres :
∞
(2.56)
𝑊(𝑥 ) = ∑ 𝑐𝑖 𝑊𝑖 (𝑥 )
𝑖=1
Où la constante 𝑐𝑖 est définie par :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 36
Dynamique des Structures
𝑙 (2.57)
𝑐𝑖 = ∫ 𝜌𝐴 (𝑥 )𝑊𝑖 (𝑥)𝑊𝑗 (𝑥 )𝑑𝑥 𝑖 = 1,2, …
0
Et
𝑙
𝑑2 𝑑 2 𝑊 (𝑥 ) (2.58)
𝜔2𝑖 𝑐𝑖 = ∫ 𝑊𝑖 (𝑥 ) [𝐸𝐼( 𝑥 ) ] 𝑑𝑥 𝑖 = 1,2, ….
0 𝑑𝑥 2 𝑑𝑥 2
2.5 Vibrations forcées :
L’équation du mouvement d’une poutre soumise à l’action d’une
distribution transversale de force 𝑓(𝑥, 𝑡) est la suivante (voir équation (34)) :
𝜕2 𝜕2 𝑤 𝜕2 𝑤 (2.59)
2
[𝐸𝐼( 𝑥 ) 2
( 𝑥, 𝑡) ] + 𝜌𝐴 ( 𝑥 ) 2
(𝑥, 𝑡) = 𝑓(𝑥, 𝑡)
𝜕𝑥 𝜕𝑥 𝜕𝑡
En utilisant l’approche du mode normal (analyse modale), la solution de
l’équation (2.59) est supposée être une combinaison linéaire des modes
normaux de la poutre :
∞
(2.60)
( ) ( ) (
𝑤 𝑥, 𝑡 = ∑ 𝑊𝑖 𝑥 𝜂𝑖 𝑡)
𝑖 =1
où 𝑊𝑖 (𝑥) sont les modes normaux trouvées en résolvant l’équation (en
utilisant les quatre conditions aux limites de la poutre)
𝜕2 𝜕2 𝑊𝑖 (2.61)
[𝐸𝐼 ( 𝑥 ) (𝑥, 𝑡) ] + 𝜌𝐴(𝑥 )𝜔2𝑖 𝑊𝑖 (𝑥 ) = 0
𝜕𝑥 2 𝜕𝑥 2
Et 𝜂𝑖 (𝑡) sont les coordonnées généralisées.
En utilisant l’équation (2.60) dans l’équation (2.59), on obtient :
∞ ∞
𝜕2 𝜕2 𝑤 𝑑 2 𝜂𝑖 (𝑡) (2.62)
∑ 2 [𝐸𝐼 (𝑥 ) 2 (𝑥, 𝑡) ] 𝜂𝑖 (𝑡) + 𝜌𝐴(𝑥 ) ∑ 𝑊𝑖 (𝑥 ) = 𝑓 (𝑥, 𝑡)
𝜕𝑥 𝜕𝑥 𝑑𝑡 2
𝑖=1 𝑖 =1
En utilisant l’équation (2.45) l’équation (2.46) nous donne :
∞ ∞
𝑑 2 𝜂𝑖 (𝑡) (2.63)
𝜌𝐴(𝑥 ) ∑ 𝜔2𝑖 𝑊𝑖 (𝑥 )𝜂𝑖 (𝑡) + 𝜌𝐴(𝑥 ) ∑ 𝑊𝑖 (𝑥 ) = 𝑓 (𝑥, 𝑡)
𝑑𝑡 2
𝑖=1 𝑖 =1
Multiplions l’équation (2.63) par 𝑊𝑗 (𝑥) et integrons entre 0 et 𝑙, on a :
∞ 𝑙 ∞
𝑑 2 𝜂𝑖(𝑡) 𝑙 𝑙
∑𝜂𝑖 (𝑡) ∫ 𝜌𝐴 (𝑥)𝜔2𝑖 𝑊𝑗 (𝑥)𝑊𝑖 (𝑥)𝑑𝑥 + ∑ ∫ 𝜌𝐴 (𝑥 )𝑊𝑗 (𝑥 )𝑊𝑖 (𝑥 )𝑑𝑥 = ∫ 𝑊𝑗 (𝑥)𝑓 (𝑥, 𝑡)𝑑𝑥
0 𝑑𝑡 2 0 0
𝑖=1 𝑖=1
Au vu de la condition d’orthogonalité, tous les termes dans chacune des
sommations du premier membre disparaissent à l’exception du terme pour
lequel 𝑖 = 𝑗, en laissant :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 37
Dynamique des Structures
𝑑 2 𝜂𝑖 (𝑡) (2.64)
+ 𝜔2𝑖 𝜂𝑖 (𝑡) = 𝑄𝑖 (𝑡), 𝑖 = 1,2, …
𝑑𝑡 2
où 𝑄𝑖 (𝑡) est la force généralisée correspondant au i-ème mode, donnée par :
𝑙
𝑄𝑖 (𝑡) = ∫ 𝑊𝑖 (𝑥 )𝑓(𝑥, 𝑡)𝑑𝑥
0
La solution complète de l’équation (2.64) peut être exprimée par :
1 𝑡 (2.65)
𝜂𝑖 (𝑡) = 𝐴𝑖 𝑐𝑜𝑠 𝜔𝑖 𝑡 + 𝐵𝑖 𝑠𝑖𝑛 𝜔𝑖 𝑡 + ∫ 𝑄𝑖 (𝑡) sin 𝜔𝑖 (𝑡 − 𝜏) 𝑑𝑡
𝜔𝑖 0
Ainsi, la solution de l’équation (2.59) est donnée par l’équation (2.60) et par
l’équation (2.65) :
∞
1 𝑡 (2.66)
𝑤 (𝑥, 𝑡) = ∑ [𝐴𝑖 𝑐𝑜𝑠 𝜔𝑖 𝑡 + 𝐵𝑖 𝑠𝑖𝑛 𝜔𝑖 𝑡 + ∫ 𝑄𝑖 (𝜏) sin 𝜔𝑖 (𝑡 − 𝜏)𝑑𝜏 ] 𝑊𝑖 (𝑥 )
𝜔𝑖 0
𝑖 =1
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 38
Dynamique des Structures
III. Vibrations des membranes :
Une membrane est une fine plaque ou lame parfaitement flexible qui est
soumise à une tension.
Il a une résistance négligeable aux forces de cisaillement ou de flexion, et les
forces de rappel surviennent exclusivement des forces d'étirement ou de
traction dans le plan. La peau et les diaphragmes des microphones à
condensateur sont des exemples de membranes.
3.1 Equation du mouvement :
Considérons une membrane homogène et parfaitement flexible délimitée par
une courbe plane 𝐶 dans le plan 𝑥𝑦 à l'état non déformé. Il est soumis à une
charge de pression d'intensité
𝑓 (𝑥, 𝑦, 𝑡) par unité de surface dans la direction transversale ou 𝑧 et une
tension d'amplitude 𝑇 par unité de longueur le long du bord comme dans le
cas d'une peau de tambour.
Chaque point de la membrane est supposé se déplacer uniquement dans la
direction 𝑧, et le déplacement, 𝑤(𝑥, 𝑦, 𝑡), est supposé être très petit par rapport
aux dimensions de la membrane. Considérons un élément de surface de la
membrane, 𝑑𝑥𝑑𝑦, avec des forces de traction d'amplitude 𝑇𝑑𝑥 et 𝑇𝑑𝑦 agissant
sur les côtés parallèles aux axes 𝑥 et 𝑦, respectivement, comme illustré à la
Figure 4.
Figure 11 : (a) membrane non déformée dans le plan xy (b) membrane déformée vue
dans le plan xz.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 39
Dynamique des Structures
Après déformation, les forces nettes agissant sur l'élément de la membrane le
long de la direction 𝑧 en raison des forces 𝑇𝑑𝑥 et 𝑇𝑑𝑦 sera :
𝜕2 𝑤 𝜕2 𝑤
𝑇 𝜕𝑦2 𝑑𝑥𝑑𝑦 et 𝑇 𝜕𝑥 2 𝑑𝑥𝑑𝑦
La force de pression qui agit sur l’élément de surface de la membrane dans
la direction 𝑧 est 𝑓 (𝑥, 𝑦)𝑑𝑥𝑑𝑦. La force d’inertie sur l’élément de surface est
donnée par :
𝜕2 𝑤
(
𝜌 𝑥, 𝑦 ) 𝑑𝑥𝑑𝑦
𝜕𝑡 2
Où 𝜌(𝑥, 𝑦) est la masse par unité de surface. L’application de la seconde loi
de Newton donne l’équation du mouvement de la vibration transversale
forcée de la membrane comme suit :
𝜕2 𝑤 𝜕2 𝑤 𝜕2 𝑤 (3.1)
𝑇( 2 + ) + 𝑓 = 𝜌
𝜕𝑥 𝜕𝑦 2 𝜕𝑡 2
Lorsque la force extérieure 𝑓(𝑥, 𝑦) = 0, l’équation du mouvement de la
vibration libre de la membrane est donnée par :
2
𝜕 2 𝑤 𝜕2 𝑤 𝜕2 𝑤 (3.2)
𝑐 ( 2+ ) =
𝜕𝑥 𝜕𝑦 2 𝜕𝑡 2
Avec 𝑐 = √𝑇⁄𝜌
Condition initiales et aux limites :
Les équations (3.1) 𝑒𝑡 (3.2) sont des équations aux dérivées partielles du
second ordre de chacune des variables 𝑥, 𝑦 𝑒𝑡 𝑡. Nous avons alors besoin de
deux conditions initiales et quatre conditions aux limites pour trouver une
solution unique au problème.
Habituellement, le déplacement et la vitesse de la membrane à 𝑡 = 0 sont
spécifiés par 𝑤0 (𝑥 ) et 𝑤̇0 (𝑥 ) respectivement. Alors les conditions initiales
sont données par
𝑤(𝑥, 𝑦, 0) = 𝑤0 (𝑥, 𝑦)
𝜕𝑤
(𝑥, 𝑦, 0) = 𝑤̇0 (𝑥, 𝑦)
𝜕𝑡
Les conditions aux limites de la membrane peuvent être déterminées comme
suit :
i. Si la membrane est fixée à n’importe quel point (𝑥1 , 𝑦1 ) de bordure
alors la flèche doit être nulle :
𝑤 (𝑥1 , 𝑦1 , 𝑡) = 0, 𝑡≥ 0
ii. Si la membrane est libre en déplacement transversal (dans la
direction z) en n’importe que point (𝑥 2 , 𝑦2 ) de la bordure, alors il ne
peut y avoir aucune force en ce point dans la direction z :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 40
Dynamique des Structures
𝜕𝑤
𝑇 (𝑥 , 𝑦 , 𝑡) = 0, 𝑡 ≥ 0
𝜕𝑛 2 2
𝜕𝑤
Où indique la dérivée de 𝑤 par rapport la direction normal n à la bordure
𝜕𝑛
de la membrane en ce point (𝑥 2 , 𝑦2 ).
3.2 Solution Ondulatoire :
Les fonctions :
𝑤1 (𝑥; 𝑦; 𝑡) = 𝑓(𝑥 − 𝑐𝑡) (3.3)
𝑤2 (𝑥; 𝑦; 𝑡) = 𝑓(𝑥 + 𝑐𝑡) (3.4)
𝑤3 (𝑥; 𝑦; 𝑡) = 𝑓(𝑥 cos 𝜃 + 𝑦 sin 𝜃 − 𝑐𝑡 ) (3.5)
Peuvent être vérifiées si elles sont solutions de l’équation (3.2). Par exemple
considérons la fonction 𝑤3 , les dérivées partielles par rapport à 𝑥, 𝑦 et 𝑡 sont :
𝜕𝑤3 𝜕2 𝑤3
= 𝑓 ′ cos 𝜃 ; = 𝑓 ′′ cos 2 𝜃
𝜕𝑥 𝜕𝑥 2
𝜕𝑤3 ′
𝜕2 𝑤3
= 𝑓 sin 𝜃 ; 2
= 𝑓 ′′ sin2 𝜃
𝜕𝑦 𝜕𝑦
𝜕𝑤3 ′
𝜕2 𝑤3
= −𝑐𝑓 ; 2
= 𝑐 2 𝑓 ′′
𝜕𝑥 𝜕𝑥
En substituant 𝑤3 à 𝑤 dans l’équation (3.2) et en utilisant les dérivées
partielles de 𝑤3 , alors on trouve que l’équation (3.2) est vérifiée.
Les solutions données par les équations (3.3) et (3.4) sont les mêmes que pour
une corde.
L’équation (3.3) exprime une onde progressive dans le sens des 𝑥 positifs à la
vitesse 𝑐. La forme de l’onde est indépendante de 𝑦 et la membrane se
comporte comme si elle était constituée d’un nombre infini de bandes tous
parallèles à l’axe 𝑥.
De façon similaire l’équation (3.4) exprime une onde progressive dans la
direction des 𝑥 négatifs. L’équation (3.5) exprime une onde se propageant
dans une direction faisant un angle 𝜃 à la vitesse 𝑐.
3.3 Vibrations libres des membranes rectangulaires :
La vibration libre d’une membrane rectangulaire de cotés 𝑎 et 𝑏 (voir Figure
5) peut être déterminée en utilisant la méthode de séparation des variables.
Figure 12 : Membrane rectangulaire
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 41
Dynamique des Structures
Alors le déplacement 𝑤(𝑥, 𝑦, 𝑡) est exprimé sous la forme d’un produit de trois
fonctions comme suit :
𝑤(𝑥, 𝑦, 𝑡) = 𝑊 (𝑥, 𝑦)𝑇 (𝑡) = 𝑋(𝑥 )𝑌(𝑦)𝑇 (𝑡) (3.6)
En substituant l’expression (3.6) dans l’équation (3.2) de la vibration libre et
en divisant l’expression obtenue par : 𝑋 (𝑥 )𝑌(𝑦)𝑇(𝑡), on obtient :
2
1 𝑑 2 𝑋 (𝑥) 1 𝑑 2 𝑌 ( 𝑦) 1 𝑑 2 𝑌 ( 𝑦) (3.7)
𝑐 [ + ] =
𝑋(𝑥 ) 𝑑𝑥 2 𝑌(𝑦) 𝑑𝑦 2 𝑌(𝑦) 𝑑𝑦 2
Le premier membre de l’équation (3.7) est une fonction des seules variables
𝑥 et 𝑦; le second membre est une fonction de la seule variable 𝑡, alors les deux
membres sont nécessairement égaux à une constante, 𝑘:
2
1 𝑑 2 𝑋 (𝑥 ) 1 𝑑 2 𝑌 ( 𝑦) 1 𝑑 2 𝑇 (𝑡) 2
(3.8)
𝑐 [ + ] = = 𝑘 = −𝜔
𝑋 (𝑥 ) 𝑑𝑥 2 𝑌(𝑦) 𝑑𝑦 2 𝑇 (𝑡) 𝑑𝑡 2
L’équation (3.8) nous donne deux équations séparées :
1 𝑑 2 𝑋 (𝑥 ) 𝜔2 1 𝑑 2 𝑌 ( 𝑦) (3.9)
+ = −
𝑋 (𝑥 ) 𝑑𝑥 2 𝑐2 𝑌 ( 𝑦) 𝑑 𝑦 2
2 ( ) (3.10)
𝑑 𝑇 𝑡 2 ( )
+ 𝜔 𝑇 𝑡 = 0
𝑑𝑡 2
Pour l’équation (3.9) on peut appliquer un raisonnement similaire à celui de
l’équation (3.7).
Soit
1 𝑑 2 𝑋 (𝑥 ) 𝜔2 1 𝑑 2 𝑌 ( 𝑦)
+ 2 =− = 𝛽2
𝑋 (𝑥 ) 𝑑𝑥 2 𝑐 𝑌(𝑦) 𝑑𝑦 2
Et posons :
𝜔2 (3.11)
𝛼 = 2 − 𝛽2
2
𝑐
On obtient les équations (3.11) et (3.12) qui suivent :
𝑑 2 𝑋 (𝑥 ) (3.12)
+ 𝛼2 𝑋 (𝑥 ) = 0
𝑑𝑥 2
𝑑 2 𝑌 ( 𝑦) (3.13)
+ 𝛽 2 𝑌 ( 𝑦) = 0
𝑑𝑦 2
Donc le problème qui était de résoudre une équation différentielle du second
ordre de trois variables revient maintenant à la résolution de trois équations
différentielles du second ordre à une variable chacune, en l’occurrence les
équations (3.10) , (3.12) et (3.13). Les solutions peuvent donc s’exprimer sous
la forme suivante :
𝑇 (𝑡) = 𝐴 cos 𝜔𝑡 + 𝐵 sin 𝜔𝑡 (3.14)
𝑋(𝑥 ) = 𝐶1 cos 𝛼𝑥 + 𝐶2 sin 𝛼𝑥 (3.15)
( )
𝑌 𝑦 = 𝐶3 cos 𝛽𝑦 + 𝐶4 sin 𝛽𝑦 (3.16)
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 42
Dynamique des Structures
Où les constantes 𝐴 et 𝐵 sont déterminées par les conditions initiales et les
constantes 𝐶1 , 𝐶2 , 𝐶3 et 𝐶4 sont déterminées par les conditions aux limites.
On détermine finalement les solutions de l’équation (3.2) pour une membrane
rectangulaire à l’aide des équations (3.6), (3.14), (3.15) et (3.16).
3.3.1 Membrane fixée aux bordures :
Si une membrane est fixée de tous les côtés, les conditions aux limites
s’expriment par :
𝑤(0, 𝑦, 𝑡) = 0, 0 ≤ 𝑦 ≤ 𝑏, 𝑡 ≥ 0 (3.17)
𝑤 (𝑎, 𝑦, 𝑡) = 0, 0 ≤ 𝑦 ≤ 𝑏, 𝑡 ≥ 0 (3.18)
𝑤(𝑥, 0, 𝑡) = 0, 0 ≤ 𝑥 ≤ 𝑎, 𝑡 ≥ 0 (3.19)
𝑤(𝑥, 𝑏, 𝑡) = 0, 0 ≤ 𝑥 ≤ 𝑎, 𝑡 ≥ 0 (3.20)
L’équation (3.6) permet d’affirmer que les conditions aux limites deviennent :
𝑋(0) = 0 (3.21)
𝑋 ( 𝑎) = 0 (3.22)
𝑌 ( 0) = 0 (3.23)
𝑌(𝑏) = 0 (3.24)
Les équations (3.21) et (3.23) nécessitent que 𝐶1 = 0 et 𝐶3 = 0. Les fonctions
𝑋(𝑥 ) et 𝑌(𝑦) deviennent alors :
𝑋(𝑥 ) = 𝐶2 sin 𝛼𝑥
𝑌(𝑦) = 𝐶4 sin 𝛽𝑦
Pour avoir des solutions non triviales des fonctions 𝑋(𝑥 ) et 𝑌(𝑦), les
conditions 𝑋(0) = 0 et 𝑌(0) = 0 requièrent que :
sin 𝛼𝑥 = 0 (3.25)
sin 𝛽𝑦 = 0 (3.26)
Les équations (3.25) et (3.26) définissent les valeurs propres de la membrane
à travers l’équation (3.11). Les solutions des équations (3.25) et (3.26) sont
données par :
𝛼𝑚 = 𝑚𝜋 ; 𝑚 = 1, 2, … (3.27)
𝛽𝑛 = 𝑛𝜋 ; 𝑛 = 1, 2, … (3.28)
Les fréquences naturelles 𝜔𝑚𝑛 de la membrane peuvent être déterminées en
utilisant l’équation (3.11) par :
𝜔𝑚𝑛 2 = 𝑐 2 (𝛼𝑚 2 + 𝛽𝑛 2 )
𝑚 2 𝑛 2 1/2 (3.29)
𝜔𝑚𝑛 = 𝜋𝑐 √(𝛼𝑚 2 + 𝛽𝑛 2 ) = 𝜋𝑐 [( ) + ( ) ] , 𝑚 = 1, 2, … ; 𝑛 = 1, 2, …
𝑎 𝑏
La forme de mode 𝑊𝑚𝑛 de la membrane correspondant à la valeur propre
𝜔𝑚𝑛 est donné par :
𝑚𝜋𝑥 𝑛𝜋𝑦 𝑚𝜋𝑥 𝑛𝜋𝑦
𝑊𝑚𝑛 (𝑥, 𝑦) = 𝑋𝑚 (𝑥 )𝑌𝑛 (𝑦) = 𝐶2𝑚 sin 𝐶4𝑛 sin = 𝐶𝑚𝑛 sin sin
𝑎 𝑏 𝑎 𝑏
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 43
Dynamique des Structures
𝑚𝜋𝑥 𝑛𝜋𝑦 (3.30)
𝑊𝑚𝑛 (𝑥, 𝑦) = 𝐶𝑚𝑛 sin sin , 𝑚 = 1, 2, … ; 𝑛 = 1, 2, …
𝑎 𝑏
Le mode naturel de vibration correspondant à 𝜔𝑚𝑛 peut être exprimé par :
𝑚𝜋𝑥 𝑛𝜋𝑦
𝑤𝑚𝑛 (𝑥, 𝑦, 𝑡) = 𝐶𝑚𝑛 sin sin (𝐴 cos 𝜔𝑚𝑛 𝑡 + 𝐵 sin 𝜔𝑚𝑛 𝑡 )
𝑎 𝑏
Posons : 𝐴𝑚𝑛 = 𝐴. 𝐶𝑚𝑛 𝐵𝑚𝑛 = 𝐵. 𝐶𝑚𝑛
𝑚𝜋𝑥 𝑛𝜋𝑦
𝑤𝑚𝑛 (𝑥, 𝑦, 𝑡) = sin sin (𝐴𝑚𝑛 cos 𝜔𝑚𝑛 𝑡 + 𝐵𝑚𝑛 sin 𝜔𝑚𝑛 𝑡 ) (3.31)
𝑎 𝑏
La solution générale de l’équation (3.31) est donnée par la somme de tous
les modes naturels par :
∞ ∞
𝑚𝜋𝑥 𝑛𝜋𝑦 (3.32)
𝑤𝑚𝑛 (𝑥, 𝑦, 𝑡) = ∑ ∑ sin sin (𝐴𝑚𝑛 cos 𝜔𝑚𝑛 𝑡 + 𝐵𝑚𝑛 sin 𝜔𝑚𝑛 𝑡 )
𝑎 𝑏
𝑚=1 𝑛=1
Les constantes 𝐴𝑚𝑛 et 𝐵𝑚𝑛 peuvent être déterminées en utilisant les
conditions initiales.
∞ ∞
𝑚𝜋𝑥 𝑛𝜋𝑦 (3.33)
∑ ∑ 𝐴𝑚𝑛 sin sin = 𝑤0 (𝑥, 𝑦)
𝑎 𝑏
𝑚=1 𝑛=1
∞ ∞
𝑚𝜋𝑥 𝑛𝜋𝑦 (3.34)
∑ ∑ 𝐵𝑚𝑛 𝜔𝑚𝑛 sin sin = 𝑤̇0 (𝑥, 𝑦)
𝑎 𝑏
𝑚=1 𝑛=1
Les équations (3.33) et (3.34) expriment les doubles développements en série
sinusoïdales de Fourier des fonctions 𝑤0 (𝑥, 𝑦) et 𝑤̇0 (𝑥, 𝑦) respectivement. En
𝑚𝜋𝑥 𝑛𝜋𝑦
multipliant les équations (3.33) et (3.34) par sin 𝑎 sin 𝑏 et en intégrant sur la
surface de la membrane, on obtient les relations suivantes :
𝑎 𝑏
4 𝑚𝜋𝑥 𝑛𝜋𝑦
𝐴𝑚𝑛 = ∫ ∫ 𝑤0 (𝑥, 𝑦) sin sin 𝑑𝑥𝑑𝑦
𝑎𝑏 0 0 𝑎 𝑏
𝑎 𝑏
4 𝑚𝜋𝑥 𝑛𝜋𝑦
𝐵𝑚𝑛 = ∫ ∫ 𝑤̇0 (𝑥, 𝑦) sin sin 𝑑𝑥𝑑𝑦 .
𝑎𝑏𝜔𝑚𝑛 0 0 𝑎 𝑏
3.4 Vibrations forcées des membranes rectangulaires :
L’équation régissant la vibration forcée d’une membrane rectangulaire est
donnée par l’équation (3.1) qui est la suivante :
𝜕2 𝑤(𝑥, 𝑦, 𝑡) 𝜕2 𝑤(𝑥, 𝑦, 𝑡) 𝜕2 𝑤(𝑥, 𝑦, 𝑡) (3.35)
−𝑇 [ 2
+ 2
] + 𝜌 2
= 𝑓(𝑥, 𝑦, 𝑡)
𝜕𝑥 𝜕𝑦 𝜕𝑡
Nous pouvons trouver la solution de l’équation (3.35) en utilisant la
procédure de l’analyse modale. En conséquence nous supposons la solution
sous la forme :
∞ ∞
(3.36)
𝑤(𝑥, 𝑦, 𝑡) = ∑ ∑ 𝑊𝑚𝑛 (𝑥, 𝑦)𝜂𝑚𝑛 (𝑡)
𝑛=1 𝑚=1
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 44
Dynamique des Structures
où 𝑊𝑚𝑛 (𝑥, 𝑦) sont les modes naturels de la vibration et 𝜂𝑚𝑛 (𝑡) les
coordonnées généralisées correspondantes. Traitons le cas particulier d’une
membrane fixée aux bordures. Pour ce cas les fonctions propres sont :
𝑚𝜋𝑥 𝑛𝜋𝑦
𝑊𝑚𝑛 (𝑥, 𝑦) = 𝐶𝑚𝑛 sin sin , 𝑚 = 1, 2, … ; 𝑛 = 1, 2, …
𝑎 𝑏
Les fonctions propres ou modes normaux peuvent être normalisées par :
𝑎 𝑏 𝑎 𝑏
𝑚𝜋𝑥 2 𝑛𝜋𝑦
∫ ∫ 𝜌𝑊𝑚𝑛 2 (𝑥, 𝑦)𝑑𝑥𝑑𝑦 = 1 𝑠𝑜𝑖𝑡 ∫ ∫ 𝜌𝐶𝑚𝑛 2 sin2 sin 𝑑𝑥𝑑𝑦 = 1
0 0 0 0 𝑎 𝑏
La simplification donne :
2
𝐶𝑚𝑛 =
√𝜌𝑎𝑏
Alors les modes normaux prennent la forme suivante :
2 𝑚𝜋𝑥 𝑛𝜋𝑦 (3.37)
𝑊𝑚𝑛 (𝑥, 𝑦) = sin sin , 𝑚 = 1, 2, … ; 𝑛 = 1, 2, …
√𝜌𝑎𝑏 𝑎 𝑏
En substituant l’équation (3.36) dans l’équation (3.35) puis en multipliant le
résultat par 𝑊𝑚𝑛 (𝑥 ) et en intégrant sur la surface de la membrane, on
obtient :
𝜂̈ 𝑚𝑛 (𝑡) + 𝜔𝑚𝑛 2 𝜂𝑚𝑛 (𝑡) = 𝑁𝑚𝑛 (𝑡), 𝑚, 𝑛 = 1, 2, … (3.38)
Où :
𝑇 𝑚 2 𝑛 2
𝜔𝑚𝑛 = 𝜋√ [( ) + ( ) ]
𝜌 𝑎 𝑏
𝑎 𝑏 𝑎 𝑏
2 𝑚𝜋𝑥 𝑛𝜋𝑦
𝑁𝑚𝑛 (𝑡) = ∫ ∫ 𝑊𝑚𝑛 (𝑥, 𝑦)𝑓(𝑥, 𝑦, 𝑡) 𝑑𝑥𝑑𝑦 = ∫ ∫ 𝑓(𝑥, 𝑦, 𝑡) sin sin 𝑑𝑥𝑑𝑦
0 0 √𝜌𝑎𝑏 0 0 𝑎 𝑏
La solution de l’équation (3.38) peut être donnée par :
𝑡
1 𝜂̇ (0)
( )
𝜂𝑚𝑛 𝑡 = ∫ 𝑁𝑚𝑛 (𝜏) sin 𝜔𝑚𝑛 (𝑡 − 𝜏) 𝑑𝜏 + 𝜂𝑚𝑛 (0) cos 𝜔𝑚𝑛 𝑡 + 𝑚𝑛 sin 𝜔𝑚𝑛 𝑡
𝜔𝑚𝑛 0 𝜔𝑚𝑛
L’équation (3.36) devient alors :
∞ ∞
2 𝑚𝜋𝑥 𝑛𝜋𝑦 𝜂̇ (0)
𝑤 (𝑥, 𝑦, 𝑡) = { ∑ ∑ (sin sin ) [𝜂𝑚𝑛 (0) cos 𝜔𝑚𝑛 𝑡 + 𝑚𝑛 sin 𝜔𝑚𝑛 𝑡]}
√𝜌𝑎𝑏 𝑛=1 𝑚=1 𝑎 𝑏 𝜔𝑚𝑛
∞ ∞
𝑡
2 1 𝑚𝜋𝑥 𝑛𝜋𝑦
+[ ∑∑ (sin sin ) × ∫ 𝑁𝑚𝑛 (𝜏) sin 𝜔𝑚𝑛 (𝑡 − 𝜏) 𝑑𝜏 +]
√𝜌𝑎𝑏 𝑛=1 𝑚=1 𝜔𝑚𝑛 𝑎 𝑏 0
3.5 Vibration des membranes circulaires :
En notant le Laplacien en coordonnées cartésiennes est définie par :
𝜕2 𝜕2
∇= +
𝜕𝑥 2 𝜕𝑦 2
L’équation (3.1) du mouvement d’une membrane rectangulaire peut être
réécrite comme suit :
𝜕2 𝑤(𝑥, 𝑦, 𝑡)
𝑇 ∇ 𝑤(𝑥, 𝑦, 𝑡) + 𝑓(𝑥, 𝑦, 𝑡) = 𝜌
𝜕𝑡 2
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 45
Dynamique des Structures
Pour une membrane circulaire l’équation régissant le mouvement peut être
déduit en utilisant une approche d’équilibre en considérant un élément
différentiel en coordonnées polaire 𝑟 et 𝜃.
Les relations de transformation de coordonnées peuvent être utilisées pour
exprimer le Laplacien en coordonnées polaires :
𝑥 = 𝑟 cos 𝜃 et 𝑦 = 𝑟 sin 𝜃
Le Laplacien est alors :
1 𝜕 𝜕 1 𝜕² 𝜕2 1𝜕 1 𝜕²
∇= (𝑟 ) + 2 2 = 2 + + 2 2
𝑟 𝜕𝑟 𝜕𝑟 𝑟 𝜕𝜃 𝜕𝑟 𝑟 𝜕𝑟 𝑟 𝜕𝜃
Alors l’équation du mouvement pour la vibration forcée d’une membrane
circulaire peut etre exprimée par :
𝜕2 𝑤(𝑟, 𝜃, 𝑡) 1 𝜕𝑤(𝑟, 𝜃, 𝑡) 1 𝜕²𝑤(𝑟, 𝜃, 𝑡) 𝜕²𝑤(𝑟, 𝜃, 𝑡)
𝑇 [ + + ] + 𝑓 ( 𝑟, 𝜃, 𝑡) = 𝜌
𝜕𝑟2 𝑟 𝜕𝑟 𝑟2 𝜕𝜃 2 𝜕𝑡 2
Pour une vibration libre l’équation devient :
2
𝜕2 𝑤(𝑟, 𝜃, 𝑡) 1 𝜕𝑤(𝑟, 𝜃, 𝑡) 1 𝜕²𝑤(𝑟, 𝜃, 𝑡) 𝜕²𝑤(𝑟, 𝜃, 𝑡)
𝑐 [ + + ] + 𝑓 ( 𝑟, 𝜃, 𝑡 ) = 𝜌
𝜕𝑟2 𝑟 𝜕𝑟 𝑟2 𝜕𝜃 2 𝜕𝑡 2
Où 𝑐 = √𝑡/𝜌.
Le déplacement 𝑤 est maintenant une fonction de 𝑟, 𝜃 et 𝑡 on peut utiliser la
méthode séparation des variables pour la résolution de l’équation, en posant :
𝑤(𝑟, 𝜃, 𝑡) = 𝑅(𝑟)Θ(𝜃)𝑇 (𝑡)
Où 𝑅(𝑟), Θ(𝜃) et 𝑇 (𝑡) sont des fonctions des seules variables
𝑟, 𝜃 et 𝑡 respectivement.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 46
Dynamique des Structures
Chapitre 3 : Méthodes Variationnelle de
Caractérisation des valeurs propres
Introduction
Nous avons vu que, dans beaucoup de problèmes d’études de systèmes
dynamiques, la connaissance des valeurs propres nous permet de déterminer
les fréquences de vibration de ces systèmes. Cependant la détermination
analytique de ces valeurs est dans de nombreux cas fastidieux voire
impossible. Face à ces cas, l’ingénieur dispose méthodes numériques
destinées à déterminer des valeurs approximatives des valeurs propres et
vecteurs propres du système. Ces méthodes donnent de bonnes
approximations de ces valeurs propres. Dans cette leçon, nous présentons
quelques algorithmes de valeurs propres classiques. Individuellement,
ces outils sont utiles dans certaines situations, en particulier l'itération
inverse, qui est la méthode standard pour déterminer un vecteur propre
lorsque la valeur propre correspondante est connue.
Soit 𝐴 une matrice carrée d’ordre 𝑛.
Le problème de calcul des valeurs propres consiste à trouver un scalaire 𝜆 tel
qu’il existe un vecteur 𝑦, 𝑦 ≠ 0, tel que : 𝑨𝒚 = 𝝀𝒚.
Un tel vecteur 𝑦 est appelé vecteur propre associé à la valeur propre 𝜆.
Dans la suite on notera 𝜆1 , 𝜆 2 , … , 𝜆𝑛 les valeurs propres de la matrice 𝐴.
1. Quotient de Rayleigh :
Nous savons que l’équation du mouvement d’un système dynamique linéaire
s’écrit sous la forme suivante :
[𝑀]{𝑋̈ } + [𝐶]{𝑋̇ } + [𝐾]{𝑋} = {𝐹 }
où [𝑀], [𝐶] et [𝐾] sont respectivement les matrices de masse,
d’amortissement et de raideur et {𝑋} et {𝐹 } les vecteurs déplacement et forces
extérieures appliquées au système respectivement.
Pour un système dynamique conservatif et sans amortissement, cette
équation devient :
[𝑀]{𝑋̈ } + [𝐾]{𝑋} = 0
L’équation aux valeurs propres associée à cette équation est :
(𝐾 − 𝜔2 𝑀)𝑈 = 0 ou bien (𝐾 − 𝜆𝑀)𝑈 = 0
où 𝜔2 est une valeur propre asociee au vecteur propre 𝑈.
Pour un tel système et pour un vecteur {𝑋} quelconque, l’énergie cinétique de
déformation 𝐾 et l’énergie potentielle de déformation 𝑇 sont definies par :
𝐾 = {𝑋} 𝑇 [𝑀]{𝑋} et 𝑇 = {𝑋} 𝑇 [𝐾]{𝑋}
On définit, pour tout vecteur {𝑋} non nul, le quotient de Rayleigh 𝑅𝑄({𝑋}) du
système par le ratio de l’énergie potentielle par rapport à l’énergie cinétique :
{𝑋} 𝑇 [𝐾]{𝑋}
𝑅𝑄 ({𝑋}) =
{𝑋} 𝑇 [𝑀]{𝑋}
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 47
Dynamique des Structures
L’intérêt principal du quotient de Rayleigh est que pour une valeur propre
𝜔2 du système associée au vecteur propre 𝑈, on a :
𝑅𝑄({𝑈}) = 𝜔2 = 𝜆
Donc une fois un vecteur propre identifié, on détermine immédiatement la
valeur propre correspondante en calculant le quotient de Rayleigh.
Le principe classique de Rayleigh stipule que dans l’espace réel, le quotient
de Rayleigh atteint une valeur stationnaire à chaque mode propre et est
délimité par le carré des fréquences naturelles les plus basses et les plus
élevées du système.
Pour un système linéaire amorti, l’équation du mouvement est la suivante :
[𝑀]{𝑋̈ } + [𝐶]{𝑋̇ } + [𝐾]{𝑋} = 0
Dans ce cas les valeurs propres sont les solution de l’équation quadratique :
([𝑀]𝜆2 + [𝐶 ]𝜆 + [𝐾])𝑈 = 0
avec 𝜆 une valeur propre et 𝑈 un vecteur propre qui sont tous de dimension
𝑛, et généralement complexes. La solution de l'équation donne 2𝑛 valeurs
propres 𝜆𝑗 (𝑗 = 1,2, . . . ,2𝑛) et leurs vecteurs propres 𝑈𝑗 correspondants, qui ne
peuvent être déterminés qu’à un scalaire multiplicateur arbitraire près. Le
spectre d'un système linéaire amorti est composé de valeurs propres réelles
négatives, associées à une décroissance exponentielle de la réponse, et de
valeurs propres conjuguées complexes, dont la partie réelle négative est un
taux de décroissance et dont la partie imaginaire est une fréquence de
vibration amortie. Les valeurs propres réelles peuvent toujours avoir des
vecteurs propres réels qui leur sont associés.
Le système dissipatif est associé à trois type de quotients de Rayleigh :
{𝑍} ∗ [𝐴]{𝑍}
𝑅𝑄𝐴𝐵 ({𝑍}) =
{𝑍} ∗ [𝐵 ]{𝑍}
Calcule en remplaçant les matrices réelles, symétriques et définies positives
[𝐴] et [𝐵 ] par deux quelconques des matrices [𝑀], [𝐶 ] et [𝐾].
Pour simplifier, les quotients de la forme de l'équation (17) seront appelés
quotients de Rayleigh. Comme les vecteurs propres de l'équation (16) sont
généralement complexes, les quotients de l'équation (17) sont fonctions de
vecteurs complexes à 𝑛 dimensions, où {𝑍} ∗ désigne la transposée
conjuguée complexe : {𝑍} ∗ = ̅̅̅̅ {𝑍} 𝑇 .
Bien que {𝑍} puisse être complexe, le quotient 𝑅𝑄𝐴𝐵 ({𝑍}) est réel.
A chacun des trois types de quotient de Rayleigh est associée la valeur propre
généralisée :
([𝐴] + 𝜇𝐴𝐵 [𝐵 ]){𝑊 }𝐴𝐵 = 0
avec la valeur propre 𝜇𝐴𝐵 associée au vecteur propre de dimension 𝑛 {𝑊 }𝐴𝐵 .
Toutes les valeurs propres sont des nombres réels négatifs alors elles peuvent
être rangées dans un ordre de valeur absolue croissante :
0 ≤ |𝜇𝐴𝐵 𝐴𝐵 𝐴𝐵
1 | ≤ |𝜇2 | ≤ |𝜇3 | ≤ ⋯ ≤ |𝜇𝑛 |
𝐴𝐵
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 48
Dynamique des Structures
De plus les vecteurs propres correspondants {𝑊 }𝑗𝐴𝐵 sont linéairement
indépendant et orthogonal par rapport à [𝐴] et par rapport à [𝐵 ].
2. Méthode de la Puissance itérée :
Soit 𝐴 une matrice symétrique définie positive, supposons la matrice
𝐴 diagonalisable, c’est-à-dire pouvant se mettre sous la forme :
𝐴 = 𝑃 −1 𝐷𝑃
Avec 𝐷 une matrice diagonale et 𝑃 une matrice inversible d’ordre 𝑛.
Supposons les valeurs propres de 𝐴, 𝜆 1 , 𝜆 2 , … , 𝜆 𝑛 , réelles et distinctes en
module :
|𝜆 1 | > |𝜆 2 | > ⋯ > |𝜆 𝑛 |
Valeur propre dominante-Vecteur propre dominant :
Soit 𝜆 1 , 𝜆 2 , … , 𝜆 𝑛 les valeurs propres d’une matrice carrée 𝐴 d’ordre 𝑛.
𝜆 1 est appelé valeur propre dominante si |𝜆 1 | > |𝜆 𝑖 | pour 𝑖 = 2, 3, … . , 𝑛
Le vecteur propre correspondant à 𝜆 1 est appelé vecteur propre dominant.
La méthode la puissance est une méthode itérative. Sa première hypothèse est
que la matrice 𝐴 admet une valeur propre dominante avec des vecteurs
propres dominants correspondants.
On détermine alors une approximation initiale 𝒙𝟎 de l’un des vecteurs propres
dominants. Cette approximation initiale doit être un vecteur non nul de ℝ𝑛 .
Finalement on forme la séquence suivante :
𝑥 1 = 𝐴𝑥 0
𝑥 2 = 𝐴𝑥 1 = 𝐴2 𝑥 0
𝑥 3 = 𝐴𝑥 2 = 𝐴(𝐴2 𝑥 0 ) = 𝐴3 𝑥 0
……..
𝑥 𝑘 = 𝐴𝑥 𝑘−1 = 𝐴(𝐴𝑘−1 𝑥 0 ) = 𝐴𝑘 𝑥 0
Dès que l’on détermine une bonne approximation d’un vecteur propre
dominant, le quotient de Rayleigh nous permet de déterminer la valeur propre
dominante correspondante.
Condition suffisante de convergence de la méthode la puissance :
Si 𝐴 est une matrice carrée d’ordre 𝑛 diagonalisable avec une valeur propre
dominante alors il existe un vecteur non nul 𝒙𝟎 tel que la suite
𝐴𝑥 0 ; 𝐴2 𝑥 0 ; 𝐴3 𝑥 0 ; … ; 𝐴𝑘 𝑥 0 ; … converge vers un multiple du vecteur propre
dominant de 𝐴.
Exemple d’Application : Déterminer la plus grande valeur propre de la
matrice suivante :
1 6 1
𝐴 = ( 1 2 0)
0 0 3
1
Posons 𝑋0 = (0) le vecteur initial.
0
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 49
Dynamique des Structures
1
1 1 6 1 1
𝑋1 = 𝐴𝑋0 = (1
0) (0) = (1) = 1 (1)
2
0
3 0 0 0 0
1
1 6 1 1 7 3 1
𝑋2 = 𝐴𝑋1 = (1 2 0) (1) = (3) = 7 ( ) = 7 (0.42)
0 7
0 0 3 0 0
0
1 6 1 1 3.52 1
𝑋3 = 𝐴𝑋2 = (1 2 0) (0.42) = (1.54) = 3.52 (0.52)
0 0 3 0 0 0
1 6 1 1 4.12 1
𝑋4 = 𝐴𝑋3 = (1 2 0) (0.52) = (2.04) = 4.12 (0.49)
0 0 3 0 0 0
1 6 1 1 3.94 1
𝑋5 = 𝐴𝑋4 = (1 2 0) (0.49) = (1.98) = 3.94 (0.50)
0 0 3 0 0 0
1 6 1 1 4 1
𝑋6 = 𝐴𝑋5 = (1 2 0) (0.50) = (2) = 4 (0.50)
0 0 3 0 0 0
On voit que 𝑋6 − 𝑋5 = 0, alors la plus grande valeur propre est 4 et le vecteur
1
propre correspondant est (0.50) .
0
Cette méthode de la puissance permet également de trouver la plus petite
valeur propre, dans ce cas il faut appliquer la méthode avec la matrice 𝐴−1 .
I. Méthodes d’itération vectorielle
1. Concepts de base :
Considérons le problème aux valeurs propres suivant : [𝑨]{𝑼} = 𝝀{𝑼}.
Soient 𝜆1 , 𝜆 2 , … … , 𝜆 𝑛 les valeurs associees respectivement aux vecteurs
propres {𝑈1 }, {𝑈2 }, … … , {𝑈𝑛 } . Nous savons que ces vecteurs propres constituent
une base, alors pour un vecteur quelconque, on peut écrire (𝑐1 , 𝑐2 , … , 𝑐𝑛 étant
des scalaires) :
{ 𝑋 1 } = 𝑐1 {𝑈1 } + 𝑐2 {𝑈2 } + ⋯ + 𝑐𝑛 {𝑈𝑛 }.
Calculons le vecteur {𝑋 2 } = [𝐴]{𝑋 1 }, on a :
{𝑋 2 } = [𝐴]{𝑋 1 } = 𝑐1 [𝐴]{𝑈1 } + 𝑐2 [𝐴] {𝑈2 } + ⋯ + 𝑐𝑛 [𝐴]{𝑈𝑛 }
{𝑋 2 } = 𝑐1 𝜆1 {𝑈1 } + 𝑐2 𝜆 2 {𝑈2 } + ⋯ + 𝑐𝑛 𝜆 𝑛 {𝑈𝑛 }
En multipliant, à gauche, par [𝐴] 𝑚 fois, on obtient :
{𝑋 𝑚+1 } = 𝑐1 𝜆𝑚 𝑚 𝑚
1 {𝑈1 } + 𝑐2 𝜆 2 {𝑈2 } + ⋯ + 𝑐𝑛 𝜆 𝑛 {𝑈𝑛 }
𝜆2 𝑚 𝜆𝑛 𝑚
{𝑋 𝑚+1 } = 𝑐1 𝜆𝑚 1 [{𝑈1 } + 𝑐 2 ( ) { 𝑈2 } + ⋯ + 𝑐𝑛 ( ) {𝑈𝑛 }]
𝜆1 𝜆1
En classant les valeurs propres par ordre décroissant tel que : 𝜆 1 > 𝜆 2 >
𝜆 𝑚
… … > 𝜆 𝑛 alors après plusieurs itérations, on a : (𝜆 𝑖 ) → 0 pour 𝑖 =
1
1, 2, … , 𝑛, i. e. on a :
{𝑋 𝑚+1 } ≈ 𝑐1 𝜆𝑚
1 {𝑈1 }
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 50
Dynamique des Structures
1
Par exemple, pour une valeur propre égale a , apres 10 itérations, on a :
5
1 10
(5 ) = 0,0000001024. Ce qui est sensiblement négligeable. Dons si les valeurs
propres sont bien séparées, le vecteur {𝑋 𝑚+1 } contiendra, essentiellement, les
composantes du mode propre {𝑈1 }. On dit que {𝑋 𝑚+1 } converge vers {𝑈1 }.
A l’aide du quotient de Rayleigh nous obtenons une estimation de 𝜆 1 .
A noter que, indépendamment de notre choix de {𝑋 1 } la convergence est
toujours vers la plus grande valeur propre et le vecteur propre correspondant.
La rapidité de convergence dépend de comment les valeurs propres sont bien
séparées. Si une valeur propre est proche de la plus grande valeur propre,
alors même après plusieurs itérations, sa contribution ne sera pas
négligeable. Pour accélérer la convergence nous devons séparer les valeurs
propres et pour cela nous pouvons utiliser la propriété de substitution des
valeurs propres.
2. Itération directe et inverse :
Nous pouvons appliquer ce que nous venons de voir à notre équation aux
valeurs propres typique :
([𝐾] − 𝜔2 [𝑀]){𝑈} = 0 soit [𝐾]{𝑈} = 𝜔2 [𝑀]{𝑈}
- Itération directe :
[𝐾]{𝑈} = 𝜔2 [𝑀]{𝑈} soit [𝑀]−1 [𝐾]{𝑈} = 𝜔2 {𝑈}
- Itération inverse :
1
[𝐾]{𝑈} = 𝜔2 [𝑀]{𝑈} soit {𝑈} = [𝐾]−1 [𝑀]{𝑈}
𝜔2
Dans chaque cas une estimation de la valeur absolue peut être obtenue en
calculant le quotient de Rayleigh 𝑅𝑄 défini par :
{𝑋} 𝑇 [𝐾]{𝑋}
𝑅𝑄 =
{𝑋} 𝑇 [𝑀]{𝑋}
Le quotient de Rayleigh est base sur le fait que :
[𝐾]{𝑈} = 𝜔2 [𝑀]{𝑈} ⟹ {𝑈} 𝑇 [𝐾]{𝑈} = 𝜔2 {𝑈} 𝑇 [𝑀]{𝑈}
Soit :
{𝑈} 𝑇 [𝐾]{𝑈}
𝜔2 =
{𝑈} 𝑇 [𝑀]{𝑈}
Ceci est, en réalité, une mesure du rapport de l’énergie potentielle sur
l’énergie cinétique. Nous savons que l’amplitude du vecteur {𝑈} n’est pas
unique. Pour fixer l’amplitude du vecteur {𝑈} , on effectue la normalisation
{𝑈}𝑇 [𝑀]{𝑈} = 1. Alors chaque itération calcule simplement le 𝑅𝑄, et continue
le processus d’itération jusqu’à ce que la valeur absolue converge dans la
limite de tolérance prédéfinie.
Schéma de calcul avec la méthode de l’itération inverse :
Etape 1 : Formuler les matrices [𝑀] 𝑒𝑡 [𝐾] de la structure
Etape 2 : Choisir un vecteur initial {𝑋 1 }
Etape 3 : calculer {𝑅} = [𝑀]{𝑋 1 }
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 51
Dynamique des Structures
Etape 4 : Résoudre [𝐾]{𝑋̅} = {𝑅}
Etape 5 : calculer {𝑋 2 } à partir de {𝑋̅} tel que {𝑋 2 }𝑇 [𝑀]{𝑋 2 } = 1
i.e.
1
{𝑋 2 } = {𝑋̅}
√{𝑋̅}𝑇 [𝑀]{𝑋̅}
Etape 6 : calculer 𝜆 = {𝑋} 𝑇 [𝐾]{𝑋}
Répéter les étapes 3 à 6 jusqu’à ce que 𝜆 converge dans la tolérance
prédéfinie.
Le schéma de l’itération inverse est un outil très efficace pour la
détermination de la première paire propre d’un système.
3. Déflation de Gram-Schmidt
Nous avons : [𝐾]{𝑈} = 𝜔2 [𝑀]{𝑈}
Soit 𝜆1 , { 𝑈1 } le plus petit couple propre, i.e. 𝜆1 est la plus petite valeur propre
puis à partir d’un vecteur quelconque {𝑋2 }que nous choisissons pour obtenir
le vecteur propre {𝑈2 }, soustraire 𝛽1 {𝑈1 }, où :
𝛽1 = {𝑈1 }𝑇 [𝑀]{𝑋2 }
{𝑋̅2 } = {𝑋2 } − 𝛽1 {𝑈1 }
Alors, après une itération appropriée, {𝑋̅2 } convergera vers {𝑈2 }.
En effet, posons :
{𝑋2 } = 𝑐1 {𝑈1 } + 𝑐2 {𝑈2 } + 𝑐3 {𝑈3 } + ⋯ + 𝑐𝑛 { 𝑈𝑛 }
Alors :
{𝑈1 }𝑇 [𝑀]{𝑋2 } = 𝑐1 {𝑈1 }𝑇 [𝑀]{𝑈1 } + 𝑐2 {𝑈1 }𝑇 [𝑀] {𝑈2 } + 𝑐3 {𝑈1 }𝑇 [𝑀] {𝑈3 } + ⋯ + 𝑐𝑛 {𝑈1 }𝑇 [𝑀]{𝑈𝑛 }
La propriété d’orthogonalité des vecteurs propres implique que :
{𝑈1 }𝑇 [𝑀]{𝑈𝑗 } = 0 , pour 𝑗 = 2, 3, … , 𝑛
{𝑈1 }𝑇 [𝑀]{𝑈1 } = 1
Alors
𝛽1 = {𝑈1 }𝑇 [𝑀]{𝑋2 } = 𝑐1
Alors :
{𝑋̅2 } = {𝑋2 } − 𝛽1 {𝑈1 } = 𝑐1 {𝑈1 } + 𝑐2 {𝑈2 } + 𝑐3 {𝑈3 } + ⋯ + 𝑐𝑛 {𝑈𝑛 } − 𝑐1 {𝑈1 }
= 𝑐2 {𝑈2 } + 𝑐3 {𝑈3 } + ⋯ + 𝑐𝑛 {𝑈𝑛 }
Nous avons ainsi éliminé la composante {𝑈1 } du vecteur initial {𝑋̅2 }.
Dans sa compilation dans l’itération inverse, nous devons effectue cette
déflation de Gram-Schmidt dans chaque itération. Notre vecteur initial
convergera alors vers le couple propre 𝜆 2 , {𝑈2 }. Apres avoir déterminer
𝜆 2 , {𝑈2 } , pour déterminer 𝜆 3 , { 𝑈3 } nous utiliserons alors le vecteur déflaté
donne par :
{𝑋̅} = {𝑋} − 𝛽1 {𝑈1 } − 𝛽2 {𝑈2 }
Où 𝛽2 est déterminé de façon similaire a 𝛽1 .
Séquences de l’itération inverse avec la déflation de Gram-Schmidt :
Etape 1 : Formuler les matrices [𝑀] 𝑒𝑡 [𝐾] du système
Etape 2 : Choisir un vecteur initial {𝑋1 }
Etape 3 : Effectuer la déflation, i.e.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 52
Dynamique des Structures
𝑚
{̅̅̅
𝑋1 } = {𝑋1 } − ∑ 𝛽𝑖 {𝑈𝑖 }
𝑖=1
où 𝛽𝑖 = { 𝑈𝑖 }𝑇 [𝑀]{𝑋1 } et m le nombre de pairs propres déjà trouvés.
Etape 4 : Calculer {𝑅} = [𝑀]{̅̅̅̅̅̅
𝑋1 }
Etape 5 : Résoudre l’équation [𝑀]{̿̿̿̿̿̿ 𝑋1 } = {𝑅}
Etape 6 : Obtenir {𝑋2 } à partir de {̿̿̿̿̿̿
𝑋1 } tel que :
{𝑋2 }𝑇 [𝑀]{𝑋2 } = 1
1
{𝑋2 } = √
̿̿̿̿̿̿
{𝑋2 }𝑇 [𝑀]{̿̿̿̿̿̿
𝑋2 }
Etape 7 : calculer la (𝑚 + 1)𝑒 valeur propre 𝜆𝑚+1 : 𝜆𝑚+1 = { 𝑋2 }𝑇 [𝐾]{𝑋2 }.
Répéter les 3 à 7 jusqu’à ce qu’une convergence satisfaisante soit atteinte.
4. Méthode de l’itération simultanée :
Dans la technique de la déflation de Gram-Schmidt chaque vecteur initia est
redu orthogonal a tous les vecteurs propres obtenus précédemment. Ainsi
nous nous assurons de la convergence vers le couple propre était obtenue de
manière séquentielle, en commençant par le premier jusqu’au 𝑚𝑒 vecteur
désiré. Dans la présente méthode nous commençons un ensemble de vecteurs
initiaux et faire l’itération sur eux simultanément.
Ce schéma d’itération est similaire à la méthode d’itération sur un vecteur de
base, mais avec plusieurs vecteurs initiaux. L’idée de base est d’assurer que
les vecteurs initiaux sont mutuellement orthogonaux dans chaque itération
pour qu’ils convergent vers des vecteurs propres distincts. Nous avons :
[𝐾]{𝑈} = 𝜔2 [𝑀]{𝑈}
Les étapes de base dans cette méthode sont les suivantes :
Etape 1 : Choisir les vecteurs {̅̅̅ 𝑋1 }𝑛× 𝑚 . Ces vecteurs doivent être deux a deux
orthogonaux. Supposons que d’une manière quelconque nous puissions extraire
{𝑋1 }𝑛×𝑚 , i.e. :
{𝑋1 }𝑇 [𝑀]{𝑋1 } = [𝐼 ]
Etape 2 : Calculer {𝑅}𝑛×𝑚 = [𝑀]𝑛×𝑛 { 𝑋1 }𝑛×𝑚 .
Etape 3 : Déterminer {̅̅̅ 𝑋2 }𝑛×𝑚 en resolvant l’équation [𝐾]𝑛×𝑛 {̅̅̅
𝑋2 }𝑛×𝑚 = {𝑅} 𝑛×𝑚
Etape 4 : Calculer [𝜆]𝑚×𝑚 = {𝑋1 } [𝐾]{𝑋1 } . [𝜆] est maintenant une matrice
𝑇
carrée d’ordre m dont les éléments diagonaux sont les approximations de nos
valeurs propres. On continue les itérations jusqu’à ce que [𝜆] converge
complètement, [𝜆] sera de la forme :
𝜆1 0
[𝜆] = [ 𝜆 2 ]
⋱
0 𝜆𝑛
Et les {𝑋} correspondant auront les vecteurs propres.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 53
Dynamique des Structures
Nous voyons que les étapes de base sont identiques a celles de l’itération
inverse, excepté que nous opérons sur plusieurs vecteurs initiaux.
Cette différence significative dans la procédure pour obtenir des vecteurs
orthonormalisés {𝑋} à partir des vecteurs originaux {𝑋̅} est maintenant
décrite.
Pour avoir {𝑋} à partir {𝑋̅}, nous utilisons la condition suivante :
{𝑋} 𝑇𝑚×𝑛 [𝑀]𝑛×𝑛 {𝑋} 𝑛×𝑚 = [𝐼 ]𝑛×𝑛
Soit
{𝑋} 𝑛×𝑚 = {𝑋̅}𝑛× 𝑚 [𝑆]𝑚×𝑚
Alors on a (en omettant les dimensions) :
{𝑋} 𝑇 [𝑀]{𝑋} = [𝑆]𝑇 {𝑋̅}𝑇 [𝑀]{𝑋̅}[𝑆] = [𝑆]𝑇 [𝐷][𝑆] = [𝐼 ]
où
[𝐷] = {𝑋̅}𝑇 [𝑀]{𝑋̅}
Effectuons la factorisation Cholesky de [𝐷], posons :
[𝐷] = [𝐿 ][𝐿]𝑇
Alors, on obtient que :
[𝑆] = [𝐿 ]−𝑇
Ainsi, notre méthode pour trouver les vecteurs orthonormés {𝑋} à partir
{𝑋̅} implique l’utilisation de la matrice triangulaire inferieure [𝐿] obtenue par
la factorisation de Cholesky de {𝑋̅} 𝑇 [𝑀]{𝑋̅}. Nous pouvons ainsi écrire :
{𝑋} = {𝑋̅}[𝐿 ]−𝑇 ou {𝑋} [𝐿]𝑇 = {𝑋̅}
En évitant le calcul explicite de [𝐿]−1 nous préférons obtenir {𝑋} comme
solution en remontant.
Séquences de la méthode de l’itération simultanée :
Etape 1 : Formuler les matrices globales assemblées de raideur et de masse
[𝐾] et [𝑀] .
Etape 2 : Choisir les vecteurs {̅̅̅
𝑋1 }𝑛× 𝑚 .
Etape 3 : Calculer [𝐷]𝑚×𝑚 = {̅̅̅𝑋1 }𝑇𝑛×𝑚 [𝑀]{̅̅̅
𝑋1 }𝑛×𝑚
Etape 4 : Déterminer la factorisation de Cholesky [𝐷]𝑚×𝑚 = [𝐿 ]𝑇𝑚×𝑚 [𝐿 ]𝑚×𝑚
Etape 5 : Déterminer {𝑋1 }𝑛×𝑚 en résolvant l’équation {𝑋1 }𝑛× 𝑚 = [𝐿 ]𝑇𝑚×𝑚 {̅̅̅
𝑋1 }𝑛×𝑚
Etape 6 : Calculer 𝑅 𝑛×𝑚 = 𝑀 𝑛×𝑛 𝑋1 𝑛×𝑚
[ ] [ ] { }
Etape 7 : Déterminer {̅̅̅
𝑋2 }𝑛×𝑚 en résolvant l’équation [𝐾]𝑛×𝑛 {̅̅̅
𝑋2 }𝑛×𝑚 = [𝑅]𝑛×𝑚
Etape 8 : Calculer la matrice quotient de Rayleigh [𝜆] = {𝑋1 }𝑇𝑛×𝑚 [𝐾]𝑛×𝑛 {𝑋1 }𝑛×𝑚
Répéter les étapes 3 à 8 jusqu’à ce qu’une convergence satisfaisante soit
obtenue en [λ], i.e. les éléments non diagonaux soit très négligeables et que
les changements dans les valeurs d’une itération à une autre soit très
négligeables.
Factorisation de Cholesky :
On se propose, pour une matrice [𝐴] donnée, d’obtenir une factorisation sous
la forme :
[𝐴] = [𝐿 ][𝐿 ]𝑇
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 54
Dynamique des Structures
Où [𝐿 ] est une matrice triangulaire inférieure, i.e. [𝐿 ] est de forme :
𝑙 11 0
𝑙 22
𝑙
[𝐿 ] = 12 ⋱
⋮
[𝑙1𝑛 𝑙 𝑛𝑛 ]
Pour ce faire, on effectue la multiplication [𝐿 ][𝐿 ]𝑇 et on compare aux
coefficients de [𝐴]. Par identification, on obtient les équations permettant de
déduire les coefficients 𝑙𝑖𝑗 . On a :
𝑙 211 = 𝑎11 ⟹ 𝑙 11 = √𝑎11
𝑎1𝑗
𝑙11 . 𝑙 1𝑗 = 𝑎1𝑗 ⟹ 𝑙 1𝑗 = pour 𝑗 = 2, … , 𝑛
𝑙 11
𝑙 212 + 𝑙 222 = 𝑎22 2
⟹ 𝑙 22 = √𝑎22 − 𝑙12
𝑎2𝑗 − 𝑙 12 . 𝑙1𝑗
𝑙 12 . 𝑙 1𝑗 + 𝑙 22 . 𝑙 2𝑗 = 𝑎2𝑗 ⟹ 𝑙 2𝑗 = , 𝑗 = 3, … , 𝑛
𝑙 22
2
𝑙13 + 𝑙 223 + 𝑙 233 = 𝑎33 ⟹ 𝑙 33 = √𝑎33 − 𝑙 213 − 𝑙 223
𝑎3𝑗 − 𝑙13 . 𝑙 1𝑗 − 𝑙 23 . 𝑙 2𝑗
𝑙 13 . 𝑙1𝑗 + 𝑙 23 . 𝑙 2𝑗 + 𝑙 33 . 𝑙 3𝑗 = 𝑎3𝑗 ⟹ 𝑙 3𝑗 = , 𝑗 = 4, … , 𝑛
𝑙 33
Ainsi de suite, tous les coefficients de [𝐿 ] peuvent être déterminés. Nous
avons :
[𝐴]{𝑥 } = {𝑏}
[𝐿 ][𝐿 ]𝑇 {𝑥} = {𝑏}
Soit
[𝐿 ]𝑇 {𝑥} = {𝑦}
{
[𝐿 ]{𝑦} = {𝑏}
𝑙 11 0 𝑦1 𝑏1
𝑙 𝑦2 𝑏2
𝑙 12 21 ⋮
⋱ = ⋮
⋮
[𝑙1𝑛 𝑙 𝑛𝑛 ] {𝑦𝑛 } {𝑏 }
𝑛
La première colonne donne :
𝑏
𝑙 11 𝑦1 = 𝑏1 ⟹ 𝑦1 = 1
𝑙11
La deuxième colonne donne :
𝑙12 𝑦1 + 𝑙 22 𝑦2 = 𝑏2
En remplaçant 𝑦1 par son expression, on trouve 𝑦2 et ainsi de suite jusqu’à
𝑦𝑛 . Une fois le vecteur {𝑦} déterminé, on utilise :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 55
Dynamique des Structures
𝑇 𝑥1 𝑦1
𝑙 11 0
𝑙 21 𝑥2 𝑦2
𝑙 12 ⋮ = ⋮
⋱
⋮
[𝑙 1𝑛 𝑙 𝑛𝑛 ] {𝑥 𝑛 } {𝑦𝑛 }
On détermine 𝑥 𝑛 , 𝑝𝑢𝑖𝑠 𝑥 𝑛−1 , ainsi 𝑥 2 , 𝑥 1 par remontée.
Séquences de la factorisation de Cholesky :
Etape 1 : Factorisation de Cholesky [𝐴] = [𝐿 ][𝐿 ]𝑇
Etape 2 : [𝐿]{𝑦} = {𝑏}
Etape 3 : Substitution en remontée [𝐿]𝑇 {𝑥 } = {𝑦}.
Exemple d’application : vibrations axiales d’une barre encastrée a une
extrémité et libre à l’autre
Considérons une barre uniforme de longueur 𝐿 = 1 𝑚, de section 𝐴 =
3𝑜 × 10−6 𝑚2 , de module de Young 𝐸 = 2 × 1011 𝑁⁄𝑚2 𝑒𝑡 𝜌 = 7800 𝑘𝑔⁄𝑚3 .
Figure 13 : Barre encastrée
En utilisant les matrices de raideur et de masse globales suivantes :
0.360 0.180 0 0.052 0.013 0
[𝐾] = 108 [−0.180 0.360 −0.180] et [𝑀] = [ 0.013 0.052 0.013]
0 −0.180 0.180 0 0.013 0.026
Déterminons les valeurs propres et vecteurs propres de ce système en
utilisant la méthode de l’itération simultanée :
Soit :
1 0 0
{̅̅̅
𝑋1 }3× 3 = [0 1 0]
0 0 1
Alors :
1 0 0 0.052 0.013 0 1 0 0
𝑋1 }𝑇3×3 [𝑀]{̅̅̅
[𝐷]3×3 = {̅̅̅ 𝑋1 }3×3 = [ 0 1 0] [ 0.013 0.052 0.013] [0 1 0]
0 0 1 0 0.013 0.026 0 0 1
0.052 0.013 0
𝑋1 }𝑇3×3 [𝑀]{̅̅̅
[𝐷]3×3 = {̅̅̅ 𝑋1 }3×3 = [ 0.013 0.052 0.013]
0 0.013 0.026
Il faut ensuite déterminer la factorisation de Cholesky : [𝐷]3×3 = [𝐿 ]𝑇3×3 [𝐿 ]3×3
0.228 0 0 0.228 0.057 0
On obtient : [𝐿 ]3×3 [ 0.057 0.220 0 ] et [𝐿]𝑇3×3 = [ 0 0.220 0.059]
0 0.059 0.150 0 0 0.150
On détermine ensuite { 𝑋1 }3×3 en résolvant l’équation {𝑋1 }3×3 = [𝐿 ]𝑇3×3 {̅̅̅ 𝑋1 }3×3
4.39 −1.13 0.444
On trouve : {𝑋1 }3×3 = [ 0 4.53 −1.78]
0 0 6.66
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 56
Dynamique des Structures
Il faut ensuite calculer [𝑅]3×3 avec l’équation : [𝑅]3×3 = [𝑀]3×3 { 𝑋1 }3×3
0.228 0 0
[𝑅] = [ 0.057 0.221 0 ]
0 0.0589 0.15
Il s’agit ensuite de déterminer {̅̅̅ X2 } en résolvant l’équation [K] {̅̅̅
X2 } = [R]
0.158 0.155 0.0834
On obtient : {̅̅̅ X2 } = 10−7 [0.190 0.311 0.167 ]
0.190 0.343 0.250
On calculer ensuite la matrice quotient de Rayleigh [λ] = {X1 }T [K] {X1 }
On répète les étapes 3 à 8 pour les itération suivantes. Après cinq itérations,
on obtient notre approximation des valeurs propres et des vecteurs propres,
qui sont donnes respectivement par les éléments diagonaux de la matrice de
Rayleigh et les colonnes de la matrice [𝜙]:
0.0049 0 0 1.50 −3.59 2.37
[ λ] = [ 0 0.0015 0 ] et [𝜙] = [2.59 0.008 −4.12]
0 0 0.0008 2.99 3.57 4.76
𝜆1 = 0.0049 𝑟𝑎𝑑⁄𝑠 ; 𝜆 2 = 0.0015 𝑟𝑎𝑑⁄𝑠 ; 𝜆 3 = 0.0008 𝑟𝑎𝑑⁄𝑠
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 57
Dynamique des Structures
Chapitre 4 : Méthodes Numériques Appliquées à
la Dynamique des Structures
Méthode des Eléments Finis en Dynamique des Poutres
Introduction
Comme nous le savons la plupart des phénomènes physiques sont gouvernés
par des équations aux dérivées partielles (EDP). La résolution de ces équations
sous forme analytique est souvent très ardue voire impossible. C’est ainsi que
l’on a développé des méthodes numériques de résolution de ces équations.
En physique, l’objectif des méthodes numériques est de fournir des solutions
approchées du comportement réel d’un phénomène physique qui est
généralement décrit par des équations différentielles ou intégrales. Ce sont
des méthodes d’approximation qui ont un même objectif : remplacer un
problème mathématique défini sur un milieu continu (équations
différentielles ou intégrales) par un problème mathématique discret
(équations matricielles), problème de dimension finie que l’on sait résoudre.
La méthode des éléments finis (MEF) consiste à remplacer le système continu
par un système discret. La structure réelle est alors subdivisée en certains
nombre de parties appelés éléments dont l’assemblage permet de reconstituer
la géométrie initiale. Chacun des éléments est relié à ses voisins par des
nœuds dont les degrés de liberté constituent les inconnus du problème. On
vient ainsi de définir le maillage d’une structure. Le maillage, c’est l’ensemble
des éléments et nœuds utilisés pour la discrétisation. Cette discrétisation
peut se faire, selon la nature du problème, par des éléments
monodimensionnels (barre, poutre), par des éléments bidimensionnels
(plaques) ou par des éléments tridimensionnels.
Sur chacun des éléments finis, il est possible de linéariser l'EDP, c'est-à-dire
de remplacer l'équation aux dérivées partielles par un système d'équations
linéaires, par approximation. Ce système d'équations linéaires peut se décrire
par une matrice ; il y a donc une matrice par élément fini. Cependant, les
conditions aux frontières sont définies sur les frontières du système global et
pas sur les frontières de chaque élément fini ; il est donc impossible de
résoudre indépendamment chaque système. Les matrices sont donc réunies
au sein d'une matrice globale. Le système d'équations linéaires global est
résolu par ordinateur (des systèmes simples peuvent être résolus à la main et
constituent en général des exercices d'apprentissage).
L'EDP est résolue aux nœuds du maillage, c'est-à-dire que la solution est
calculée en des points donnés (résolution discrète) et non en chaque point du
domaine Ω. Cela nécessite de pouvoir interpoler, c'est-à-dire déterminer les
valeurs en tout point à partir des valeurs connues en certains points. On
utilise en général des fonctions polynomiales. Un élément fini est la donnée
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 58
Dynamique des Structures
d'une cellule élémentaire et de fonctions de base de l'espace d'approximation
dont le support est l'élément, et définies de manière à être interpolantes.
I. Présentation de la MEF :
1. Procédure de mise en œuvre de la MEF :
Les étapes suivantes peuvent être suivies pour un calcul par éléments finis.
Etape 1 : Discrétisation du système réel et sélection des types d’éléments
La première étape consiste à diviser une région de solution en éléments finis.
Le maillage d'éléments finis est généralement généré par un programme de
préprocesseur. La description du maillage se compose de plusieurs tableaux
dont les principaux sont les coordonnées nodales et les connectivités des
éléments.
Dans l’exemple de la figure 1, le solide continu (figure (a)) est subdivisé en un
nombre fini d’éléments de formes géométriques spécifique interconnectés
aux nœuds (figure (b)).
Le modèle discrétisé offre une approximation dans la géométrie du solide
continu original. Une différence notable est que le solide continu avec son
contour est maintenant représenté par des cordes de cotés droits. On associe
un système de coordonnées fixes au modèle discrétisé.
Figure 14
La structure continue est subdivisée en un certain nombre d'éléments par des
lignes ou des surfaces imaginaires.
- Les éléments interconnectés peuvent avoir des tailles et des formes
différentes.
- Établir le maillage EF avec des coordonnées définies, des numéros
d'éléments et des numéros de nœuds.
- Le modèle EF discrétisé doit être situé avec un système de coordonnées
- Les éléments et les nœuds du modèle EF discrétisé doivent être identifiés
par des numéros d'éléments et des numéros de nœuds.
- Les nœuds sont identifiés par les numéros de nœuds attribués et leurs
coordonnées correspondantes.
Etape 2 : Identifier les principales quantités inconnues pour l'analyse
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 59
Dynamique des Structures
Il s’agit de déterminer la première et principale quantité inconnue à obtenir
par la méthode des éléments finis.
Exemple : En analyse des contraintes : déplacement {𝑢} aux nœuds.
Dans l'analyse des contraintes, les principales inconnues sont les
déplacements nodaux, mais les quantités inconnues secondaires
comprennent : les déformations dans les éléments qui peuvent être obtenues
par les relations déformations-déplacements et les contraintes inconnues dans
les éléments qui peuvent être obtenues par les relations contraintes-
déformations (loi de Hooke).
Etape 3 : Interpolation des fonctions de déplacements ou fonctions
d’interpolation ou encore fonctions d’approximation
La fonction de déplacement est le point de départ de l'analyse mathématique.
Elle représente la variation du déplacement au sein de l'élément.
La fonction de déplacement peut être approximée sous la forme d'une
fonction linéaire ou d'une fonction d'ordre supérieur.
Une façon pratique de l'exprimer est d'utiliser des expressions polynomiales,
le degré du polynôme dépend du nombre de nœuds attribués à l'élément.
La forme ou la géométrie de l'élément peut également être approximée.
Les fonctions d’interpolation relient les quantités primaires dans les éléments
et celles associées aux nœuds des mêmes éléments.
Les fonctions d’interpolation vérifient la relation :
{u(𝑥, 𝑦, 𝑧)} = [𝑁(𝑥, 𝑦, 𝑧)]{𝑞}
où {u(𝑥, 𝑦, 𝑧)} représente les quantités primaires inconnues dans l’élément ;
[𝑁(𝑥, 𝑦, 𝑧)]: Est la fonction d’interpolation sous forme matricielle ;
{𝑞}: Représente les mêmes quantités primaires inconnues aux nœuds associés
au même élément.
Les fonctions d’interpolation auront les formes suivantes pour les différentes
formes d’éléments :
[𝑁(𝑥 )]𝑇 Fonctions 1𝐷 pour les éléments de barres et de poutres.
[𝑁(𝑥, 𝑦)] Fonctions 2𝐷 pour les éléments de plaques triangulaires et
quadrilatérales
[𝑁(𝑥, 𝑦, 𝑧)] Fonctions 3𝐷 pour les éléments tétraédriques et hexaédriques
Etape 4 : Définir les relations déformations/déplacements et
contraintes/déformations
Dans le cas de la méthode des éléments finis appliquée à la mécanique des
solides, il s’agit des lois de l’Elasticité.
Etape 5 : déterminer les matrices de rigidité élémentaires
Une fois la structure discrétisée avec les formes d'éléments souhaitées, les
matrices de rigidité des éléments individuels sont formulées. Ces sont
exprimées dans un premier temps en repère local (lié à l’élément) puis dans
le repère global.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 60
Dynamique des Structures
Etape 6 : Assembler la matrice de rigidité globale et l’équation globale
Une fois les matrices de rigidité des éléments formées dans le repère global,
elles sont assemblées pour former la matrice de rigidité globale.
- L'assemblage se fait via les nœuds communs aux éléments adjacents.
- La matrice de rigidité globale est symétrique et en bandes.
Etape 7 : Introduire des conditions aux limites
Les conditions de contrainte aux limites doivent être imposées dans la matrice
de rigidité. Il existe différentes techniques pour satisfaire les conditions aux
limites. La taille de la matrice de rigidité global peut être réduite ou condensée
dans sa forme finale. Pour faciliter l'aspect programmation informatique et
intégrer élégamment les conditions aux limites, la taille de la matrice globale
est conservée à la même valeur.
Etape 8 : Résoudre la matrice globale (pour obtenir les déplacements) ou
résolution de l’équation globale
Les déplacements nodaux inconnus sont calculés par la multiplication du
vecteur de force avec l'inverse de la matrice de rigidité.
Le système global peut être linéaire ou non. Il définit soit un problème
d'équilibre qui concerne un cas stationnaire ou statique soit un problème de
valeurs critiques où il faut déterminer les valeurs et vecteurs propres du
système qui correspondent généralement aux fréquences et modes propres
d'un système physique. Un problème de propagation qui concerne le cas
transitoire (non stationnaire) dans lequel il faut déterminer les variations dans
le temps des variables physiques et la propagation d'une valeur initiale.
Etape 9 : détermination des déformations et contraintes (à partir des
déplacements)
Les déplacements nodaux sont utilisés pour le calcul des contraintes ou des
déformations. Cela peut être fait pour tous les éléments de la structure
continue ou peut être limité à certains éléments prédéterminés.
2. Avantages et inconvénients de la méthode des éléments finis
a. Avantages de la méthode des éléments finis :
- La modélisation de géométries complexes et de formes irrégulières est
plus facile car une variété d'éléments finis est disponible pour la
discrétisation du domaine.
- Les conditions aux limites peuvent être facilement incorporées dans la
méthode des éléments finis.
- Différents types de propriétés de matériaux peuvent être facilement pris
en compte dans la modélisation d'un élément à l'autre ou même au sein
d'un élément.
- Des éléments d'ordre supérieur peuvent être mis en œuvre.
- La méthode des éléments finis est simple, compacte et axée sur les
résultats et donc très populaire au sein de la communauté des ingénieurs.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 61
Dynamique des Structures
- La disponibilité d'un grand nombre de progiciels informatiques et de
littérature fait de la méthode des éléments finis une méthode numérique
polyvalente et puissante.
b. Inconvénients de la méthode des éléments finis
- Une grande quantité de données est requise en entrée pour le maillage
utilisé en termes de connectivité nodale et d'autres paramètres en fonction
du problème.
- Elle nécessite un ordinateur numérique et est assez étendue
- Elle nécessite un temps d'exécution plus long que la méthode des éléments
finis.
- Le résultat de sortie variera considérablement.
c. Limitations de la MEF
i. Il faut faire preuve de jugement technique pour interpréter les résultats.
ii. Il faut une mémoire informatique importante et du temps de calcul pour
obtenir les résultats escomptés.
iii. Il existe certaines catégories de problèmes pour lesquels d'autres
méthodes sont plus efficaces, par exemple, les problèmes de fluides
ayant des limites à l'infini sont mieux traités par la méthode des
éléments limites.
iv. Pour certains problèmes, il peut y avoir une quantité considérable de
données d'entrée. Des erreurs peuvent se glisser dans leur préparation
et les résultats ainsi obtenus peuvent également sembler acceptables,
ce qui indique un état de fait trompeur. Il est toujours souhaitable de
faire une vérification visuelle des données d'entrée.
v. Dans la méthode des éléments finis, de nombreux problèmes
conduisent à des erreurs d'arrondi. L'ordinateur fonctionne avec un
nombre limité de chiffres et la résolution du problème avec un nombre
restreint de chiffres peut ne pas donner le degré de précision souhaité
ou peut donner des résultats totalement erronés dans certains cas. Pour
de nombreux problèmes, l'augmentation du nombre de chiffres à des
fins de calcul améliore la précision.
II. Mise en œuvre de la MEF :
1. Discrétisation :
La méthode des éléments finis repose sur l’idée de remplacer de physique (et
de mathématique) continu par un problème de physique (et de mathématique)
discret.
Cette discrétisation consiste à subdiviser donc le domaine réelle, de
géométrie souvent très complexe, en un grand nombre petits domaines de
géométrie bien choisie, appelé éléments ou mailles, facilitant l’étude du
problème. Les points de jonction entre les éléments sont appelés nœuds.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 62
Dynamique des Structures
Figure 15 : Discrétisation d’un domaine réel
Types d’Eléments finis :
La discrétisation peut se faire, selon la nature du problème, par des éléments
monodimensionnels (barre, poutre), par des éléments bidimensionnels
(plaques) ou par des éléments tridimensionnels.
Les éléments unidimensionnels 1D : barre, poutre rectiligne ou courbe.
Un solide peut être considéré comme étant 1D lorsqu’une des dimensions est
très large par rapport aux deux autres.
Les éléments bidimensionnels 2D : élasticité plane (contraintes et
déformations planes), plaque en flexion, coque courbe.
Dans le cas du 2D, deux grandeurs sont très larges par rapport à la troisième.
On a par exemple les plaques, les coques minces, les membranes.
Les éléments tridimensionnels 3D : Ce sont les solides ou les trois dimensions
sont dans le même ordre de grandeur.
On a les éléments de volume ou coques épaisses.
Eléments classiques utilisés pour la discrétisation 1D
Barre à 2 nœuds Barre à 3 nœuds Poutre à 2 nœuds
Eléments classiques utilisés pour la discrétisation 2D
Triangle de Triangle de quadrilatère de quadrilatère de
degré 1 degré 2 degré 1 degré 2
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 63
Dynamique des Structures
3 nœuds 6 nœuds 4 nœuds 8 nœuds
Eléments classique utilisés pour la discrétisation en 3D
Tétraèdres de Tétraèdres de hexaèdres de hexaèdres de
degré 1 degré 2 degré 1 degré 2
4 nœuds 10 nœuds 8 nœuds 20 nœuds
Eléments classique utilisés pour la discrétisation en 3D
Pentaèdre (prisme) Pentaèdre (prisme) Pyramide de Pyramide de
de degré 1 de degré 2 degré 1 degré 1
6 nœuds 15 nœuds 4 nœuds 10 nœuds
2. Identification des principales quantités inconnues pour l'analyse
Il s’agit de déterminer la première et principale quantité inconn ue à obtenir
par la méthode des éléments finis.
En analyse de contrainte, les principales inconnues sont les déplacements aux
nœuds {𝑢} , {𝑢} 𝑇 = {𝑢 𝑥 𝑢 𝑦 𝑢 𝑧 } avec 𝑢 𝑥 𝑢 𝑦 𝑢 𝑧 les déplacements respectifs
suivant les axes 𝑥, 𝑦, et 𝑧.
Les quantités inconnues secondaires comprennent : les déformations dans les
éléments et peuvent être obtenues par les « relations déformation-
déplacement » et les contraintes inconnues dans les éléments par les relations
contrainte-déformation (loi de Hooke).
3. Interpolation nodale ou des fonctions de déplacements :
La méthode consiste le plus souvent à approximer la fonction de déplacement
au moyen d’une fonction polynomiale formée d’un nombre fini de paramètres
et à l’exprimer en fonction des déplacements nodaux (les déplacements
associés aux degrés de liberté). On aura donc :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 64
Dynamique des Structures
𝑛
𝑢(𝑥, 𝑦, 𝑧) = ∑ 𝑁𝑖 (𝑥, 𝑦, 𝑧) ∙ 𝑢 𝑖 = 𝑁1 (𝑥, 𝑦, 𝑧) ∙ 𝑢 1 + 𝑁2 (𝑥, 𝑦, 𝑧) ∙ 𝑢 2 + ⋯ + 𝑁𝑛 (𝑥, 𝑦, 𝑧) ∙ 𝑢 𝑛
𝑖=1
où :
𝑢 𝑖 désigne les déplacements aux degrés de liberté, i.e. aux nœuds ;
𝑁𝑖 (𝑥, 𝑦, 𝑧) désigne les fonctions de forme de l’élément permettant d’obtenir les
déplacements de celui-ci en n’importe quel point de sa géométrie et toujours
à partir des déplacements nodaux.
En regroupant toutes les fonctions de forme de l’élément au sein de la matrice
[𝑁] et les déplacements nodaux dans le vecteur {𝑞𝑒 }, le champ de déplacement
s’exprime alors :
{𝑢} = [𝑁]{𝑞𝑒 }
Le degré du polynôme retenu pour l’élément est donc intimement lié aux
nombres de nœuds et de degrés de liberté par nœud. Par exemple, un élément
unidimensionnel à trois nœuds aura une fonction d’approximation
parabolique.
Généralement, les fonctions d’interpolation retenues sont linéaires,
quadratiques et plus rarement cubiques ce qui ne se révèle pas toujours
conforme aux théories visées telles que celles des poutres, des plaques ou des
coques. En effet, les ligne et surface élastique découlant de ces théories sont
plutôt du troisième voire du quatrième degré. Une solution revient à
discrétiser le modèle de telle manière à reconstituer le champ réel de
déplacement à partir de fonctions de degré inférieur. La qualité et l’efficacité
d’une telle reconstitution dépendront donc d’une part du type d’élément
choisi, en d’autres termes de ses capacités, et d’autre part de la densité du
maillage. Dans le cas des éléments courants et toujours dans une hypothèse
de linéarité, ce choix aura pour conséquence de ne pas garantir la continuité
du champ de contrainte d’où la recherche d’une discrétisation adaptée
permettant sa bonne reconstitution. De ce fait, la taille des éléments utilisés
sera généralement petite.
Les fonctions de déplacement sont intimement liées aux choix des éléments
du maillage.
4. Relations déformations/déplacements et contraintes/déformations
En élasticité le champ de déformations et de contraintes se déduisent
directement du champ des déplacement. On cherchera d’abord les
déplacements des nœuds à partir desquels on détermine les autres champs.
Les hypothèses des petits déplacements et des petites déformations restent
applicables.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 65
Dynamique des Structures
Déplacements Déformations Contraintes
𝑢1 𝜀𝑥𝑥 𝜀𝑥𝑥 𝜎𝑥𝑥
𝜀𝑦𝑦 𝜀 𝑦𝑦 𝜎𝑦𝑦
{𝑢} = 𝑢 2 𝜀 𝜀𝑧𝑧 𝜎
{𝜀 } = 𝛾𝑧𝑧 = 2𝜀 {𝜎} = 𝜏 𝑧𝑧
𝑥𝑦 𝑥𝑦 𝑥𝑦
𝛾𝑥𝑧 𝜏𝑥𝑧
(𝑢 3 ) 2𝜀𝑥𝑧
( 𝛾𝑦𝑧 ) (2𝜀𝑦𝑧 ) 𝜏
( 𝑦𝑧 )
Les champ de déplacements {𝑢} , de déformations {𝜀 } et de contraintes {𝜎} sont
liés par les relations :
{𝜀 } = [𝜕]{𝑢} puis {𝜎} = [𝐻 ] ∙ {𝜀 }
𝜕
0 0
𝜕𝑥
𝜕
0 0
𝜀𝑥𝑥 𝜕𝑦 𝑢1
𝜀𝑦𝑦 𝜕
0 0
𝜀 𝜕𝑧
{𝜀 } = 𝛾𝑧𝑧 = 𝑢2 ;
𝑥𝑦 𝜕
𝛾𝑥𝑧 0 0
𝜕𝑦
(𝛾𝑦𝑧 ) (𝑢 3 )
𝜕 𝜕
0
𝜕𝑧 𝜕𝑦
𝜕 𝜕
( 𝜕𝑧 0
𝜕𝑥 )
E est le module d’élasticité de Young ;
𝑣 est le coefficient de Poisson.
Energie de déformation :
L’énergie de déformation 𝑊 est déterminée à partir des tenseurs de
déformations et de contraintes.
1
𝑊 = ∫ {𝜀 }𝑇 {𝜎}𝑑𝑉
2 𝑉
5. Détermination des matrices de rigidité élémentaires
On se base sur la formulation de l’énergie potentielle totale minimale pour
déterminer les matrices de rigidité élémentaires. Son expression est (déjà
exprimée) :
1
𝑊 = ∫ {𝜀 }𝑇 {𝜎} 𝑑𝑉
2 𝑉
Comme exprimé dans l’étape précédente, le champ des déformations {𝜀 } :
{𝜀 } = [𝜕]{𝑢}
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 66
Dynamique des Structures
Suite à la discrétisation, le champ de déplacement {𝑢} devient {𝑢} = [𝑁]{𝑞𝑒 },
[𝑁] étant la matrice regroupant les fonctions de forme des éléments utilisées
pour la discrétisation et {𝑞𝑒 } représente le vecteur des déplacements nodaux.
On a alors :
{𝜀 } = [𝜕]{𝑢} = [𝜕][𝑁]{𝑞𝑒 } = [𝐵 ]{𝑞𝑒 }
Le champ des contraintes {𝜎} :
{𝜎} = [𝐻] ∙ {𝜀 }
Alors l’énergie de déformation associée à l’élément 𝑒 peut être déduite par :
1 1 1
𝑊𝑒 = ∫ [𝜀 ]𝑇 [𝜎]𝑑𝑉𝑒 = ∫ [𝑞𝑒 ]𝑇 [𝐵 ]𝑇 [𝐻][𝐵]{𝑞𝑒 }𝑑𝑉𝑒 = {𝑞𝑒 }𝑇 [𝑘𝑒 ]{𝑞𝑒 }
2 𝑉𝑒 2 𝑉𝑒 2
On définit la matrice de rigidité élémentaire par l’intégrale volumique sur le
volume 𝑉𝑒 :
[𝑘𝑒 ] = ∫ [𝐵 ]𝑇 [𝐻][𝐵]𝑑𝑉𝑒
𝑉𝑒
a. Matrice de rigidité élémentaire pour l’élément poutre en repère local
:
Soit un élément de poutre plane de longueur 𝐿, de section 𝑆 et d’inertie 𝐼
constantes et dont le matériau a un module d’élasticité longitudinal 𝐸, voir
figure 16.
Figure 16 : Elément de poutre plane
Cet élément à deux nœuds destiné au calcul des réseaux de poutres chargés
dans leur plan fait appel à la théorie des poutres qui permet de ramener le
problème tridimensionnel à un problème unidimensionnel en condensant ses
caractéristiques au niveau de sa fibre moyenne. Chacun de ses nœuds possède
𝑑𝑣𝑖
trois degrés de liberté 𝑢 𝑖 , 𝑣𝑖 , qui permettent de reconstituer les fonctions
𝑑𝑡
de déplacements axial.
Champ de déplacement axial 𝑢(𝑥 ) ∶
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 67
Dynamique des Structures
En l’absence de flexion, l’élément poutre travaille uniquement en traction-
compression, le matériau utilise étant supposé linéaire, la déformation
constante sur la hauteur de la section. Le déplacement axial 𝑢(𝑥 ) est alors une
fonction linéaire de la forme 𝑢(𝑥 ) = 𝑎0 + 𝑎1 𝑥 . La déformation associée est
homogène à la section et vaut :
𝑑𝑢
𝜀𝑥𝑥 = = 𝑐𝑠𝑡𝑒 ⟹ 𝑢(𝑥 ) = 𝑎0 + 𝑎1 𝑥
𝑑𝑥
De plus, si aux nœuds 𝑖 et 𝑗, on a : 𝑢 (0) = 𝑢 𝑖 et 𝑢(𝐿 ) = 𝑢𝑗 , on peut determiner
les constantes 𝑎0 et 𝑎1 . On aura :
𝑥 𝑥 𝑥 𝑥 𝑢𝑖
𝑢(𝑥 ) = (1 − ) 𝑢 𝑖 + 𝑢𝑗 = {1 − }{ }
𝐿 𝐿 𝐿 𝐿 𝑢𝑗
𝑥 𝑥
𝑢(𝑥 ) = (1 − ) 𝑢 𝑖 + 𝑢𝑗 = 𝑢(𝑥 ) = 𝑁1 𝑢 𝑖 + 𝑁2 𝑢𝑗 = [𝑁]{𝑞𝑒 } .
𝐿 𝐿
Champ de déplacement transversal 𝑣(𝑥 ) ∶
En flexion, chacun des nœuds d’un élément possède deux degrés de liberté,
un en translation dans la direction transversale 𝑣 et un associé à la rotation 𝛽
autour de l’axe 𝑧 perpendiculaire au plan 𝑥𝑦. Quatre conditions aux limites
𝑣𝑖 , 𝛽𝑖 , 𝑣𝑗 , 𝛽𝑗 peuvent donc être utilisées pour définir sa fonction
d’approximation ce qui explique l’expression de 𝑣 (𝑥 ) sous la forme d’un
polynôme de degré 3.
Le champ de déplacements suppose pour la poutre devrait être tel qu’il prend
les valeurs de déflection et de pente à chaque extrémité comme données par
les valeurs nodales 𝑣𝑖 , 𝛽𝑖 , 𝑣𝑗 , 𝛽𝑗 . Soit 𝑣(𝑥 ) donné par :
𝑐0
𝑐
𝑣 (𝑥 ) = 𝑐0 + 𝑐1 𝑥 + 𝑐2 𝑥 2 + 𝑐3 𝑥 3 = { 1 𝑥 𝑥 2 𝑥 3 } { 𝑐1 }
2
𝑐3
Alors, on a :
𝑑𝑣⁄𝑑𝑥 = 𝑐1 + 2𝑐2 𝑥 + 3𝑐3 𝑥 2
A partir des conditions aux limites, 𝑥 = 0 et 𝐿, le vecteur des déplacements
nodaux s’écrit :
𝑣𝑖 = 𝑐0 , 𝛽𝑖 = 𝑐1 , 𝑣𝑗 = 𝑐0 + 𝑐1 𝐿 + 𝑐2 𝐿2 + 𝑐3 𝐿3 , 𝛽𝑗 = 𝑐1 + 2𝑐2 𝐿 + 3𝑐3 𝐿2
𝑣𝑖 1 0 0 0 𝑐0
𝜃𝑖 0 1 0 0 𝑐
{𝑣𝑗 } = [ 2 3 ] { 𝑐1 }
1 𝐿 𝐿 𝐿 2
𝜃𝑗 0 1 2𝐿 3𝐿2 𝑐3
1 0 0 0 𝑣𝑖 𝑐0
0 1 0 0 𝛽 𝑖 𝑐
Posons [𝑅] = [ 2 3 ], alors on a : {𝑣𝑗 } = [𝑅] { 𝑐1 }
1 𝐿 𝐿 𝐿 2
0 1 2𝐿 3𝐿2 𝛽𝑗 𝑐3
𝑐0 𝑣𝑖 1 0 0 0 𝑣𝑖
𝑐 𝛽𝑖 0 1 0 0 𝛽𝑖
Soit {𝑐1 } = [𝑅]−1 { 𝑣𝑗 } = [−3⁄𝐿2 −2⁄𝐿 3⁄𝐿2 −1⁄𝐿 ] { 𝑣𝑗 }
2
𝑐3 𝛽𝑗 2⁄𝐿3 1⁄𝐿2 −2⁄𝐿3 1⁄𝐿2 𝛽𝑗
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 68
Dynamique des Structures
Finalement, on a :
𝑣𝑖
2 3 𝛽𝑖
3𝑥 2𝑥 2
2𝑥 𝑥 3
3𝑥2
2𝑥3
𝑥2
𝑥 3
𝑣(𝑥 ) = {1 − 2
− 3 𝑥− − 2 − − + 2 } {𝑣𝑗 }
𝐿 𝐿 𝐿 𝐿 𝐿2 𝐿3 𝐿 𝐿
𝛽𝑗
On pose :
3𝑥 2 2𝑥 3 2𝑥 2 𝑥 3 3𝑥 2 2𝑥 3 𝑥2 𝑥3
𝑁1 = 1 − − ; 𝑁2 = 𝑥 − − 2; 𝑁3 = 2 − 3 ; 𝑁4 = − + 2
𝐿2 𝐿3 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿
𝑣𝑖
𝛽
𝑣(𝑥 ) = {𝑁1 𝑁2 𝑁3 𝑁4 } {𝑣𝑖 }
𝑗
𝛽𝑗
On vient ainsi d’exprimer le champ de déplacement en fonction des
déplacements aux nœuds et des fonctions 𝑁1 ; 𝑁2 ; 𝑁3 ; 𝑁4 qui sont appelés
fonctions de forme.
Champ des déformations :
Le champ de déplacement complet s’écrit en tout point de la poutre par :
𝑢𝑖
𝑥 𝑥
1− 0 0 0 0 𝑣𝑖
𝑢(𝑥 ) 𝐿 𝐿 𝛽𝑖
{ }= 𝑢𝑗
3𝑥 2 2𝑥 3 2𝑥 2
𝑥 3
3𝑥 2
2𝑥 3
𝑥 2
𝑥 3
𝑣(𝑥 ) 0 1− 2 − 3 𝑥− − 2 0 − 3 − + 2 𝑣𝑗
{ 𝐿 𝐿 𝐿 𝐿 𝐿2 𝐿 𝐿 𝐿 } {𝛽 }
𝑗
𝑢(𝑥 )
{ } = [𝑁]{𝑞𝑒 }
𝑣(𝑥 )
La déformation de la poutre résultant de concomitance de la flexion et des
variations de longueur, on a :
𝑑𝑢 𝑑 𝑑𝑣 𝑑𝑢 𝑑2𝑣
𝜀𝑥𝑥 = − (−𝑦 ) = −𝑦
𝑑𝑥 𝑑𝑥 𝑑𝑥 𝑑𝑥 𝑑𝑥 2
Le champ de déformation devient alors :
𝑣𝑖
2 2
𝑑 𝑁 𝑑 𝑁2 𝑑 2 𝑁3 𝑑 2 𝑁4 𝜃𝑖
𝜀𝑥𝑥 = −𝑧 [ 21 ] { 𝑣𝑗 }
𝑑𝑥 𝑑𝑥 2 𝑑𝑥 2 𝑑𝑥 2
𝜃𝑗
𝑢𝑖
𝑣𝑖
1 −6 12𝑥 4 6𝑥 1 −6 12𝑥 2 6𝑥 𝛽𝑖
𝜀𝑥 = [− ( 2 + 3 ) ∙ 𝑦 ( − 2) ∙ 𝑦 ( 2 + 3 ) ∙ 𝑦 ( − 2 ) ∙ 𝑦 ] 𝑢𝑗
𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿
𝑣𝑗
{𝛽𝑗 }
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 69
Dynamique des Structures
𝑢𝑖
𝑣𝑖
𝛽𝑖
En posant {𝑞𝑒 } = 𝑢 , on a :
𝑗
𝑣𝑗
{𝛽𝑗 }
1 −6 12𝑥 4 6𝑥 1 −6 12𝑥 2 6𝑥
𝜀𝑥𝑥 = [− ( 2 + 3 ) ∙ 𝑦 ( − 2) ∙ 𝑦 ( 2 + 3 ) ∙ 𝑦 ( − 2 ) ∙ 𝑦] {𝑞𝑒 }
𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿
Soit la matrice [𝐵] définie par :
1 −6 12𝑥 4 6𝑥 1 −6 12𝑥 2 6𝑥
[𝐵 ] = [− ( 2 + 3 ) ∙ 𝑦 ( − 2) ∙ 𝑦 ( 2 + 3 ) ∙ 𝑦 ( − 2 ) ∙ 𝑦]
𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿 𝐿
D’où : 𝜀𝑥𝑥 = [𝐵 ]{𝑞𝑒 }
L’expression de l’énergie de déformation élémentaire est :
1 1 1
𝑊𝑒 = ∫ [𝜀 ]𝑇 [𝜎]𝑑𝑉𝑒 = ∫ [𝑞𝑒 ]𝑇 [𝐵]𝑇 [𝐻][𝐵 ]{𝑞𝑒 }𝑑𝑉𝑒 = [𝑞𝑒 ]𝑇 [𝑘𝑒 ]{𝑞𝑒 }
2 𝑉𝑒 2 𝑉𝑒 2
Sachant que la matrice [𝐻], les vecteurs contraintes et déformation se
réduisent respectivement aux seuls termes 𝐸, 𝜎𝑥𝑥 , 𝜀𝑥𝑥 , 𝑊𝑒 s’écrit :
1 1 1
𝑊𝑒 = ∫ 𝜎𝑥𝑥 ∙ 𝜀𝑥𝑥 𝑑𝑉𝑒 = ∫ {𝑞𝑒 }𝑇 [𝐵]𝑇 [𝐻][𝐵 ]{𝑞𝑒 }𝑑𝑉𝑒 = {𝑞𝑒 }𝑇 (𝐸 ∫ [𝐵 ]𝑇 [𝐵 ]𝑑𝑉𝑒 ) { 𝑞𝑒 }
2 𝑉𝑒 2 𝑉𝑒 2 𝑉𝑒
La matrice de rigidité de l’élément 𝑒 est alors définie par :
𝐿
[𝑘𝑒 ] = ∫ [𝐵 ]𝑇 [𝐵 ]𝐸𝑑𝑣 = 𝐸 ∫ (∫ 𝑑𝑆𝑒 ) ∙ [𝐵 ]𝑇 ∙ [𝐵 ] ∙ 𝑑𝑥
𝑉𝑒 0 𝑆𝑒
1
−
𝐿
−6 12𝑥
( + )∙𝑦
𝐿2 𝐿3
4 6𝑥
𝐿 ( − )∙𝑦 1 −6 12𝑥 4 6𝑥 1 −6 12𝑥 2 6𝑥
[ 𝑘𝑒 ] = 𝐸 ∫ (∫ 𝑑𝑆𝑒 ) 𝐿 𝐿2 [− ( − )∙𝑦 ( − ) ∙ 𝑦] 𝑑𝑥
1 ( + )∙𝑦 2 ( + )∙𝑦
0 𝑆𝑒 𝐿 2
𝐿 𝐿3 𝐿 𝐿 𝐿 2
𝐿 𝐿3 𝐿 𝐿2
𝐿
−6 12𝑥
( + )∙𝑦
𝐿2 𝐿3
2 6𝑥
{ ( 𝐿 − 𝐿2 ) ∙ 𝑦 }
Avec ∫𝑆 𝑑𝑆𝑒 = 𝑆 (section constante) ; ∫𝑆 𝑦 2 𝑑𝑆𝑒 = 𝐼 (inertie constante)
𝑒 𝑒
∫𝑆 𝑦𝑑𝑆𝑒 = 0 (moment statique nul par rapport à la fibre moyenne).
𝑒
On notera que les termes en 𝐸𝐼 sont bien liés aux degrés de liberté associes à
la flexion de l’élément poutre.
On obtient :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 70
Dynamique des Structures
Dans le cas où la poutre ne travaille qu’en flexion, les composantes de
l’effort axial s’annulent. Le vecteur {𝑞𝑒 } se réduit et le déplacement
devient :
0 0 0 0 𝑣𝑖
𝛽𝑖
{𝑣 ( 𝑥 ) } = { 3𝑥 2 2𝑥 3 2𝑥 2 𝑥 3 }
3𝑥 2 2𝑥 3
𝑥2 𝑥3 𝑣𝑗
1− 2 − 3 𝑥− − 2 − 3
− + 2
𝐿 𝐿 𝐿 𝐿 𝐿2 𝐿𝐿 𝐿 {𝛽 }
𝑗
En effectuant les calcul du champ de déformation et des contraintes, la
matrice de rigidité élémentaire devient :
12 6𝐿 −12 6𝐿
𝐸𝐼 6𝐿 4𝐿 2
−6𝐿 2𝐿2 ]
[ 𝑘𝑒 ] = 3 [
𝐿 −12 −6𝐿 12 6𝐿
2
6𝐿 2𝐿 −6𝐿 4𝐿2
b. Matrice de rigidité élémentaire pour l’élément poutre en repère global
:
L’élément poutre plane possédant deux nœuds à trois degrés de liberté, la
matrice de changement de base s’écrit :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 71
Dynamique des Structures
On notera l’ajout du 1 associe aux rotations autour des axes 𝑧 et 𝑍 qui sont
colinéaires dans le cas plan.
La matrice de rigidité globale obtenue est la suivante :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 72
Dynamique des Structures
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 73
Dynamique des Structures
6. Assemblage de la matrice rigidité globale :
Le calcul de la matrice de rigidité globale de la structure nécessite une
opération d’assemblage de toutes les matrices de rigidité des éléments
constituant le maillage. Pour cela, toutes les matrices élémentaires doivent
être exprimées dans la base globale. L’assemblage consiste à ranger chaque
terme dans la matrice globale, cela se fait par le biais des matrices booléennes
d’assemblage [𝐴𝑒 ]. Les matrices booléennes d’assemblage permettent de
placer les degrés de liberté associés à un élément donné dans le vecteur global
des déplacements {Q} de la structure, Elles sont telles que :
{Qe }=[Ae ]{Q}
Exemple : Soient deux éléments de connectivite 1, 2 𝑒𝑡 3.
Figure 17 : Assemblage de deux éléments
Des vecteurs de déplacement en repère global des deux éléments :
𝑈1
𝑉1
𝑈1 𝑈2
𝐵1
𝑉1 𝑉2
𝑈2
𝐵1 𝐵2
{𝑄1 } = , {𝑄2 } = et du vecteur global {𝑄} = 𝑉2 , on déduit pour
𝑈2 𝑈3
𝐵2
𝑉2 𝑉3
𝑈3
{𝐵2 } {𝐵3 }
𝑉3
{𝐵3 }
[A1 ] 𝑒𝑡 [A2 ]:
1 0 0 0 0 0 0 0 0 𝑈1
𝑈1 0 1 0 0 0 0 0 0 0 𝑉1
𝑉1 0 0 1 0 0 0 0 0 0 𝐵1
0 0 0 1 0 0 0 0 0 𝑈2
𝐵1
{ 𝑄1 } = = 0 0 0 0 1 0 0 0 0 𝑉2
𝑈2
0 0 0 0 0 1 0 0 0 𝐵2
𝑉2 0 0 0 0 0 0 0 0 0 𝑈3
{𝐵2 } 0 0 0 0 0 0 0 0 0 𝑉3
[0 0 0 0 0 0 0 0 0] { 𝐵 }
3
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 74
Dynamique des Structures
0 0 0 0 0 0 0 0 0 𝑈1
𝑈2 0 0 0 0 0 0 0 0 0 𝑉1
𝑉2 0 0 0 0 0 0 0 0 0 𝐵1
𝐵2 0 0 0 1 0 0 0 0 0 𝑈2
{𝑄2 } = = 0 0 0 0 1 0 0 0 0 𝑉2
𝑈3
0 0 0 0 0 1 0 0 0 𝐵2
𝑉3 0 0 0 0 0 0 1 0 0 𝑈3
{𝐵3 } 0 0 0 0 0 0 0 1 0 𝑉3
[0 0 0 0 0 0 0 0 1] { 𝐵3 }
En introduisant {Qe }=[Ae ]{Q} l’expression dans l’expression de l’énergie
potentielle, celle-ci devient :
1
𝐸 = 𝑊 − 𝑇 = ∑ ( {𝑄} 𝑇 [Ae ]𝑇 [Ke ][Ae ]{𝑄} − {𝑄} 𝑇 [Ae ]𝑇 ({𝐹𝑒𝑣 } + {𝐹𝑒𝑠 } + {𝑃𝑒 }))
2
𝑒
1 𝑇 𝑇 𝑣 𝑠
𝐸 = {𝑄}𝑇 (∑[Ae ] [Ke ][Ae ]) {𝑄} − {𝑄}𝑇 ∑[Ae ] ({𝐹𝑒 } + {𝐹𝑒 } + {𝑃𝑒 })
2
𝑒 𝑒
Sachant que les termes [Ae ] [Ke ][Ae ] 𝑒𝑡 [Ae ]𝑇 ({𝐹𝑒𝑣 } + {𝐹𝑒𝑠 } + {𝑃𝑒 })
𝑇
correspondent respectivement aux positionnements de [Ke ] dans [K] et des
{𝐹𝑒 } dans {𝐹 }. On obtient alors :
1
𝐸 = 𝑊 − 𝑇 = ∑ 𝑊𝑒 − ∑ 𝑇𝑒 = {𝑄} 𝑇 [K]{𝑄} − {𝑄} 𝑇 {𝐹 }
2
𝑒 𝑒
Avec :
∑[Ae ]𝑇 [Ke ][Ae ] = [𝐾] ∶ matrice de rigidité de la structure
𝑒
∑[ Ae ]𝑇 {𝑃𝑒 } = {𝑃} ∶ vecteur des charges nodales concentrées définies dans le repère global.
𝑒
Comme l’équilibre est caractérisé par une variation nulle d’énergie potentielle
par rapport au déplacement, le système final s’écrit :
𝑑𝐸
= {𝑑𝑄}𝑇 ∙ [K] ∙ {𝑄} − {𝑑𝑄}𝑇 ∙ {𝐹}
𝑑𝑄
𝑑𝐸
= 0 ⟺ {𝑑𝑄}𝑇 ∙ [K] ∙ {𝑄} − {𝑑𝑄}𝑇 ∙ {𝐹} ⟺ {𝑑𝑄}𝑇 ∙ ([𝐾] ∙ {𝑄} − {𝐹}) = 0
𝑑𝑄
[𝐾] ∙ {𝑄} − {𝐹 } = 0
De par sa construction, la matrice [𝐾] est symétrique. Elle est encore
singulière (det[𝐾] = 0) du fait de l’absence des conditions d’appui. Nous allons
maintenant introduire les conditions d’appui.
7. Introduction des conditions aux limites :
Pour fixer la structure au sol, les conditions d’appui sont introduites de
manière à assurer sa stabilité d’ensemble.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 75
Dynamique des Structures
Suivant la nature des degrés de liberté bloquées, ces supports peuvent
reconstituer des appuis ponctuels ou des encastrements, mais également
correspondre à des déplacements nuls ou imposés.
Regroupant les déplacements connus lies aux appuis au sein d’un vecteur
{𝑄𝐵 } et ceux connus dans {𝑄𝐿 }, le système [𝐾] ∙ {𝑄} − {𝐹 } = 0 peut-être
transforme de telle manière à séparer les degrés de liberté libres et bloques,
soit :
[𝐾 ] [𝐾𝐿𝐵 ] {𝑄𝐿 } {𝐹 }
[𝐾] ∙ {𝑄} − {𝐹 } = 0 ⟺ [ 𝐿𝐿 ]{ } = { 𝑒𝑥𝑡 }
[𝐾𝐵𝐿 ] [𝐾𝐵𝐵 ] {𝑄 } {𝑅}
𝐵
Avec :
{𝑅} ∶ vecteur des réactions d’appui,
{𝐹𝑒𝑥𝑡 } ∶ vecteur des charges associés aux degrés de liberté libres.
Les déplacements inconnues {𝑄𝐿 } peuvent donc être déterminés en résolvant
le système linéaire :
[𝐾𝐿𝐿 ] ∙ {𝑄𝐿 } = {𝐹𝑒𝑥𝑡 } − [𝐾𝐿𝐵 ] ∙ {𝑄𝐵 }
Pour être inversible, donc ne pas être singulière, [𝐾𝐿𝐿 ] ne devra pas comporter
des termes diagonaux nuls. Ceux-ci ne pourront pas être non plus négatifs du
fait de son caractère définie positive. Dans le cas contraire on parle de pivot
négatif ou nul, terme qui est d’ailleurs largement repris par les logiciels
lorsque le calcul est interrompu pour ces raisons.
Une fois {𝑄𝐿 } déterminé, les réactions d’appui peuvent être calculer à partir
de :
{𝑅 } = [𝐾𝐵𝐿 ] ∙ { 𝑄𝐿 } + [𝐾𝐵𝐵 ] ∙ { 𝑄𝐵 }
Calcul des effort internes en repère local :
Une fois les déplacements en repère global connus et pour déterminer les
efforts internes aux nœuds, il suffit de calculer pour chacun des éléments, les
déplacements associés en repère local puis de les injecter dans la relation de
rigidité élémentaire. On a donc :
{𝑓𝑒 }𝑖𝑛𝑡𝑒𝑟𝑛𝑒𝑠 = [𝑘𝑒 ] ∙ {𝑞𝑒 } = [𝑘𝑒 ] ∙ [𝑅𝑒 ] ∙ {𝑄𝑒 }
Vecteur charges exprimées en repère local :
Considérant une poutre plane chargée uniformément, le vecteur charge {𝑓𝑒 }
se réduit à celui des forces nodales de surface, soit
{𝑓𝑒𝑠 } = ∫ [𝑁]𝑇 ∙ { 𝑓 𝑠 } ∙ 𝑑𝑆𝑒
𝑆𝑒
𝐿
0
{𝑓𝑒 } = {𝑓𝑒𝑠 } = ∫ [𝑁 ]𝑇 ∙ {𝑓 𝑠}
∙ 𝑑𝑆𝑒 = ∫ [𝑁]𝑇 ∙ {−𝑞 } 𝑑𝑥
𝑆𝑒 0 𝑦
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 76
Dynamique des Structures
3𝑥 2 2𝑥 3 𝑞𝑦 𝐿
1− 2
− 3 −
𝐿 𝐿 2
𝐹𝑖𝑦 2𝑥 2
𝑥 3
𝑞𝑦 𝐿2
𝑀𝑖𝑧 𝐿 𝑥 − − −
= ∫ 𝐿 𝐿2 ∙ { 0 } ∙ 𝑑𝑥 = 12
𝐹𝑗𝑦 3𝑥 2
2𝑥 3 −𝑞𝑦 𝑞 𝑦𝐿
0 −
−
{𝑀𝑗𝑧 } 𝐿2 𝐿3 2
𝑥 2
𝑥 3 𝑞𝑦 𝐿2
{ − + { 12 }
𝐿 𝐿2 }
D’une manière générale et dans le cas de forces orientées dans le sens 𝑦
négatif, le vecteur charges d’un élément de poutre plane soumis à un système
de charges transversales est égal à l’inverse des réactions et moments
d’encastrement de la poutre bi-encastrée subissant le même chargement. En
d’autres termes, ce vecteur traduit les actions nodales équivalentes au
chargement appliqué sur la poutre. Le tableau suivant donne les charges
nodales équivalentes pour différents type de chargements.
Calcul des efforts internes en repère local :
Une fois les déplacements en repère global connus, il suffit de calculer pour
chacun des éléments, les déplacements associés en repère local puis de les
injecter dans la relation de rigidité élémentaire pour obtenir les efforts aux
nœuds correspondants. On a donc :
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 77
Dynamique des Structures
Avec :
{𝑞𝑒 } = [𝑅𝑒 ]{𝑄𝑒 }
{𝑓𝑒 } ∶ vecteur des charges nodales équivalente au chargement appliqué sur la
poutre.
8. Résolution :
La méthode de résolution s’inspire généralement de celle de Rayleigh – Ritz
qui consiste via une minimisation de l’énergie potentielle, en une recherche
des termes d’une fonction approximant le champ inconnu et satisfaisant aux
conditions aux limites. On définit l’énergie potentielle 𝐸 comme étant la
différence entre l’énergie de déformation 𝑊 et le travail 𝑇 des forces de
volume et de surface :
1
𝐸 = 𝑊 − 𝑇 = ∫ {𝜀 }𝑇 ∙ {𝜎} ∙ 𝑑𝑉 − ∫ {𝑢} 𝑇 ∙ { 𝑓 𝑣 } ∙ 𝑑𝑉 − ∫ {𝑢} 𝑇 ∙ {𝑓 𝑠 } ∙ 𝑑𝑆
2 𝑉 𝑉 𝑆
Le problème revient alors à rechercher un champ de déplacement minimisant
𝐸 de manière à caractériser son équilibre.
Ecole Nationale d’Ingénieurs Abderhamane Baba Touré 78