# Etude de stabilité d'un solide sous frottement de glissement par optimisation de la surface de contact

Xavier Lorang, Marc Bonnet, Qoc Son Nguyen

► To cite this version:

Xavier Lorang, Marc Bonnet, Qoc Son Nguyen. Etude de stabilité d'un solide sous frottement de glissement par optimisation de la surface de contact. 7e colloque national en calcul des structures, CSMA, May 2005, Giens, France. hal-01813052

HAL Id: hal-01813052

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

Submitted on 12 Jun 2018

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 CC0 1.0 - Universal - International License
# Etude de stabilité d'un solide sous frottement de glissement par optimisation de la surface de contact.

Xavier Lorang, Marc Bonnet, Qoc Son Nguyen

Laboratoire de Mécanique des Solides - CNRS UMR 7649
Ecole Polytechnique, 91128 Palaiseau cedex, France
lorang@lms.polytechnique.fr, bonnet@lms.polytechnique.fr,
son@lms.polytechnique.fr

RÉSUMÉ. On s'intéresse à la dynamique d'un disque élastique en glissement stationnaire par rapport à un solide rigide. Il a été démontré pour le disque en rotation uniforme et pour un coefficient de frottement donné que l'équilibre glissant est instable. L'étude numérique de la stabilité se traduit par la résolution d'un problème aux valeurs propres généralisées. L'objet de cette communication est de présenter une méthode permettant d'évaluer la sensibilité de la valeur propre par rapport à la géométrie de la surface de contact et de définir un algorithme permettant de trouver la surface de contact optimale afin de stabiliser un mode précis. La méthode de calcul de la sensibilité est illustrée sur un exemple simple afin de la valider. Ce travail s'inscrit à long terme dans une optique de réduction du phénomène de crissement des freins à disque.

**ABSTRACT.** This communication deals with the steady-state sliding of a rotating elastic disk on a rigid body. Previous investigations have shown that, if the friction coefficient is constant, the sliding equilibrium is unstable. The numerical study of stability is based on a generalized eigenvalue problem. This work aims at adjusting the shape of the contact surface so as to preclude the occurrence of the unstable eigenmodes. Towards this goal, a methodology for computing the sensitivity of a given eigenvalue of that problem with respect to the shape of the contact surface, based on Lagrangian differentiation with respect to shape parameters, is proposed and validated on a numerical example. This work is expected to ultimately help reducing disk brake squeal via an optimal design of the contact surface.

**MOTS-CLÉS :** Stabilité dynamique, frottement, optimisation.

**KEYWORDS:** Dynamical stability, friction, optimization.
# 1. Contexte et hypothèses

On considère un disque élastique $D$ en rotation autour de son axe de révolution. Ce disque est en contact frottant avec un solide indéformable fixe sur une surface $S_c$ (frottement de Coulomb et coefficient de frottement $f$ constant). On s'intéresse ici à l'évolution d'une perturbation lorsque l'on est proche de l'équilibre stationnaire. Les équations du contact unilatéral frottant sont telles que l'on ne peut se permettre de considérer une perturbation quelconque, le caractère non différentiable de ses équations rendant la linéarisation impossible. On choisit donc une perturbation qui conserve la condition de non décollement et on ne tient pas compte du phénomène d'adhérence. Dans ces conditions, la perturbation $\hat{u}$, dite glissante, du déplacement $u$ est définie par $u = u_e + \hat{u}$ avec $\hat{u}_n = 0$ sur $S_c^e$ ($u_e$ déplacement à l'équilibre, $S_c^e$ zone de contact à l'équilibre). Le contact frottant devient bilatéral avec glissement imposé : $T = fNw/|w|$ avec $w$ la vitesse de glissement et $N$ la réaction normale de contact. On choisit d'étudier la stabilité de l'équilibre par la méthode de linéarisation en négligeant l'effet du frottement plan (étudié numériquement). On a alors $\hat{T} = f\hat{N}t$ avec $t = v/|v|$, $v$ étant la vitesse du disque à l'équilibre (orthoradiale ici).

# 2. Equation de la dynamique du disque et instabilité glissante

En écrivant les équations de la dynamique du disque et en négligeant les termes convectifs [MOIR00], $(\hat{u}, \hat{N})$ est régi par le problème variationnel suivante :

