Résolution approchée de problèmes de dynamique

en régime transitoire par superposition modale

F. Louf

Dans cette fiche, on montre comment calculer une solution approchée à un problème de
dynamique transitoire par une technique de superposition modale et à l’aide de la méthode
des éléments finis. L’étude sera menée sur une barre en traction puis on montrera comment
faire cela dans CATIA.
1 Démarche générale
La démarche ressemble fortement à celle proposée au chapitre ??. Puisque, dans le cas gé-
néral, il n’est pas possible de calculer les modes propres et fréquences propres exactes, nous
commençons par calculer une base modale approchée par la méthode des éléments finis. On
obtient alors, contrairement au modèle continu, un nombre fini Nddl de modes propres. En pra-
tique, pour une structure complexe, on se limite à nombre de modes bien inférieur au nombre
de degrés de liberté.
Le calcul de la solution par superposition modale, consiste donc à chercher la solution sous la
forme (dans le cas unidimensionel) :
N
X ddl

u(x, t) = φn (x)gn (t)

n=1

Cela nécessite ensuite de résoudre des équations différentielles découplées en gn (t), ici en
nombre fini, où apparaissent au second membre les efforts imposés projetés dans la base mo-
dale. Nous verrons qu’il est possible de résoudre exactement ces équations différentielles sous
l’hypothèse, peu contraignante en pratique, d’efforts imposés linéaires par morceaux.
Dans la suite, nous illustrons cette méthode sur l’exemple de la barre en traction, soumise à un
effort en bout, de type « échelon ».

2 Rappel du problème continu

x=0 x=L
~x

F IG . 1 – La barre et les conditions aux limites

La structure étudiée dans ce chapitre est unidimensionnelle. Il s’agit d’une barre soumise à
un chargement de traction compression. Ses caractéristiques géométriques et matériau sont
notées :
– L : longueur ;
– S : section ;
– E : module d’Young ;
– ρ : masse volumique.
On note u(x, t) le déplacement d’un point d’abscisse x de la barre à l’instant t. La figure 1
présente la structure dans son environnement. Le déplacement, en x = 0, est imposé nul :

u(0, t) = 0 ∀t ∈ [0, T ]

où [0, T ] représente l’intervalle temporel d’étude. Le chargement est appliqué en bout de

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 2

poutre. On le note F (t). On obtient alors :
∂2u ∂2u
ES − ρS = 0 ∀x ∈ [0, L] ∀t ∈ [0, T ] (1)
∂x2 ∂t2
∂u
u(0, t) = 0 et ES (L, t) = F (t) ∀t ∈ [0, T ] (2)
∂x
∂u
u(x, 0) = 0 et (x, 0) = 0 ∀x ∈ [0, L] (3)
∂t

3 Formulation variationnelle
On suit la démarche classique. Pour éviter toute confusion avec le champ de vitesse, on note
le champ test u∗ . Le problème posé est alors équivalent à :
Z L 2 Z L 2
∂ u ∗ ∂ u ∗
ES 2
u dx = ρS 2
u dx ∀u∗ (4)
0 ∂x 0 ∂t
∂u
u(0, t) = 0 et ES (L, t) = F (t) ∀t ∈ [0, T ] (5)
∂x
∂u
u(x, 0) = 0 et (x, 0) = 0 ∀x ∈ [0, L] (6)
∂t
En intégrant par partie, et en utilisant la relation entre déplacement et effort imposé en x = L,
on obtient la formulation variationnelle suivante :

∀u∗ /u∗ (0, t) = 0

Z L Z L
∗ ∂u ∂u∗ ∂2u ∗
F (t)u (L) = ES dx + ρS u dx (7)
0 ∂x ∂x 0 ∂t2
u(0, t) = 0 ∀t ∈ [0, T ] (8)
∂u
u(x, 0) = 0 et (x, 0) = 0 ∀x ∈ [0, L] (9)
∂t

4 Discrétisation par éléments finis

4.1 Choix du maillage
On divise le domaine [0, L] en N segments à deux nœuds de taille identique. On exprime alors
le champ de déplacement recherché à l’aide des fonctions de base éléments finis :

u(x, t) = [N(x)]{u(t)} (10)

Le vecteur {u(t)} représente les valeurs aux nœuds du champ de déplacement à l’instant t. On
choisit u∗ (x, t) de la même forme :

u∗ (x, t) = [N(x)]{u∗ (t)} (11)

4.2 Discrétisation des équations d’équilibre

En poursuivant la démarche classiquement utilisée pour des problèmes ne faisant pas intervenir
le temps, on introduit u(x, t) et u∗ (x, t) de la forme choisie dans la formulation variationnelle.

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 3

Avant toute chose, il faut définir :
∂u ′ ∂u∗
(x, t) = [N (x)]{u(t)} et (x, t) = [N ′ (x)]{u∗ (t)} (12)
∂x ∂x
∂2u
(x, t) = [N(x)]{ü(t)} (13)
∂t2
L’équation (7) devient alors :
Z L
T T T
∗
F (t){u (t)} [N(L)] ={u (t)} ES∗
[N ′ (x)]T [N ′ (x)]dx{u(t)}
0
Z L
T
∗
+ {u (t)} ρS [N(x)]T [N(x)]dx{ü(t)} (14)
0

