# Etude du crissement d'un frein à disque ferroviaire. Réduction de modèle en dynamiques transitoire et stationnaire

Andréa Loyer, Jean-Jacques Sinou, Olivier Chiello, Xavier Lorang

► To cite this version:

Andréa Loyer, Jean-Jacques Sinou, Olivier Chiello, Xavier Lorang. Etude du crissement d'un frein à disque ferroviaire. Réduction de modèle en dynamiques transitoire et stationnaire. 10e colloque national en calcul des structures, May 2011, Giens, France. hal-00592943

HAL Id: hal-00592943

https://hal.science/hal-00592943v1

Submitted on 3 May 2011

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Étude du crissement d'un frein à disque ferroviaire
## Réduction de modèle en dynamiques transitoire et stationnaire

A. Loyer<fnref n="1" />,<fnref n="2" />,<fnref n="3" />, J.-J. Sinou<fnref n="1" />, O. Chiello<fnref n="2" />,
X. Lorang<fnref n="3" />

**Résumé** — Ce travail s’intéresse à la simulation numérique du crissement des freins à disque dans le domaine ferroviaire. Le modèle retenu considère une déstabilisation du système par couplage de modes induite par le frottement entre le disque et les garnitures. Compte tenu des temps de calculs nécessaires pour obtenir une solution dans le domaine temporel, la mise en place de technique de réduction est impérative. L’influence de la réduction sur les réponses vibratoires non-linéaires obtenues est étudiée sur un modèle éléments finis simplifié d’un frein TGV.

**Mots clés** — crissement, contact, frottement, instabilité.

## 1 Introduction

Le crissement des freins à disque est un phénomène complexe qui préoccupe de plus en plus les industriels tant dans l'automobile que le ferroviaire ou encore l'aéronautique. Ce phénomène a fait l'objet de nombreuses études, plus particulièrement dans les domaines automobiles [5, 3, 4] et aéronautiques [2] et n'a n'a été abordé que récemment l'objet d'études dans le domaine ferroviaire [9]. Le bruit de crissement est lié au frottement et peut être expliqué par deux types d'approches principales. La première se fonde sur des aspects tribologiques tels que des variations de coefficients de frottement [13]. La deuxième considère une déstabilisation de la structure par le couplage de modes induit par le frottement [14, 6]. C'est cette seconde explication qui est retenue dans cette étude. La description du contact frottant se limitera donc au contact Signorini-Coulomb à coefficient constant.

Afin de réduire les coûts de calculs liés à la simulation temporelle en dynamique non-linéaire, le calcul sur base réduite de la solution sera étudié. En particulier, cette étude s’intéresse à montrer les effets d’une telle réduction sur le comportement non-linéaire du système. Dans un travail précédent [8], Lorang a montré qu’une solution temporelle d’un système de frein TGV simplifié était convenablement représentée a posteriori par projection sur la base des modes complexes issus de l’étude de stabilité du système. Ces modes seront donc utilisés pour la réduction du système. Dans ce travail, le comportement d’un modèle de frein TGV simplifié est étudié. Dans un premier temps, une étude de stabilité sur ce système est menée, puis la réponse temporelle est calculée avec et sans réduction. L’impact de la réduction sur le calcul de la solution de systèmes déstabilisés par le frottement est mis en évidence. Cet article est divisé en deux parties. La première expose les choix de modélisation du contact et les algorithmes de résolution utilisés, la seconde présente le modèle de frein simplifié utilisé et les résultats concernant la réduction pour le calcul de sa réponse temporelle.

## 2 Modélisation du problème

### 2.1 Traitement du contact

On considère un milieu élastique compris dans un domaine $\Omega$ de frontière $\partial\Omega$ admettant la partition $\{\partial\Omega_f, \partial\Omega_u, S_c\}$, où $\partial\Omega_f$ est la partie sur laquelle un effort surfacique $f_s$ est appliqué, $\partial\Omega_u$ celle soumise à un déplacement imposé $u_d$ et $S_c$ une interface de contact frottant potentiel. Un jeu à l'interface de contact $g_0$ est pris en compte, ainsi qu'une vitesse de glissement Eulérienne $v_g$. Dans le cas du contact

<footnote n="1">¹ LTDS, École Centrale Lyon, France, jean-jacques.sinou@ec-lyon.fr</footnote>

<footnote n="2">² LTE, IFSTTAR Bron, France, olivier.chiello@ifsttar.fr</footnote>