$$ \left\{ \begin{array}{l} \displaystyle \int_D \rho \hat{u}_{,tt} \cdot u^* d\Omega + \int_D \varepsilon(\hat{u}) : L : \varepsilon(u^*) d\Omega - \int_{S_c} \hat{N}(n+ft) \cdot u^* dS = 0, \forall u^* \\[1em] \displaystyle \int_{S_c} \hat{u}_n N^* dS = 0, \forall N^* \end{array} \right. \qquad (1) $$

Poser $\hat{u}(x, t) = \hat{u}_0(x) \exp(\lambda t)$ et $\hat{N}(x, t) = \hat{N}_0(x) \exp(\lambda t)$, conduit alors au problème aux valeurs propres généralisé :

$$ \left\{ \begin{array}{l} \displaystyle \lambda^2 \int_D \rho \hat{u}_0 \cdot u^* d\Omega + \int_D \varepsilon(\hat{u}_0) : L : \varepsilon(u^*) d\Omega - \int_{S_c} \hat{N}_0(n+ft) \cdot u^* dS = 0, \forall u^* \\[1em] \displaystyle \int_{S_c} \hat{u}_{0n} N^* dS = 0, \forall N^* \end{array} \right. \qquad (2) $$

qui gouverne les modes d'instabilité dynamique du régime glissant. Une écriture plus compacte du problème [2] est :

$$ \lambda^2 M(u, u^*) + K(u, u^*) + L(N, u_n^*) + L(u_n, N^*) + F(N, u_t^*) = 0, \forall (u^*, N^*) \quad (3) $$

ayant utilisé les notations simplifiées $u$ et $N$ pour $\hat{u}_0$ et $\hat{N}_0$, et avec

$$ \begin{aligned} M(u, u^*) &= \int_D \rho u \cdot u^* d\Omega & K(u, u^*) &= \int_D \nabla u : L : \nabla u^* d\Omega \\ L(N, u_n^*) &= - \int_{S_c} N u_n^* dS & F(N, u_t^*) &= - \int_{S_c} f N u_t^* dS \end{aligned} $$
### 3. Problème d'optimisation

L'évolution de la structure suivant tout mode solution de [3] est

$$ \mathrm{Re}(u(x)e^{\lambda t}) = e^{\mathrm{Re}(\lambda)t}[\mathrm{Re}(u(x)) \cos(\mathrm{Im}(\lambda)t) - \mathrm{Im}(u(x)) \sin(\mathrm{Im}(\lambda)t)] $$

On note que si $\mathrm{Re}(\lambda) > 0$, alors le mode est instable. L'objectif poursuivi ici est d'éviter la déstabilisation d'un mode choisi en modifiant la géométrie de la zone de contact $S_c$. La stratégie de calcul est de déterminer cette zone par un algorithme de minimisation de la fonction coût

$$ \mathcal{J}(S_c) = 2\mathrm{Re}(\lambda_i(S_c)) = \lambda_i(S_c) + \bar{\lambda}_i(S_c) $$

i étant fixé. Pour cela, la surface plane ajustable $S_c$ est décrite à l'aide d'une transformation géométrique $\Phi$ telle que $x = \Phi(X, p)$ dépendant de paramètres géométriques $p$. Il s'agit de chercher $p$, et donc $S_c$, qui minimise la fonction coût $J(p) = \mathcal{J}(S_c(p))$. Pour gagner en efficacité, on se place dans l'optique de l'utilisation d'une méthode d'optimisation employant le gradient de $J(p)$. Le développement qui suit a pour objet de proposer une procédure d'évaluation numériquement peu coûteuse des dérivées $J_{,p_k} = \lambda_{i,p_k} + \bar{\lambda}_{i,p_k}$ (dans la suite, on utilisera la notation abrégée ($\cdot$) pour $(\cdot)_{,p_k}$). Afin de calculer $\dot{\lambda}_i$, on effectue une dérivation par rapport à $p_k$ de l'équation [3] :

$$ (\lambda^2) \dot{M}(u, u^*) + \lambda^2 \dot{M}(u, u^*) + \dot{K}(u, u^*) + \dot{L}(N, u_n^*) + \dot{L}(u_n, N^*) + \dot{F}(N, u_t^*) = 0 \quad (4) $$

Le calcul de $\dot{M}$, $\dot{K}$, $\dot{L}$ et $\dot{F}$ necessite la dérivation d'intégrales de volume et de surface. Or la singularité de certains intégrands sur le bord de la surface de contact conduit à utiliser des formules de dérivation lagrangienne, en analogie avec les procédures de calcul de taux de restitution d'énergie de type « $G - \theta$ » en mécanique linéaire de la rupture. On est donc amené à introduire, compte tenu du caractère bilinéaire par rapport aux champs des composantes du problème variationnel [3], les expressions suivantes :

$$ \begin{aligned} \dot{M}(u, u^*) &= M(\dot{u}, u^*) + M(u, \dot{u}^*) + M'(u, u^*) \\ \dot{K}(u, u^*) &= K(\dot{u}, u^*) + K(u, \dot{u}^*) + K'(u, u^*) \\ \dot{L}(N, u_n^*) &= L(\dot{N}, u_n^*) + L(N, \dot{u}_n^*) + L'(N, u_n^*) \\ \dot{F}(N, u_t^*) &= F(\dot{N}, u_t^*) + F(N, \dot{u}_t^*) + F'(N, u_t^*) \end{aligned} $$

dans lesquelles $M'(u, u^*)$, $K'(u, u^*)$, $L'(N, u_n^*)$ et $F'(N, u_t^*)$ correspondent aux valeurs prises par les dérivées lagrangienes des intégrales correspondantes à dérivée lagrangienes des champs nulles :

$$ \begin{aligned} M'(u, u^*) &= \int_D \rho u \cdot u^* \operatorname{div} \theta_k d\Omega \\ K'(u, u^*) &= \int_D (\nabla u : L : \nabla u^*) \operatorname{div} \theta_k - (\nabla u : \nabla \theta_k) : L : \nabla u^* - \nabla u : L : (\nabla u^* \nabla \theta_k) d\Omega \end{aligned} $$
$$L'(N, u_n^*) = - \int_{S_c} N u_n^* \mathrm{div}_S \theta_k dS$$

$$F'(N, u_t^*) = - \int_{S_c} f N u_t^* \mathrm{div}_S \theta_k dS.$$

$\theta_k(X, p) = \Phi_{p_k}(X, p)$ désignant le champ de vitesse de transformation associé à $\Phi$. On notera que $\rho$ est constant et que $\dot{t} = \dot{n} = 0$. En remarquant que $(\dot{u}^*, \dot{N}^*)$ est un champ test admissible, on obtient :

$$(\lambda^2)' M(u, u^*) + \lambda^2 (M(\dot{u}, u^*) + M'(u, u^*)) + K(\dot{u}, u^*) + K'(u, u^*) + L(\dot{N}, u_n^*) \\ + L'(N, u_n^*) + L(\dot{u}_n, N^*) + L'(u_n, N^*) + F(\dot{N}, u_t^*) + F'(N, u_t^*) = 0 \quad (5)$$

Afin d'éliminer les termes contenant $(\dot{u}, \dot{N})$ de [5], on considère le problème aux valeurs propres suivant, problème adjoint au problème [3] :

$$\lambda^2 M(v, v^*) + K(v, v^*) + L(v_n, Q^*) + L(Q, v_n^*) + F(v_t, Q^*) = 0, \forall (v^*, Q^*) \quad (6)$$

Les choix de fonction-test $(u^*, N^*) = (v, Q)$ dans [5] et $(v^*, Q^*) = (\dot{u}, \dot{N})$ dans [6] permettent, par combinaison des deux égalités obtenues, d'éliminer $(\dot{u}, \dot{N})$ et d'obtenir directement l'expression de la dérivée de la valeur propre :

$$\dot{\lambda} = - \frac{\lambda^2 M'(u, v) + K'(u, v) + L'(N, v_n) + F'(N, v_t)}{2\lambda M(u, v)} \quad (7)$$

L'évaluation pratique de [7] a pu être simplifiée en remarquant que, si $G_\theta(u)$ désigne le taux de restitution de l'énergie en régime quasistatique pour le champ $u$ associé à une fissure fictive dont le front est le bord de la surface de contact, alors

$$K'(u, v) = -2G_\theta(u)$$

Compte tenu par ailleurs du caractère bilinéaire en $(u, v)$ de $K'(u, v)$, on obtient

$$K'(u, v) = \frac{1}{2} (G_{\theta}(u - v) - G_{\theta}(u + v))$$

L'intérêt pratique de cette remarque réside dans la possibilité d'évaluer la contribution $K'(u, v)$ à [7] au moyen d'un opérateur de calcul de taux de restitution d'énergie, une possibilité offerte par exemple par CAST3M. Un traitement similaire n'étant pas disponible pour $M'(u, v)$, le calcul de ce terme a été programmé directement. La procédure mise en œuvre consiste alors à calculer les modes $u$ et $v$ dans MATLAB, puis évaluer $\dot{\lambda}$ par [7] au moyen d'un post-traitement des champs $(u, v)$, cette fois dans l'environnement CAST3M afin d'utiliser l'opérateur donnant $G_\theta$.

Ainsi, le calcul de $\lambda_{i,p_k}$ pour tous les paramètres $p_k$ permet d'évaluer le gradient $\nabla_p J(p)$ de $J(p)$. L'optimisation de $J(p)$ sera réalisée à l'aide d'une procédure du optimization toolbox de MATLAB, à laquelle l'information $\nabla_p J(p)$ est transmise.
<figure><img src="image_2.png" /><figcaption>**Figure 1.** *Configuration géométrique de l'exemple.*</figcaption></figure>

**4. Exemple**

La méthode n'ayant pas encore été appliquée au mécanisme industriel, un exemple simple de validation est présenté ici. On considère, sous l'hypothèse des déformations planes, un solide élastique rectangulaire en contact frottant glissant sur $S_c = S_{c1} \cup S_{c2}$ (Figure 1), avec pour conditions aux limites $u|_{S_1} = u|_{S_2}$ (condition de périodicité), $u|_{S_3} = 0$. Les caractéristiques du matériau sont : $E = 2.5 \times 10^{11}$ Pa, $\nu = 0.3$, $\rho = 7800$ kg/m$^3$ et $f = 1$. La surface de contact est paramétrée par le seul paramètre $\alpha$. Le premier mode instable est le mode 3.

Afin de valider la démarche, on calcule en fonction de $\alpha$ la valeur propre $\lambda_3$ associée au mode 3 et la dérivée $\lambda_3^\alpha$ donnée par [7]. Il est alors possible de vérifier la cohérence du résultat en comparant la dérivée calculée par [7] avec la dérivée numérique de $\lambda(\alpha)$ (Figures 2 à 5). Les résultats font apparaître une bonne concordance entre les dérivées analytique et numérique. On remarque cependant que cette dernière présente des oscillations, vraisemblablement dues au fait que dans la procédure de dérivation numérique les valeurs de $\lambda(\alpha)$ et $\lambda(\alpha + \Delta\alpha)$ sont obtenues sur des maillages différents.

<figure><img src="image_3.png" /><figcaption>**Figure 2.** $Re(\lambda)$ en fonction d'{$\alpha$}</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 3.** $Im(\lambda)$ en fonction d'{$\alpha$}</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>Figure 4. Re(dλ/dα) en fonction d'α</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5. Im(dλ/dα) en fonction d'α</figcaption></figure>

**5. Conclusion et perspectives**

A ce stade de l'étude, une procédure à base de dérivation lagrangienne permettant le calcul de la dérivée d'une valeur propre du problème de perturbation glissante associé à un solide en rotation a été formulée et validé. Les étapes suivantes consistent à intégrer cette procédure à l'algorithme d'optimisation, puis à traiter le cas industriel. Les premiers résultats de l'algorithme d'optimisation proposé seront présentés à partir du modèle numérique élément fini d'un mécanisme simplifié de frein à disque. On montrera le contour optimal de la surface défini par les paramètres p optimaux.

A terme, les résultats de la méthode permettront de suggérer une nouvelle géométrie des garnitures du frein, pour laquelle il faudra alors vérifier par des mesures acoustiques qu'il y a effectivement disparition de la fréquence associée au mode optimisé. Par ailleurs, pour ce qui concerne la simulation numérique, on pourra chercher à vérifier que la géométrie optimale obtenue par l'algorithme proposé ici conduit également à la stabilisation de la solution dynamique du problème non linéaire.

**6. Bibliographie**

[MOIR00] F. MOIROT, Q.S NGUYEN, 2000, Brake squeal : a problem of flutter instability of the steady sliding solution ? *Arch. Mech.* **52**,645-661.

[MOIR98] F. MOIROT, 1998, Etude de la stabilité d'un équilibre en présence du frottement de coulomb, *Thèse, Ecole polytechnique*.

[MAR95] J.A.C MARTINS, J. GUIMARAES, L.O. FARIA, 1995, Dynamic surface solutions in linear elasticity and viscoelasticity with frictional boundary conditions, *J. Vibration and Acoustics*, **117**, 445-451