On nomme respectivement K, M, {F (t)}, les matrices de raideurs, de masse, et le vecteur des

forces généralisées définis par :
Z L
K = ES [N ′ (x)]T [N ′ (x)]dx {F (t)} = [N(L)]T F (t) (15)
0
Z L
M = ρS [N(x)]T [N(x)]dx (16)
0

Le problème est donc de trouver le vecteur {u(t)} avec [N(0)]{u(t)} = 0, vérifiant, ∀{u∗ (t)} :
h i
{u∗ (t)}T K{u(t)} + M{ü(t)} − {F (t)} = 0 (17)
[N(0)]{u(t)} = 0 ∀t ∈ [0, T ] (18)
{u(0)} = {0} et {u̇(0)} = {0} (19)

En conséquence, cela nous permet de nous ramener à la résolution de N + 1 équations diffé-

rentielles couplées, sur tout l’intervalle d’étude :

M{ü(t)} + K{u(t)} = {F (t)}

u(0, t) = [N(0)]{u(t)} = 0 ∀t ∈ [0, T ]
{u(0)} = {0} et {u̇(0)} = {0}

5 Calcul des modes approchés

5.1 Recherche de solutions synchrones
On recherche tout d’abord les solutions de l’équation d’équilibre éléments finis homogène de
type :
u(x, t) = g(t)[N(x)]{φ} (20)
soit encore :
{u(t)} = g(t){φ} (21)
En injectant cette forme de solution dans l’équation d’équilibre éléments finis homogène on
obtient :
M{φ}g ′′ (t) + K{φ}g(t) = 0 (22)

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 4

 5.2 Obtention des pulsations propres et des modes propres

soit encore
g ′′ (t)
− M{φ} = K{φ} (23)
g(t)
Montrons que la constante λ définie par :

g ′′ (t)
λ= (24)
g(t)

est réelle et négative. Dans le cas général où le mode propre {φ} est complexe, {φ} = {a} +
i{b}, on obtient :    
−λM {a} + i({b} = K {a} + i{b} (25)

En pré-multipliant cette dernière équation par {φ} = {a} − i{b}, on trouve :

       
−λ {a}T − i{b}T M {a} + i{b} = {a}T − i{b}T K {a} + i{b} (26)
{a}T K{a} + {b}T K{b}
λ=− (27)
{a}T M{a} + {b}T M{b}

Les matrices K et M sont symétriques et positives. Par conséquent l’expression précédente

implique que la constante λ est réelle et négative. On pose donc λ = −ω 2 .

5.2 Obtention des pulsations propres et des modes propres

Les pulsations propres sont les ωn telles que les solutions du système linéaire suivant soient
non-triviales :
(K − ω 2 M){φ} = {0} (28)
où {φ} doit toujours vérifier des conditions aux limites. En effet, celles-ci n’ont pas encore été
prises en compte. Dans le cas où elles sont homogènes, il suffit d’éliminer dans les matrices
K et M les lignes et colonnes associées aux degrés de liberté bloqués. Dans le cas présent,
on éliminera la première ligne et la première colonne. Les matrices K et M de dimensions
(N + 1) × (N + 1) deviennent donc les matrices K e et M
f de dimensions N × N.
Les pulsations propres ωn sont alors telles que les solutions du système linéaire suivant soient
non-triviales :
e − ω 2 M){
(K f φ} e = {0} (29)
n

et donc telles que :

e − ω 2 M)|
|K f =0 (30)
n

Les pulsations propres ωn du modèle éléments finis sont donc les N racines réelles et positives
d’un polynôme de degré 2N en ω. Comme dans le cas continu, on les range dans un ordre
croissant.
Une fois que les ωn sont calculées, il reste à calculer les vecteurs {φen } tels que :