<footnote n="3">³ SNCF, Direction de l'Innovation et de la Recherche, France, {andrea.loyer,xavier.lorang}@sncf.fr</footnote>
entre deux milieux déformables, les champs à l'interface désignent les champs relatifs. La dynamique du milieu peut être décrite par l'équation (1) où $\mathcal{A}$ et $\mathcal{B}$ sont respectivement les tenseurs de contraintes élastiques et visqueuses du matériau et où $\mathbf{r}$ est la réaction de contact.

$$
\left\{
\begin{array}{l}
\text{trouver } \mathbf{u} \in \mathcal{U} \text{ et } \mathbf{r} \text{ tels que } \forall \dot{\mathbf{u}}^* \in \mathcal{U}^0 \\
\\
\int_{\Omega} \rho \dot{\mathbf{u}}^* \ddot{\mathbf{u}} d\Omega + \int_{\Omega} \varepsilon^* (\mathcal{B}\dot{\varepsilon} + \mathcal{A}\varepsilon) d\Omega = \int_{\partial\Omega_f} \dot{\mathbf{u}}^* \cdot \mathbf{f}_s dS + \int_{\Omega} \dot{\mathbf{u}}^* \cdot \mathbf{f}_v d\Omega \\
\qquad + \int_{S_c} \dot{\mathbf{u}}^* \cdot \mathbf{r} dS \\
\\
(\dot{\mathbf{u}}, \mathbf{r}) \text{ vérifiant les conditions de frottement sur } S_c \qquad (1) \\
\\
(\mathbf{u}, \mathbf{r}) \text{ vérifiant le contact unilatéral sur } S_c \\
\qquad \mathcal{U} = \{\mathbf{u}, \mathbf{u} = \mathbf{u}_d \text{ on } \partial\Omega_u\} \\
\qquad \mathcal{U}^0 = \{\mathbf{u}, \mathbf{u} = 0 \text{ on } \partial\Omega_u\} \\
\\
\varepsilon \text{ et } \varepsilon^* \text{ sont les gradients symétriques de } \mathbf{u} \text{ et } \mathbf{u}^*
\end{array}
\right.
$$

Si la vitesse de glissement n'est pas une translation rectiligne uniforme, des termes convectifs apparaissent dans cette équation. Ici, le choix est fait de les négliger.

En décomposant un champ **x** sous la forme **x** = $x_n$$\mathbf{n} + \Pi_t$$\mathbf{x}$, avec **n** le vecteur normal à la surface de contact et $\Pi_t$ l'opérateur de projection dans le plan tangent au contact, la modélisation du contact frottant à partir des lois de Signorini et de Coulomb peut se mettre sous la forme donnée par les équations (2,3), où $K_\mu$ est le cône de Coulomb.

$$r_n = \operatorname{proj}_{\mathbb{R}^{-}} (r_n - \rho_n^u(u_n - g)), \forall \rho_n^u > 0 \quad (2)$$

$$\mathbf{r}_t = \Pi_t \operatorname{proj}_{K_\mu} (\mathbf{r} - \rho_t (\dot{\mathbf{u}} - \mathbf{v}_g)), \forall \rho_t > 0 \quad (3)$$

Les paramètres $\rho_n^u$ et $\rho_t$ sont appelés paramètres d'augmentation.

La discrétisation du problème par éléments finis en utilisant des éléments linéaires sur la zone de contact permet de vérifier le contact en vérifiant (2,3) en chaque nœud. Le problème à résoudre est alors donnée par :

$$
\left\{
\begin{array}{l}
[M] \{\ddot{U}\} + [C] \{\dot{U}\} + [K] \{U\} = \{F\} + \{R_C\} \\
(\{U\}, \{\ddot{U}\}, \{R_C\}) \text{ vérifiant (2) et (3) à chaque nœud}
\end{array}
\right.
\qquad (4)
$$

où $[M]$, $[K]$, et $[C]$ sont les matrices de masse, raideur et amortissement du système, $\{F\}$ et $\{R_C\}$ sont les vecteurs force généralisée correspondant à l'effort extérieur $\mathbf{f}_s$ et à la réaction de contact $\mathbf{r}$.

Pour la résolution du problème de dynamique, un schéma d'intégration du premier ordre du type $\theta$-méthode avec un pas de temps $\tau$ est choisi pour ne pas avoir à traiter des accélérations potentiellement non définies aux instants de choc. On notera $\{F_{MTD}\}$ l'effort issu du schéma d'intégration comprenant la participation de l'effort extérieur et des termes d'inertie et $[J]$ la matrice de l'intégrateur. Les réactions de contact sont traitées comme des impulsions dont l'intégrale sur le pas de temps est donnée par la valeur à la fin de l'intervalle multipliée par la durée $\tau$.

Ainsi le problème à résoudre pour passer du pas de temps $i$ au pas $i+1$ est :

$$
\left\{
\begin{array}{l}
[J] \{\dot{U}\}^{i+1} = \{F_{MTD}\}^{i+1} + \{R_C\}^{i+1} \\
\\
[J] = \frac{[M]}{\tau} + \theta([C] + \tau\theta[K]) \\
\\
\{F_{MTD}\}^{i+1} = (\theta\{F\}^{i+1} + (1-\theta)\{F\}^i) \\
\qquad + ([M] - \tau(1-\theta)[C] + \tau\theta[K]) \frac{\{\dot{U}\}^i}{\tau} - [K]\{U\}^i
\end{array}
\right.
\qquad (5)
$$

Les lois de contact formulées en terme de projection permettent de définir un algorithme de résolution de type point fixe :
$$
\mathbf{r} = f(\mathbf{r}, \rho_n^u, \rho_t) \tag{6}
$$

où la fonction *f* est définie par :

$$
\begin{equation}
\begin{aligned}
f : (\mathbf{r}, \rho_n^u, \rho_t) \mapsto & \operatorname{proj}_{K_\mu} \left( \mathbf{r}_t + \mathbf{n} \operatorname{proj}_{\mathbb{R}^-} \left( r_n - \rho_n^u (u_n(\mathbf{r}) - g_0) \right) \right. \\
& \left. - \rho_t (\dot{\mathbf{u}}(\mathbf{r}) - \mathbf{v}_g) \right) + \mathbf{n} \operatorname{proj}_{\mathbb{R}^-} \left( r_n - \rho_n^u (u_n(\mathbf{r}) - g_0) \right)
\end{aligned}
\tag{7}
\end{equation}
$$

En introduisant un intégrateur dans l'expression de *f* pour exprimer le déplacement en fonction de la vitesse on obtient :

$$
\begin{equation}
\begin{gathered}
\{R\}_n^{i+1} = \operatorname{proj}_{\mathbb{R}^-} \left( \{R\}_n^{i+1} - \rho_n \left( \{\dot{U}\}_n^{i+1} + \frac{\{\tilde{B}\}_n^i}{\tau} \right) \right), \quad \forall \rho_n > 0 \qquad (8) \\
\text{avec, } \{\tilde{B}\}_n^i = (\{U\}_n^i - \{G_0\}) + \tau(1-\theta) \{\dot{U}\}_n^i
\end{gathered}
\end{equation}
$$

Plutôt que d’utiliser le même intégrateur que pour le problème global, une version modifiée de la
θ–méthode est proposée. Ce schéma, décrit dans [7], permet l’obtention de chocs inélastiques entre
deux pas de temps. Il entraîne l’expression du jeu équivalent {B̃} donnée dans (8).

La convergence de cet algorithme dépend du choix des paramètres d'augmentation : il s'agit d'assurer
la contractance de la fonction *f*. Ici le choix fait est de prendre une approximation diagonale de la matrice
de l'intégrateur [*J*] réduite aux degrés de liberté de contact. En remplaçant la matrice de l'intégrateur par
celle de raideur et en imposant un contact glissant, on obtient un algorithme permettant de calculer
l'équilibre glissant quasi-statique de la structure.

2.2 Étude de stabilité

Si les conditions aux limites sont constantes dans le temps, on peut définir une solution quasi statique ($\mathbf{u}^e, \mathbf{r}^e$) qui vérifie :

$$
\left\{
\begin{array}{l}
\displaystyle \int_{\Omega} \dot{\varepsilon}^* \mathcal{A} \varepsilon^e d\Omega = \int_{\partial\Omega_f} \dot{\mathbf{u}}^* \cdot \mathbf{f}_s dS + \int_{\Omega} \dot{\mathbf{u}}^* \cdot \mathbf{f}_v d\Omega + \int_{S_c} \dot{\mathbf{u}}^* \cdot \mathbf{r}^e dS, \forall \dot{\mathbf{u}}^* \in U^0 \\
u_n^e \le 0 \\
r_n^e \le 0 \text{ et } \mathbf{r}_t^e = -\mu r_n^e \frac{\mathbf{v}_g}{\|\mathbf{v}_g\|} = -\mu r_n^e \mathbf{t}^e \\
u_n^e r_n^e = 0 \\
\varepsilon^e \text{ et } \varepsilon^* \text{ sont les gradients symétriques de } \mathbf{u}^e \text{ et } \mathbf{u}^*
\end{array}
\right.
\tag{9}
$$

La partie de $S_c$ où il y a effectivement contact est appelée zone de contact effective $S_c^e$.

L'étude de la stabilité de cet équilibre se fait en linéarisant les équations autour de cette solution.
Compte tenu des non-régularités induites par le contact unilatéral, on restreint la stabilité aux perturba-
tions préservant l'état de contact. Ces perturbations sont appelées perturbations glissantes et n'admettent
donc pas de déplacement relatif en zone de contact effective.

On définit une base locale sur $S_c^e$, $(\mathbf{n}, \mathbf{t}, \mathbf{b})$, où $\mathbf{n}$ est la normale sortante, $\mathbf{t}$ est la direction de glissement, $\mathbf{t} = \frac{\mathbf{v}_g}{\|\mathbf{v}_g\|}$ et $\mathbf{b} = \mathbf{n} \wedge \mathbf{t}$.

En cherchant les solutions du problème de dynamique sous la forme $(\mathbf{u}, \mathbf{r}) = (\mathbf{u}^e + \hat{\mathbf{u}}, \mathbf{r}^e + \hat{\mathbf{r}})$, où $(\hat{\mathbf{u}}, \hat{\mathbf{r}})$ est une perturbation glissante harmonique, $(\hat{\mathbf{u}}, \hat{\mathbf{r}}) = (\hat{\mathbf{u}}_0, \hat{\mathbf{r}}_0 \exp(\lambda t))$, on obtient le système quadratique aux valeurs propres (QEP) suivant :

$$
\left\{
\begin{array}{l}
\text{trouver } \hat{\mathbf{u}}_0 \in \mathcal{U}, \hat{\mathbf{r}}_0 \in \mathcal{U}^0 \text{ et } \lambda \in \mathbb{R} \text{ tels que } \forall \dot{\mathbf{u}}^* \\
\qquad \lambda^2 \mathcal{M}(\hat{\mathbf{u}}_0, \dot{\mathbf{u}}^*) + \lambda \hat{\mathcal{C}}(\hat{\mathbf{u}}_0, \dot{\mathbf{u}}^*) + \mathcal{K}(\hat{\mathbf{u}}_0, \dot{\mathbf{u}}^*) = \int_{S_c^e} \hat{\mathbf{r}}_{n0} \dot{\mathbf{u}}^* \cdot (\mathbf{n} - \mu_e \mathbf{t}) \, dS \quad (10) \\
\\
\hat{\mathbf{u}}_0 \cdot \mathbf{n} = 0 \text{ on } \varepsilon S_c \\
\\
\hat{\mathcal{C}}(\hat{\mathbf{u}}_0, \mathbf{u}^*) = \mathcal{C}(\hat{\mathbf{u}}_0, \mathbf{u}^*) - \int_{\varepsilon S_c} \mu \frac{e^{\varepsilon r_n}}{\|\mathbf{v}_g\|} \hat{u}_{0b} u_b^* dS
\end{array}
\right.
$$
où $\mathcal{M}, \mathcal{C}, \mathcal{K}$ sont les opérateurs de masse, d'amortissement et de raideur de la structure.

Le terme d'amortissement additionnel qui apparaît dans l'opérateur $\hat{C}$ est lié à la linéarisation de la direction de frottement et a été pris en compte dans l'étude de stabilité pour la première fois dans [11].

On peut reformuler le QEP en déplacement uniquement :

$$
\left\{
\begin{array}{l}
\text{trouver } \hat{\mathbf{u}}_0 \in \mathcal{W}_0 \text{ et, } \lambda \in \mathbb{C} \text{ tels que } \forall \dot{\mathbf{u}}^* \in \mathcal{W}^{\mu} \\
\lambda^2 \mathcal{M}(\hat{\mathbf{u}}_0, \dot{\mathbf{u}}^*) + \lambda \hat{C}(\hat{\mathbf{u}}_0, \dot{\mathbf{u}}^*) + \mathcal{K}(\hat{\mathbf{u}}_0, \dot{\mathbf{u}}^*) = 0 \qquad (11) \\
\mathcal{W}^0 = \{\mathbf{u}, \mathbf{u} \cdot \mathbf{n} = 0 \text{ on } eS_c\} \\
\mathcal{W}^{\mu} = \{\mathbf{u}, \mathbf{u} \cdot (\mathbf{n} - \mu\mathbf{t}) = 0 \text{ on } eS_c\}
\end{array}
\right.
$$

En utilisant une discrétisation élément fini du problème, on obtient :

$$
\left\{
\begin{array}{l}
(\lambda^2 [\tilde{M}] + \lambda [\tilde{C}] + [\tilde{K}]) \{U_0\} = \{0\} \\
([\tilde{M}], [\tilde{C}], [\tilde{K}]) = [T^\mu]^T ([M], [\hat{C}], [K]) [T^N] \\
[\hat{C}] = [C] - \mu [P_B]^T [D] [P_B] \\
D_{ii} = \frac{(e r_n)_i}{(v_g)_i}
\end{array}
\right.
\qquad (12)
$$
où $[T^N]$ et $[T^\mu]$ sont respectivement des bases des espaces $\mathcal{W}^0$ et $\mathcal{W}^\mu$

où $\mathbf{v}_g$ est la norme de la vitesse de glissement à chaque nœud et $[P_B]$ la matrice de projection sur la
direction $\mathbf{b}$.

Le QEP peut être linéarisé en GEP comme suit :

$$
\left\{
\begin{array}{l}
(\Lambda [\tilde{A}] + [\tilde{B}]) \{\Phi_0\} = \{0\} \\
[\tilde{A}] = \begin{bmatrix} [\tilde{C}] & [\tilde{M}] \\ [\tilde{M}] & [0] \end{bmatrix} \quad \text{et} \quad [\tilde{B}] = \begin{bmatrix} [\tilde{K}] & [0] \\ [0] & -[\tilde{M}] \end{bmatrix}
\end{array}
\right.
\qquad (13)
$$

où {$\Phi_0$} est le mode propre dans l'espace d'état, c'est à dire avec la première moitié des lignes corres-
pondant au mode dans l'espace physique, et la seconde au mode multiplié par la valeur propre.

La résolution de ce système est faite en utilisant la méthode d'itération sur les résidus [1]. Comme il s'agit d'un problème non symétrique, les vecteurs propres peuvent être complexes et les valeurs propres peuvent avoir une partie réelle positive. Les modes associés à des valeurs propres à partie réelle positive sont appelés modes instables. On définit le taux de divergence d'un mode associé à la valeur propre $\lambda$ par $2\frac{\Re(\lambda)}{\Im(\lambda)}$, ce taux correspond à un taux d'amortissement négatif.

Les perturbations de l'équilibre se développent suivant ces modes. La croissance des vibrations ob-
tenues est limitée par des effets non-linéaires. Ici, la seule non-linéarité prise en compte étant le contact
frottant, la limitation se fera soit par décollement, soit par adhérence en zone de contact.

3 Démarche de l'étude

3.1 Modèle étudié et stabilité

Le modèle étudié est présenté sur la figure 1(a), il comporte un peu plus de 30000 degrés de liberté.
Il s'agit d'un modèle simplifié d'un frein à disque TGV. La surface supérieure (resp. inférieure) des plots
supérieurs (resp. inférieurs) est encastrée. L'amortissement dans le modèle est pris de type Rayleigh
par morceaux ([C] = α[M] + β[K]) dont les coefficients sont précisés dans le tableau 1. Les paramètres
matériau utilisés sont rappelés dans le tableau 1, ils ont été choisis de manière à être proches de propriétés
moyennes identifiées sur différentes garnitures de frein de TGV. L'équilibre statique de ce modèle ainsi
que sa stabilité sont calculés pour un coefficient de frottement μ de 0.38 et une vitesse de rotation ω de
4 rad s⁻¹. La réponse temporelle non-linéaire du système complet et réduit sont calculées et comparées.
Le choix fait ici est d'utiliser les modes complexes pour la réduction plutôt que les modes réels. Ce choix
est motivé à la fois par le fait que les modes complexes sont ceux qui correspondent à la dynamique du
système autour de l'équilibre, mais aussi par le fait que les N premiers modes réels ne suffisent pas à
correctement représenter les N premiers modes complexes calculés (le cardinal après orthonormalisation
<table><caption>TABLE 1 – Paramètres matériau</caption><thead><tr><th>Paramètre</th><th>disque</th><th>plots</th></tr></thead><tbody><tr><td>Module de Young E (MPa)</td><td>2.02 10<sup>9</sup></td><td>3 10<sup>6</sup></td></tr><tr><td>Coefficient de Poisson v</td><td>0.29</td><td>0.3</td></tr><tr><td>Masse volumique ρ (kg.m<sup>-3</sup>)</td><td>7850</td><td>5102</td></tr><tr><td>Rayleigh α (s<sup>-1</sup>)</td><td>3.75</td><td>54.375</td></tr><tr><td>Rayleigh β (s)</td><td>5 10<sup>-8</sup></td><td>72.5 10<sup>-8</sup></td></tr></tbody></table>

par rapport à la matrice de masse de la base modale composée des parties réelles et imaginaires des 200 premiers modes complexes est de 344).

<figure><img src="image_2.png" /><figcaption>FIGURE 1 – Étude d'un frein TGV simplifié</figcaption></figure>

La stabilité du système est donnée sur la figure 1(b), il comporte 13 modes instables. On peut voir que les modes instables apparaissent à proximité des modes de plots (autour de 3500 Hz et 11000 Hz). La plus forte instabilité se situe autour de 11,1 kHz avec un taux de divergence à 5% en 12 modes instables. Dans les modes basse fréquence, seul le mode à 3537 Hz est instable avec un taux de 0.1%.

## 3.2 Analyse des résultats par ondelettes

Pour contribuer à la compréhension des phénomènes en jeu dans les systèmes déstabilisés par le frottement, il peut être crucial d'avoir une vision sur le contenu temps-fréquence des champs de vitesses et de déplacement étudiés. Dans ce travail, on utilise à cette fin la transformée en ondelettes continue (CWT). Cette dernière permet de s'affranchir en partie du compromis entre résolution en temps et en fréquence que l'on doit faire avec une transformée de Fourier court-terme. Cette transformée introduite avec les travaux de Morlet [12] peut être décrite succinctement comme suit.

Premièrement, le choix d'une ondelette mère, $\phi(t)$ qui est une fonction à décroissance rapide et de moyenne nulle. Une famille de fonctions analysante est générée à partir de l'ondelette mère par translation et dilatation suivant l'équation (14) où $a$ est le facteur d'échelle – ou de dilatation – et $b$ est le facteur de translation temporel. La multiplication par $\frac{1}{\sqrt{a}}$ garantie une normalisation en énergie à toute les échelles.

La CWT est alors effectuée par projection du signal sur la famille d'ondelettes. On obtient alors les coefficients donnés par (15) où $*$ marque la conjugaison des complexes.

$$ \phi_{a,b}(t) = \frac{1}{\sqrt{a}} \phi\left(\frac{t-b}{a}\right) \qquad (14) $$

$$ \mathrm{CWT}_{\phi}^{s}(a, b) = \int_{-\infty}^{\infty} s(t) \phi_{a,b}^{*}(t) dt \qquad (15) $$

Ces coefficients d'ondelette représentent une mesure de la corrélation du signal $s(t)$ avec la famille analysante. Ainsi, $\mathrm{CWT}_{\phi}^{s}(a, b)$ correspond à l'information locale de $s(t)$ à l'échelle $a$ au temps $b$. En pratique, on représente la transformée sous forme de diagramme temps-échelle en puissance défini par $|\mathrm{CWT}_{\phi}^{s}(f, b)|^2$. En utilisant une formule de correspondance échelle-fréquence, on peut obtenir un temps-fréquence. Ici, les ondelettes de Morelt on été choisies, l'ondelette mère est donnée par :
<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Temps-fréquence vitesse relative normale sur le frein TGV.</figcaption></figure>

$$ \phi(t) = e^{jk_0t} e^{-(1/2)t^2} \xrightarrow{FT} \hat{\phi}(\omega) = \sqrt{2\pi}e^{-(1/2)(\omega-k_0)^2} \quad (16) $$

## 3.3 Réponse temporelle sur base réduite

Des études précédentes ont montré [9] que sur des systèmes simplifiés les modes instables étaient prépondérants dans la réponse temporelle, et que l'état dynamique stationnaire non-linéaire est généralement composé de seulement quelques modes complexes. Dans cette étude on se focalisera donc sur l'impact de la réduction avec une troncature modale utilisant les modes complexes plutôt que les modes réels du système. On note $FC_N$ la base de réduction, où N est le nombre de modes inclus. Le critère de troncature ici est de 1.5 fois la fréquence de la sollicitation a été retenu, en prenant comme fréquence de sollicitation la plus haute fréquence des modes instables. Ce choix entraîne la sélection des 461 premiers modes complexes dans la base. La base est enrichie avec l'équilibre statique et les modes statiques de liaison de l'interface. Pour garantir une bonne gestion du contact, l'interface physique n'est pas réduite et les degrés de liberté généralisés sont restreints à la trace des modes en dehors de la zone de contact. Pour réduire les temps de simulation, on met une condition initiale suivant un des modes instables. Ceci permet de sortir plus rapidement le système de l'équilibre statique initial.

### 3.3.1 Résultats sur le frein

La base réduite utilisée permet de réduire d'un facteur 4 le temps de calcul par rapport à la référence. La figure 2 montre une analyse temps-fréquence pour la composante normale au contact de la vitesse relative d'un point. On peut aussi constater que les fréquences apparaissant dans le calcul temporel ne sont pas exactement celles des modes instables, la contribution basse fréquence se situant à 3.42 kHz au lieu de 3.53 kHz et celle haute fréquence à 10.96 kHz au lieu de 11.1 kHz. Cet écart n'est pas surprenant étant donné que les modes sont calculés à l'équilibre en linéarisant les variations de direction de frottement et qu'il s'agit d'une analyse sur un degré de liberté de la zone de contact uniquement. De plus la condition initiale choisie est telle qu'un décollement se produise à $t = 10^{-6}$ s. Ces deux considérations peuvent expliquer les différences observées.

Le temps-fréquence de la base $FC_{461}$ semble être très proche de la référence. On peut noter une légère avance dans l'arrivée des fréquences harmoniques.

La figure 3 montre une comparaison des champs cinématiques moyennés sur la zone de contact pour les différentes bases. Là encore, on peut constater que la base $FC_{461}$ fournit des résultats corrects par rapport à la référence tant en fréquence qu'en amplitude. On peut ici voir une avance de l'ordre de 2.5 ms induit par la réduction sur cette base. Ce phénomène d'avance peut être vu comme une surestimation du taux de croissance de l'instabilité par la base réduite, comme le laisse suggérer l'étude faite sur un modèle 2D [10].

La figure 4 montre une comparaison des spectres sur les dernières 2.5 ms. On pourra remarquer que pour les deux calculs, la fréquence principale autour de 11 kHz est en fait à 10.8 kHz, ce qui correspond à l'harmonique d'ordre 3 de la basse fréquence. On peut encore une fois aussi voir le très bon accord entre la référence et le calcul réduit. Une surestimation des niveaux des pics à $f_1 = 10.8$ kHz, $2f_1$ et $3f_1$ est visible, l'ordre 2 étant plus important que le fondamental. Pour les spectres en puissance, la valeur
<figure><img src="image_4.png" /><figcaption>FIGURE 3 – Comparaison des champs de vitesse et de déplacement moyens en zone de contact</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>FIGURE 4 – FFT en puissance cinétique totale sur le frein TGV.</figcaption></figure>

de la basse fréquence émergente est plus proche du mode instable que celle trouvée dans les évolutions temps-fréquence sur un degré de liberté en zone de contact. En revanche le décalage pour la contribution haute fréquence reste assez marqué.

Finalement, la figure 5 présente la projection des composantes principales de l'état final obtenues par une décomposition en valeurs singulières et normées en masse sur les modes complexes, eux aussi normés en masse. On peut remarquer que les contributions des modes basse fréquence sont moindres dans le calcul de référence que dans le calcul réduit. De plus on notera que la deuxième déformée principale du calcul de référence se projette faiblement sur la base des modes complexes, ce qui veut dire qu'on manque une partie non-négligeable de la solution dans le calcul réduit. Les modes contribuant sont ceux fréquentielle ment proche des modes instables : du numéro 170 à 216 avec moins de 5% d'écart par rapport à la fréquence des modes instables autour de 11 kHz, du numéro 20 à 34 avec moins de 2% d'écart à la fréquence du mode instable à 3.5 kHz ou proche de leurs harmoniques d'ordre 2 (les deux

<figure><img src="image_6.png" /><figcaption>FIGURE 5 – Projection des déformées principales de l'état final sur les modes complexes</figcaption></figure>
plages de modes 107-130 et 144-160 sont à moins de 5% de fréquences comprises dans les harmoniques et combinaisons des fréquences des modes instables.)

# 4 Conclusions

Dans ce travail une démarche de simulation en temporel sur base réduite de systèmes déstabilisés par le frottement a été présentée. L'application d'une stratégie de réduction à un modèle simplifié de frein TGV montre que la base de réduction contenant les 461 premiers modes complexes arrive à capter les bons comportements. Des travaux sont en cours pour tester d'autres bases de réduction et mener les études sur des intervalles de temps plus long afin de tester si l'état final obtenu ici est un état stationnaire ou non. L'analyse des temps-fréquences de référence semblent montrer une légère croissance de la contribution basse fréquence. D'autres bases de réduction sont aussi envisagées, inspirées par la figure 5. L'étude des variations des bassins attracteurs entraînée par la réduction, phénomène observé sur un modèle 2D, sera aussi prochainement menée.

# Références

[1] Adrien Bobillot and Etienne Balmes. Iterative computation of modal sensitivities. *American Institute of Aeronautics and Astronautics Journal*, 44 :1332–1338, 2006.

[2] Fabrice Chevillot, Jean-Jacques Sinou, Nicolas Hardouin, and Louis Jézéquel. Effects of damping on the speed of increase and amplitude of limit cycle for an aircraft braking system subjected to mode-coupling instability. *Archive of Applied Mechanics*, 80 :1045–1054, 2010. 10.1007/s00419-009-0352-8.

[3] Nicolas Coudeyras, Samuel Nacivet, and Jean-Jacques Sinou. Periodic and quasi-periodic solutions for multi-instabilities involved in brake squeal. *Journal of Sound and Vibration*, 328(4-5):520–540, 2009.

[4] Guillaume Vermot des Roches. *Frequency and time simulation of squeal instabilities Application to the design of industrial automotive brakes*. PhD thesis, École Centrale Paris, 2010.

[5] Guillaume Fritz, Jean-Jacques Sinou, Jean-Marc Duffal, and Louis Jézéquel. Effects of damping on brake squeal coalescence patterns - application on a finite element model. *Mechanics Research Communications*, 34(2):181–190, 2007.

[6] R. Jarvis and B. Mills. Vibration induced by dry friction. *Proceedings of the Institution of Mechanical Engineers*, pages 847–866, 1964.

[7] Michel Jean. The non-smooth contact dynamics method. *Computer Methods in Applied Mechanics and Engineering*, 177(3-4):235–257, 1999.

[8] Xavier Lorang. Instabilité vibratoire des structures en contact frottant : Application au crissement des freins de TGV. PhD thesis, École Polytechnique, 2007.

[9] Xavier Lorang and Olivier Chiello. Stability and transient analysis in the modelling of railway disc brake squeal. In *Noise and Vibration Mitigation for Rail Transportation Systems*, volume 99 of *Notes on Numerical Fluid Mechanics and Multidisciplinary Design*, pages 447–453. Springer Berlin / Heidelberg, 2008.

[10] Andréa Loyer, Jean-Jacques Sinou, Olivier Chiello, and Xavier Lorang. Model reduction strategy for railway disc brake squeal prediction. In *6th European Conference on Braking*, 2010.

[11] Frank Moirot. Étude de la stabilité d'un équilibre de Coulomb - Application au crissement des freins à disque (Stability study of a Coulomb equilibrium - Implementation on disc brake squeal). PhD thesis, École Polytechnique, France, 1998.

[12] J. Morlet, G. Arens, E. Fourgeau, and D. Giard. Wave propagation and sampling theory—part ii : Sampling theory and complex waves. *Geophysics*, 47(2):222–236, 1982.

[13] Huajiang Ouyang, John E. Mottershead, Matthew P. Cartmell, and Michael I. Friswell. Friction-induced parametric resonances in discs : Effect of a negative friction-velocity relationship. *Journal of Sound and Vibration*, 209(2):251–264, 1998.

[14] Robb T. Spurr. A theory of brake squeal. *Proceedings of the Institution of Mechanical Engineers*, 1:33–40, 1962.