(K f φen } = {0}
e − ω 2 M){ (31)
n

Le vecteur {φ} s’obtient alors sans difficulté à partir de {φen } en tenant compte de la condition
en x = 0.

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 5

 5.3 Orthogonalité

5.3 Orthogonalité
Comme pour un système continu, deux modes distincts sont orthogonaux. Soient deux modes
distincts {φn } et {φp }. Par définition, ils vérifient :

(K − ωn2 M){φn } = {0} (32)

(K − ωp2 M){φp } = {0} (33)

En pré-multipliant (32) par {φp }T et (33) par {φn }T , et en faisant la différence des équations
obtenues on obtient :

{φp }T K{φn } − {φn }T K{φp } = ωn2 {φp }T M{φn } − ωp2 {φn }T M{φp } (34)

La transposée d’une matrice 1 × 1 étant cette même matrice, on a :

{φp }T K{φn } = {φn }T KT {φp } (35)

{φp }T M{φn } = {φn }T MT {φp } (36)

Enfin, puisque les matrices de raideur et de masse sont symétriques :

0 = (ωn2 − ωp2 ){φn }T M{φp } (37)

Les modes n et p étant par hypothèse distincts, l’équation précédente conduit à :

{φn }T M{φp } = 0 pour n 6= p (38)

ce qui montre que les modes n et p sont orthogonaux. On en déduit immédiatement que :

{φn }T K{φp } = 0 pour n 6= p (39)

5.4 Masse modale et raideur modale

On définit la masse modale Mn et la raideur modale Kn par :

Mn = {φn }T M{φn } (40)

Kn = {φn }T K{φn } (41)

Le mode {φn } vérifie :

(K − ωn2 M){φn } = {0} (42)

et donc :

{φn }T K{φn } = ωn2 {φn }T M{φn } (43)

On obtient immédiatement la relation :

Kn
ωn2 = (44)
Mn

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 6

 5.5 Normalisation

5.5 Normalisation
Comme dans le problème continu, on peut choisir φn tels que :
– le vecteur propre {φn } soit de norme unitaire ;
– la masse modale associée Mn soit unitaire ;
– la raideur modale associée Kn soit unitaire.

5.6 Mise en œuvre

On considère un maillage simple constitué de 3 segments à deux nœuds. Le nombre de degré
de liberté est donc égal à 4.

5.6.1 Calcul des matrices K et M

On commence par calculer les matrices élémentaires Ke et Me qui valent :

Z Le " #
ES ES 1 −1
Ke = e [N ′ (x)]T [N ′ (x)]dx = e (45)
L 0 L −1 1
Z Le " #
1/3 1/6
Me = ρSLe [N(x)]T [N(x)]dx = ρSLe (46)
0 1/6 1/3

Par assemblage, on obtient les matrices K et M de dimensions 4 × 4, avec Le = L/N = L/3 :

 
1 −1 0 0

3ES −1 2 −1 0 

K=   (47)
L 0 −1 2 −1
0 0 −1 1
 
1/3 1/6 0 0

ρSL 1/6 2/3 1/6 0 
M=   (48)
3  0 1/6 2/3 1/6
0 0 1/6 1/3

5.6.2 Prise en compte des conditions aux limites

La condition aux limites, en x = 0 étant homogène, il suffit d’éliminer la première ligne et la

première colonne des matrices de masse et de raideur :
 
2 −1 0
Ke = 3ES  −1 2 −1

(49)
L
0 −1 1
 
2/3 1/6 0
f = ρSL 
M 1/6 2/3 1/6

(50)
3
0 1/6 1/3

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 7

 5.6 Mise en œuvre

5.6.3 Calcul des pulsations propres

e − ω2M
Calculons le déterminant de K f:

6c2 2 2 2
L2
− 2ω9 − 3c
L2
− ω18 0
e − ω 2M|
|K f =0⇔ − 3c
2
− ω18
2 6c 2
− 92ω 2 3c2 2
− L2 − ω18 = 0 (51)
L2 L2
2 2 2 2
0 − 3c
L2
− ω18 3cL2
− ω9

On obtient donc l’équation :

√ √
a3 − 2ab2 = a(a2 − 3b2 ) = a(a − 3b)(a + 3b) = 0 (52)

en ayant posé :

6c2 2ω 2
a= − (53)
L2 9
2
3c ω2
b=− 2 − (54)
L 18
Les pulsations propres sont donc les solutions positives de :

6c2 2ω 2
a=0⇔ − =0 (55)
L2 9 √ √
√ (6 − 3 3)c2 (4 + 3)ω 2
a − 3b = 0 ⇔ 2
− =0 (56)
L√ 18
√
√ (6 + 3 3)c2 (4 − 3)ω 2
a + 3b = 0 ⇔ − =0 (57)
L2 18
On obtient finalement les trois pulsations propres données dans le tableau 1 rangées par ordre
croissant. La comparaison avec les pulsations exactes calculées via le modèle continu permet
de calculer les pourcentages d’erreur commises.

Pulsation Estimée Exacte Erreur

s √
c 6−3 3 cπ
ω1 3 2 √ 1, 1%
L 4+ 3 L2
c √ c 3π
ω2 3 3 9, 3%
sL √ L 2
c 6+3 3 c 5π
ω3 3 2 √ 16, 7%
L 4− 3 L 2

TAB . 1 – Pulsations estimées par éléments finis, à l’aide d’un maillage à 3 éléments, comparées
aux pulsations exactes obtenues.

5.6.4 Calcul des modes propres

Le mode propre {φen } associé à la valeur propre ωn vérifie :

(K f φen } = {0}
e − ω 2 M){ (58)
n

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 8

 5.6 Mise en œuvre

Calcul de {φ1 } Explicitons les termes a et b :

√ √
6c2 2ω12 6c2 18c2 6 − 3 3 c2 18 3
a= 2 − = 2 − 2 √ = 2 √ (59)
L 9 L L 4+ 3 L 4+ 3
√
3c2 ω12 3c2 c2 6 − 3 3 c2 18
b=− 2 − =− 2 − 2 √ =− 2 √ (60)
L 18 L L 4+ 3 L 4+ 3
Le système à résoudre devient alors :
√    
2
3 −1 0   φ1 
 
 0
 √  3
−1 3 −1  φ1 = 0 (61)
√   4  
  
0 −2 3 φ1 0
Bien-entendu, les équations du système ne sont pas indépendantes. On peut, par exemple fixer
φ41 = 2. Rappelons que φ11 = 0. On obtient alors :
 

 0 

 1  
{φ1 } = √ (62)

 3 


  
2

Calcul de {φ2 } Explicitons les termes a et b :

6c2 2ω22 6c2 2 27c2
a= − = − =0 (63)
L2 9 L2 9 L2
3c2 ω 2 3c2 c2 27 c2 9
b=− 2 − 2 =− 2 − 2 =− 2 (64)
L 18 L L 18 L 2
Le système à résoudre devient alors :
    
2
0 −1 0  φ2     0
  3
 −1 0 −1  φ 2 = 0 (65)
 4
     
0 −1 0 φ2 0
On peut, par exemple fixer φ42 = 2. Rappelons que φ12 = 0. On obtient alors :
 

 0 

2 
{φ2 } = (66)

 0 

  
−2

Calcul de {φ3 } Explicitons les termes a et b :

√ √
6c2 2ω32 6c2 18c2 6 + 3 3 c2 18 3
a= 2 − = 2 − 2 √ =− 2 √ (67)
L 9 L L 4− 3 L 4− 3
√
3c2 ω32 3c2 c2 6 + 3 3 c2 18
b=− 2 − =− 2 − 2 √ =− 2 √ (68)
L 18 L L 4− 3 L 4− 3
Le système à résoudre devient alors :
√    
2
3 1 0   φ3 
 
 0
 √  3
 1 3 1  φ3 = 0 (69)
√   4  
  
0 2 3 φ3 0

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 9

 5.6 Mise en œuvre

Bien-entendu, les équations du système ne sont pas indépendantes. On peut, par exemple fixer
φ43 = 1. Rappelons que φ13 = 0. On obtient alors :
 

 0 


 1  
{φ3 } = √ (70)

 − 3


 

2

Tracé des modes On a représenté les trois modes calculés précédemment sur la figure 2. Ces
modes sont à comparer aux modes exacts tracés sur la figure 3.

1.5 2 2.5 3 3.5 4

Mode 1
-1
Mode 2
Mode 3

-2

F IG . 2 – Les trois modes calculés par éléments finis

1.0

0.5

0.2 0.4 0.6 0.8 1

Mode 1
-0.5 Mode 2
Mode 3

-1.0

F IG . 3 – Les quatre premiers modes pour une barre de longueur unitaire

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 10

6 Calcul de la réponse transitoire par superposition modale
6.1 Equation d’équilibre éléments finis
Nous avons montré que l’équation d’équilibre éléments finis était de la forme :

M{ü(t)} + K{u(t)} = {F (t)} (71)

6.2 Forme de la solution recherchée

On cherche une solution vérifiant les équations d’équilibre éléments finis de la forme :
N
X
{u(t)} = gn (t){φn } (72)
n=1

où N est le nombre de degrés de liberté du modèle éléments finis. Il correspond ici au nombre
d’éléments dans le cas unidimensionnel puisque un degré de liberté est bloqué par l’encastre-
ment. Les vecteurs {φn } sont les modes propres du modèle éléments finis. En dérivant deux
fois, on obtient le vecteur des valeurs nodales de l’accélération :
N
X
{ü(t)} = gn′′ (t){φn } (73)
n=1

6.3 Equations différentielles vérifiées par gn

En injectant cette forme de solution dans l’équation d’équilibre éléments finis, on obtient :
N
X N
X
M gn′′ (t){φn } + K gn (t){φn } = {F (t)} (74)
n=1 n=1

En multipliant à gauche par un mode propre {φp }T on obtient, du fait de l’orthogonalité des
modes :
{φn }T M{φn }gn′′ (t) + {φn }T K{φn }gn (t) = {φn }T {F (t)} (75)
Il faut donc maintenant résoudre N équations différentielles découplée du second ordre en gn :

Mn gn′′ (t) + Kn gn (t) = Fn (t) (76)

où Mn et Kn sont les respectivement les masses modales et raideurs modales associées au

mode n avec :
Kn
ωn2 = (77)
Mn
ce qui conduit enfin à :
gn′′ (t) + ωn2 gn (t) = fn (t) (78)
Le terme fn est appelé facteur de participation du mode n à l’excitation.

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 11

 6.4 Résolution des équations différentielles

6.4 Résolution des équations différentielles

On fait ici l’hypothèse que l’effort imposé varie linéairement par morceaux, ce qui en pratique
n’est pas très contraignant. On introduit donc une discrétisation temporelle permettant de re-
présenter les efforts extérieurs. On note ∆t le pas de temps, tr les piquets de temps. On écrit
la linéarité du facteur de participation fn :
tr+1 − t t − tr
fn (t) = fnr + fnr+1 (79)
∆t ∆t
L’équation différentielle à résoudre sur un pas devient alors :
tr+1 − t t − tr
gn′′ (t) + ωn2 gn (t) = fnr + fnr+1 (80)
∆t ∆t
gn (tr ) = gnr (81)
gn′′ (tr ) = ġnr (82)

où les deux dernières équations représentent la continuité de gn et de sa dérivée temporelle

gn′ . Puisque les conditions initiales à t = 0 sont connues, il suffit d’intégrer explicitement le
problème (80, 81, 82) sur un pas, et par récurrence, toutes les valeurs gnr et ġnr seront connues.

Solution de l’équation homogène associée à (80) La solution générale de l’équation homo-

gène associée à (80) est du type :

gn (t) = Acos(ωn (t − tr )) + Bsin(ωn (t − tr )) (83)

Solution particulière de (80) Une solution particulière de l’équation (80) est linéaire sur le
pas de temps :
1 tr+1 − t t − tr
gn (t) = 2 (fnr + fnr+1 ) (84)
ωn ∆t ∆t

Solution complète du problème (80, 81, 82) La solution complète est la somme de la so-
lution particulière et de la solution homogène. Il reste maintenant à prendre en compte les
conditions initiales.
1 r
gn (tr ) = gnr = A + f (85)
ωn2 n
1 f r+1 − fnr
gn′ (tr ) = ġnr = Bωn + 2 n (86)
ωn ∆t

On en déduit :
1 r
A = gnr − f (87)
ωn2 n
ġ r 1 f r+1 − fnr
B= n − 3 n (88)
ωn ωn ∆t

Explicitons maintenant les relations donnant gnr+1 et ġnr+1 en fonction des efforts connus sur le
pas fnr et fnr+1 , et des conditions initiales gnr et ġnr .

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 12

 6.5 Cas d’un chargement en déplacement

On trouve, en posant θn = ωn ∆t :

sin(θn )
gnr+1 =cos(θn )gnr + ∆tġnr
θn
 sin(θ )  1 sin(θn ) r+1
n
+ − cos(θn ) fnr + 2 (1 − )fn (89)
θn ωn θn
∆tġnr+1 = − θn sin(θn )gnr + cos(θn )∆tġnr
1
+ 2 (θn sin(θn ) + cos(θn ) − 1)fnr + (1 − cos(θn ))fnr+1 (90)
ωn

La relation de récurrence donnant le vecteur {gnr+1, ∆tġnr+1}T en fonction des vecteurs connus
{gnr , ∆tġnr }T et {fnr , fnr+1}T est donc :
( ) " #( )
sin(θn )
gnr+1 cos(θn ) θn
g r
n
=
∆tġnr+1 −θn sin(θn) cos(θn) ∆tġnr
" #( )
sin(θn ) sin(θn ) r
1 θn
− cos(θn ) 1 − θn
f n
+ 2 (91)
ωn θn sin(θn) + cos(θn ) − 1 1 − cos(θn) fnr+1

On peut donc conclure que, sous l’hypothèse d’efforts imposés linéaires par morceaux, l’inté-
gration des équations différentielles est exacte et ne nécessite aucune résolution complexe.

6.5 Cas d’un chargement en déplacement

Dans le cas où est aucun effort n’est imposé et où un déplacement sur une zone donnée est
imposé variable dans le temps, par exemple ici en x = L, on obtient comme précédemment
des équations différentielles (homogènes cette fois) du second ordre en gn à résoudre. Par
contre, elles sont couplées par la condition aux limites. On a en effet la relation en x = L :
N
X
T +1
[N(L)] {u(t)} = gn (t)φN
n = ud (t) (92)
n=1

+1
où φNn désigne la valeur du mode n au dernier noeud de la barre. On ne pourra trouver de
solutions à ces équations que si la forme temporelle du déplacement peut être représentée par
une base de fonctions oscillantes de dimensions finie. Dans le cas général, il faudrait donc faire
une approximation en imposant le bon déplacement aux piquets de temps, mais par les valeurs
intermédiaires. Les solutions générales de chacune des équations différentielles sont donc :

gn (t) = An cos(ωn (t − tr )) + Bn sin(ωn (t − tr )) (93)

telles que :
N
X
+1
gn (tr )φN
n = urd (94)
n=1
XN
+1
gn (tr+1 )φN
n = ur+1
d (95)
n=1

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 13

 6.6 Troncature de la base modale

En remplaçant gn par sa forme générale, on obtient le système d’équations :

N
X
+1
An φN
n = urd (96)
n=1
XN
+1
An cos(ωn ∆t) + Bn sin(ωn ∆t)φN
n = ur+1
d (97)
n=1

Les équations de continuité entre deux pas de temps imposant déjà un nombre d’équations
suffisant pour trouver tous les An et Bn , on ne pourra pas trouver de solution.

6.6 Troncature de la base modale

Souvent, le nombre de degré de liberté est trop important et on limite le nombre de mode pris
en compte dans la base modale, c’est à dire que l’on écrit :
e
N
X
{u(t)} = e<N
gn (t){φn } avec N (98)
n=1

6.7 Mise en œuvre

Le maillage que nous utiliserons dans cette partie comporte trois segments. Les éléments sont
linéaires de sorte que les résultats établis au paragraphe 5.6 peuvent être utilisés. On connaît
donc déjà :
ω1 , ω2 , ω3 et {φ1 }, {φ2 }, {φ3} (99)
Le chargement ici considéré est un effort de type échelon en bout de poutre. Nous connaissons
donc la solution exacte du problème posé.

6.7.1 Calcul des masses modales

Puisque les participations modales fn (t) sont définies par :

{φn }T {F (t)}
fn (t) = (100)
Mn
il nous faut tout d’abord calculer les trois masses modales :

Mn = {φn }T M{φn } (101)

Rappelons que la matrice de masse M et les vecteurs {φn } vallent dans notre cas :
       
4 1 0 
 1 
 
 1  
 1 
ρSL   √ √ 
M= 1 4 1 {φ1 } = 3 {φ2 } = 0 {φ3 } = − 3 (102)
18 
   
   
 

0 1 2 2 −1 2

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 14

 6.7 Mise en œuvre

On obtient alors :
  
i 4 1 0   1  
T ρSL h √   √
M1 = {φ1 } M{φ1 } = 1 3 2 1 4 1 3 (103)
18 
  
0 1 2 2
  
i 4 1 0  1 
T ρSL h  
M2 = {φ2 } M{φ2 } = 1 0 − 1 1 4 1 0 (104)
18 
  
0 1 2 −1
  
h i 4 1 0   √ 1 

T ρSL √  
M3 = {φ3 } M{φ3 } = 1 − 3 2 1 4 1 − 3 (105)
18 
 

0 1 2 2

Les masses modales obtenues sont alors :

ρSL √
M1 = (4 + 3) (106)
3
ρSL
M2 = (107)
3
ρSL √
M3 = (4 − 3) (108)
3

6.7.2 Calcul des participations modales à l’excitation

Les participations modales sont maintenant données par l’équation (100). En tenant compte du
fait que l’effort imposé est de type échelon et donc égal à F0 pour t ≥ 0, on trouve :
 
h i 
0 
1 T 3F0 √ 6F0
f1 (t) = {φ1 } {F (t)} = 1 3 2 0 = √ (109)
M1 ρSL 
   ρSL(4 + 3)
1
 
1 3F0 h i0  3F0
T
f2 (t) = {φ2 } {F (t)} = 1 0 −1 0 =− (110)
M2 ρSL 
   ρSL
1
 
1 3F0 h √ i0  6F0
T
f3 (t) = {φ3 } {F (t)} = 1 − 3 2 0 = √ (111)
M3 ρSL 
   ρSL(4 − 3)
1

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 15

 6.7 Mise en œuvre

Préliminaires
→ Définir une discrétisation temporelle permettant de représenter les efforts imposés
Boucle sur les modes
( )
gn (0)
→ Calculer le vecteur {Grn } =
∆tġn (0)
→ Calculer θn = ωn ∆t " #
sin(θn )
cos(θn ) θn
→ Former la matrice An =
−θn sin(θn) cos(θn)
" #
sin(θn ) sin(θn )
θn
− cos(θ n ) 1 − θn
→ Former la matrice Bn =
θn sin(θn) + cos(θn ) − 1 1 − cos(θn)
→ Boucle sur les piquets de temps
( )
fnr
→ Calculer le vecteur {Fnr,r+1} =
fnr+1
→ Calculer le vecteur {Gr+1 r r,r+1
n } = An {Gn } + Bn {Fn }
r+1 r
→ Décaler les indices {Gn } → {Gn }
→ Incrémenter r
→ Incrémenter n

TAB . 2 – Algorithme associé à l’intégration explicite des équations de mouvement

6.7.3 Intégration des équations différentielles

Les participations modales fn (t) sont ici constantes du fait que l’effort imposé est un éche-
lon. On pourrait donc se contenter d’un pas de temps pour calculer la solution finale gn (T )
et ġn (T ) à partir des conditions initiales. Toutefois, comme on s’intéresse à l’évolution des
déplacements et vitesses dans la structure au cours du temps, on choisit de diminuer ce pas de
temps pour obtenir plus d’informations intermédiaires.
Le programme permettant de calculer les trois vecteurs {gnr , ġnr }T à chaque piquet de temps est
donné dans le tableau 2.

6.7.4 Expression de la réponse

Les vecteurs des valeurs nodales des champs de déplacement et de vitesse sont obtenus via la
relation (72) :
3
X
{u(tr )} = gnr {φn } (112)
n=1
X3
{u̇(tr )} = ġnr {φn } (113)
n=1

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 16

 6.7 Mise en œuvre

u(x, t) u̇(x, t)
1.5
2 2
1
1.5
t 1

1.5
t
0.5 0
0
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

(a) N = 3, Nt = 10

u(x, t) u̇(x, t)
1.5
2 1 2
1
1.5
t 1.5
t
0.5 0
0
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

(b) N = 10, Nt = 30

u(x, t) u̇(x, t)
2

1.5 2 1 2
1
1.5
t 0 1.5
t
0.5
0
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

(c) N = 20, Nt = 60

u(x, t) u̇(x, t)
2
1.5
1.5 2 1 2
1
1.5
t 0.5
1.5
t
0.5 0
0 -0.5
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

(d) N = 40, Nt = 120

F IG . 4 – Effet d’un raffinement du maillage

Les champs de déplacement et vitesse sont obtenus à l’aide des vecteurs des valeurs nodales

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 17

 6.7 Mise en œuvre

et des fonctions de formes éléments finis :

3
X
u(x, tr ) = [N(x)] gnr {φn } (114)
n=1
X3
∂u
(x, tr ) = [N(x)] ġnr {φn } (115)
∂t n=1

6.7.5 Tracé des résultats

A titre d’exemple on considère un chargement de type « échelon ». Dans un premier temps on

étudie l’effet d’un raffinement du maillage. Ensuite on s’intéresse à l’effet d’une troncature de
la base modale sur la qualité de la solution approchée, en déplacement et en vitesse. Pour cela
on définit deux erreurs mesurant un écart avec la solution de référence connue.

Etude de convergence Sur la figure 4, on a tracé les champs de déplacement et de vitesse

obtenus pour différentes discrétisations spatiales et temporelles. Nous avons choisi de conser-
ver le rapport entre le nombre d’éléments et le nombre de piquets de temps, pour tracer les
différentes figures. Le premier cas correspond à l’étude mené au paragraphe précédent.

Effet d’une troncature modale Sur la figure 5, on étudie l’influence d’une troncature de la
base modale éléments finis. On considère un maillage à 40 éléments, et on trace les champs de
déplacement et de vitesse sur le domaine [0, L] × [0, T ] pour différents nombres de modes pris
en compte.
L’allure générale du champ de déplacement ne semble que peu affectée. Par contre lorsque le
nombre de mode pris en compte augmente, des oscillations de hautes fréquences sont visibles
au niveau du champ de vitesse. Il semble donc que ne pas prendre en compte la totalité des
modes pourrait améliorer la qualité de la solution approchée en vitesse. Pour conforter cette
impression, on définit une erreur absolue en déplacement et une erreur absolue en vitesse
mesurant l’écart entre la solution exacte (uex , u̇ex) et la solution approchée (uNe , u̇Ne ) :

XN X Nt h  iL jT   iL jT i2
e) =
eu (N uex , − uNe , (116)
i=0 j=0
N Nt N Nt
XN X Nt h  iL jT   iL jT i2
e) =
eu̇ (N u̇ex , − u̇Ne , (117)
i=0 j=0
N Nt N Nt

avec :
e ≤N
N
X
u(xi , tj ) = [N(xi )] gnj {φn } (118)
n=1
e ≤N
N
∂u X
(xi , tj ) = [N(xi )] ġnj {φn } (119)
∂t n=1

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 18

 6.7 Mise en œuvre

u(x, t) u̇(x, t)
2

1.5 1
2 2
1
1.5
t 0.5
1.5
t
0.5
0
0
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

e = 10
(a) N

u(x, t) u̇(x, t)
2 1.5
1.5 2 1 2
1 t 0.5
1.5
t
0.5 1.5 0
0 -0.5
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

e = 20
(b) N

u(x, t) u̇(x, t)
2 1.5
1.5 2 1 2
1 t 0.5
1.5
t
0.5 1.5 0
0 -0.5
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

e = 30
(c) N

u(x, t) u̇(x, t)
2
1.5
1.5 2 1 2
1 t 0.5 t
0.5 1.5 0 1.5

0 -0.5
0 1 0 1
0.2 0.2
0.4 0.4
0.5 0.5
0.6 0.6
0.8 0.8
x 1 0 x 1 0

e = 40
(d) N

F IG . 5 – Effet d’une troncature de la base modale éléments finis

On constate tout d’abord que la qualité du champ de déplacement augmente régulièrement

lorsque Ne croît. Ainsi le minimum de l’erreur eu est atteint pour N e = N = 40. Par contre, en
ce qui concerne le champ de vitesses, le résultat n’est pas si clair. La perturbation de la réponse
en vitesse par les modes à haute fréquence se traduit par une fonction eu̇ dont le minimum n’est

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 19

pas atteint pour N e = N = 40, mais pour, dans ce cas N e = 14. Ce résultat peut s’expliquer
par le fait que pour un maillage donné, les premières pulsations sont très correctement estimés.
Il en est d’ailleurs de même pour les modes. Par contre la qualité des pulsations suivantes se
dégrade progressivement. La figure 7 illustre ce point pour un maillage à 40 éléments. Ainsi,
tenir compte de tous les modes dans la base modale revient à tenir compte de modes de très
médiocre qualité et peut donc, finalement dégrader la solution.
Erreur en déplacement
0.25

0.2

0.15

0.1

0.05

e
N
10 20 30 40

(a) Evolution de l’erreur en déplacement

Erreur en vitesse
50

40

30

20

10

e
N
10 20 30 40

(b) Evolution de l’erreur en vitesse

e et vitesse eu̇ (N)

F IG . 6 – Evolution des erreurs en déplacement eu (N) e en fonction du nombre
de modes pris en compte dans la base modale éléments finis

7 En pratique dans Catia

On suppose ici que le maillage, les propriétés et matériaux de la structure étudiée sont déjà
définis.

7.1 Calcul de la base modale

Pour commencer, il faut calculer la base modale :

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 20

 7.2 Définition du chargement

Pulsations ωn
140

120 Pulsations approchées

Pulsations exactes
100

80

60

40

20
Numéro du mode
10 20 30 40

F IG . 7 – Evolutions des pulsations exactes et estimées en fonction du numéro du mode

– Insérer un Cas de fréquence à partir du menu Insertion et laisser cochée la case Conditions
aux limites ;
– Préciser les conditions aux limites ;
– En double cliquant sur la Solution Modale contenue dans le Cas de fréquences, préciser le
nombre de modes à calculer, c’est-à-dire la taille de la base dans laquelle on cherchera la
réponse dynamique dans la suite ; bien évidemment, ce nombre de modes ne peut excéder
le nombre de degrés de liberté du problème éléments finis.

7.2 Définition du chargement

Ensuite, il faut définir la zone et la variation temporelle du chargement. La zone est définie dans
un Cas statique qu’il n’est pas nécessaire de calculer, la dépendance temporelle est définie dans
un fichier à deux colonnes à importer.
– Insérer un Cas statique à partir du menu Insertion ; si l’on veut, on peut faire pointer les
conditions aux limites vers celles du Cas de fréquences ;
– Créer le chargement du type souhaité : force surfacique, distribuée, ou autre, selon le type
de modèle utilisé ; les déplacements imposés ne sont pas acceptées dans la suite de la dé-
marche ; ce point est tout à fait logique étant donné la méthode utilisée (cf paragraphe 6.5).

7.3 Calcul de la réponse en régime transitoire

Enfin, il faut définir le calcul de dynamique transitoire : base modale utilisée, chargement
(direction, zone d’application, type, et dépendance au temps), amortissement.
– Insérer un Cas de réponse dynamique en temps à partir du menu Insertion ; préciser que la
base modale utilisée est celle contenue dans le Cas de fréquences ;
– Importer le fichier décrivant la dépendance temporelle de l’effort à l’aide de l’outil Modula-
tion en temps ; ce fichier doit avoir l’aspect proposé dans le tableau 3 ;
– les instants repérés ne sont pas forcément régulièrement espacés ;
– l’instant initial est donné t = 0s ;
– l’instant final de l’étude souhaitée doit correspondre au dernier instant donné.

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 21

 7.4 Post traitement

– En double cliquant sur Excitation des charges, dans la fenêtre qui apparaît :
– préciser le chargement à exciter en cliquant dans l’arbre sur le chargement défini dans le
cas statique ;
– préciser la modulation en temps en cliquant sur la modulation en temps insérée précé-
demment (elle se trouve sous le noeud Modulations de l’arbre.
– En double cliquant sur Amortissement dans l’arbre, on peut préciser le type d’amortissement
et le coefficient d’amortissement retenu ; par défaut l’amortissement est modal et la valeur
donnée est identique pour chaque mode : 1%.

Remarque. – On aurait pu insérer plusieurs chargements dans un cas statique et définir plu-
sieurs modulations associées.

t(s) F
0 0
1E-6 1
5E-6 1
6E-6 0
1E-5 0

TAB . 3 – Exemple de fichier donnant la modulation temporelle de l’effort imposé : cas d’un
effort injecté de type carré

7.4 Post traitement

Le post traitement peut se faire de différentes façons :
– il est toujours possible de générer une image (cliquer droit sur Solution de réponse dyna-
mique transitoire puis Génération d’images) présentant les déplacements, vitesses, accé-
lérations, contraintes, etc et de spécifier dans la fenêtre ouverte l’Occurence, c’est à dire
l’instant auquel on souhaite représenter le champ choisi ;
– il est ensuite possible de jouer un film en cliquant sur l’icône Animation ; dans ce cas, pour
voir l’évolution temporelle d’un champ spatial préalablement tracé, il faut spécifier (bouton
Plus) que ce sont toutes les occurences qui doivent être tracées ;
– il est également possible de tracer l’évolution temporelle d’une quantité de type déplacement
locale ; pour cela, cliquer droit sur Solution de réponse dynamique transitoire et choisir
Générer un affichage 2D ; dans la fenêtre qui s’ouvre :
– cliquer sur fin
– puis dans la fenêtre Sélection de données, choisir le(s) degré(s) de liberté à tracer : sélec-
tionner le noeud souhaité sur le maillage qui s’affiche dans la fenêtre de calcul (fichier
.CATAnalysis) et préciser la ou les directions à tracer ;
– dans l’onglet configuration, préciser si l’on souhaite tracer le déplacement, la vitesse ou
l’accélération ;

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 22

 7.5 Comparaison

– dans la fenêtre d’affichage 2D, un graphique est tracé : il est possible de modifier les axes
(échelles, titres, etc) en double cliquant dessus ;
– on peut ensuite exporter la courbe au format texte (cliquer droit, puis choisir exporter).

7.5 Comparaison
Une comparaison a été effectuée entre différents calculs éléments finis et une solution de réfé-
rence obtenue via un modèle continu non amorti. Les différents calculs éléments finis corres-
pondent à :
– un maillage à trois éléments et une base modale contenant les trois premiers modes ;
– un maillage à 50 éléments et une base modale contenant les trois premiers modes ;
– un maillage à 50 éléments et une base modale contenant les dix premiers modes.
La figure 8 montre la convergence des déplacements éléments finis vers la solution issue du
modèle continu.

F IG . 8 – Evolutions du déplacement en bout en fonction du temps : solution exacte et solutions

approchées obtenues avec différents maillages et tailles de base modale

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 23

 7.5 Comparaison

Icône Nom de l’outil Description sommaire

Force distribuée Appliquer une force en N de façon distribuée sur un support
Pression Appliquer une pression en N/m2 sur une surface
Encastrement Bloquer tous les ddl d’un support
Glissement surfacique Bloquer tous les ddl selon la normale à un support
Contraintes avancées Bloquer les ddl choisis d’un support
Modulation en temps Préciser la dépendance au temps de l’effort appliqué via un fichier
Calcul Lancer le calcul de toutes ou d’une partie des analyses
Animation Animer le tracé d’un champ dépendant du temps

TAB . 4 – Outils utilisés dans Catia et icônes correspondantes dans les ateliers utilisés

– Calcul dynamique dans Catia – Eléments de théorie - F. Louf 24