# Une méthode fréquentielle pour l'identification de structures non-linéaires

Sergio Bellizzi, Robert Bouc, Mario Defilippi

► To cite this version:

Sergio Bellizzi, Robert Bouc, Mario Defilippi. Une méthode fréquentielle pour l'identification de structures non-linéaires. Mechanics & Industry, 2004, 5, pp.61-69. 10.1051/meca:2004007. hal-00088262

HAL Id: hal-00088262

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

Submitted on 29 Sep 2020

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
# Une méthode fréquentielle pour l'identification de structures non-linéaires

SERGIO BELLIZZIª, ROBERT BOUC ET MARIO DEFILIPPI

Laboratoire de Mécanique et d'Acoustique - CNRS, 31 chemin Joseph Aiguier, 13402 Marseille Cedex 20, France

**Résumé** – Cet article décrit une procédure d'identification des composantes non-linéaires d'une structure mécanique à partir de mesures spectrales de la sollicitation et de la réponse. Elle suppose des sollicitations aléatoires stationnaires large bande. La méthode repose sur la minimisation de l'écart quadratique moyen entre la densité spectrale de puissance mesurée et une approximation semi-analytique de cette dernière obtenue par linéarisation stochastique à paramètres aléatoires de l'équation du mouvement. La notion de modes couplés non-linéaires dépendant d'une seule variable scalaire (liée à la norme des amplitudes modales de la réponse non-linéaire) est utilisée pour construire le modèle linéaire équivalent. L'aspect aléatoire de ce modèle est introduit à l'aide de la loi de probabilité de l'énergie modale de la réponse non-linéaire. La méthode est utilisée pour identifier les supports non-linéaires d'une poutre à partir de mesures.

**Mots clés** : Identification paramétrique / modes non-linéaires / linéarisation à paramètres aléatoires / approximation des DSP

**Abstract** – A frequency method to identify non-linear structures. This paper deals with non-linear mechanical structures identification using spectral measurement of the input excitation and of the corresponding structure response. Identification method needs a wide-band, stationary stochastic input process and is based on minimizing the mean-square deviation between the measured power spectral density and a semi analytic approximation of it, given by a stochastic random parameters linearization of the motion equation. The equivalent linear model is derived using non-linear coupled modes depending on a scalar variable alone. Stochasticity of this model is related to the modal energy probability law of the nonlinear response. The method is applied to identify non-linear supports of a beam.

**Key words:** Parametric identification / non-linear modes / random parameters linearization / PSD approximation

## 1 Introduction

La matrice de Densités Spectrales de Puissances (DSP) de la réponse d'un système mécanique consiste une des principales caractéristiques du comportement vibratoire. Elle fournit notamment la distribution de l'énergie par bande de fréquences, d'où le rôle important qu'elle peut jouer dans les procédures d'identification ou de mise à jour de modèles non-linéaires. Des exemples d'identification fondée sur cet aspect fréquentiel sont donnés dans [1, 2].

Pour des systèmes mécaniques non-linéaires faiblement amortis des méthodes fournissant une approximation analytique de la matrice de DSP de la réponse ont été développées ces dix dernières années. Ces méthodes restituent correctement l'énergie totale de la réponse mais

reproduisent aussi les translations des fréquences de résonance ainsi que l'étalement du spectre autour des résonances; translation et étalement sont des signatures des phénomènes non-linéaires. Ces méthodes utilisent la notion de modes non-linéaires [3] et/ou le concept de Système Linéaire à Paramètres Aléatoires (SLPA) [4, 5]. La matrice de DSP de la réponse du SLPA, soumis à la même entrée que celle du système non-linéaire, fournit l'approximation de la matrice de DSP de la réponse du système non-linéaire.

Dans l'approche proposée dans [6, 7], les matrices d'amortissement et de raideur définissant le SLPA sont construites à partir des modes couplés non-linéaires du système conservatif associé. Les modes couplés non-linéaires, fréquences et formes modales, sont fonction du niveau de la réponse par l'intermédiaire des « amplitudes modales ». Ces amplitudes sont considérées comme des

<footnote>ª Auteur correspondant : bellizzi@lma.cnrs-mrs.fr</footnote>
Nomenclature

<table>
  <tr>
    <td>a<sub>i</sub></td>
    <td>amplitude modale du i-ème mode du système projeté</td>
  </tr>
  <tr>
    <td>A<sup>2</sup> = &sum;<sub>i=1</sub><sup>p</sup> a<sub>i</sub><sup>2</sup></td>
    <td>processus « énergie modale »</td>
  </tr>
  <tr>
    <td>F</td>
    <td>non-linéarité en variables physiques</td>
  </tr>
  <tr>
    <td>F</td>
    <td>non-linéarité du système projeté sur les p modes naturels du système linéaire sous-jacent</td>
  </tr>
  <tr>
    <td>[K<sub>eq</sub>]</td>
    <td>matrice de raideur approximation de la non-linéarité par LSG</td>
  </tr>
  <tr>
    <td>J(&theta;)</td>
    <td>fonction coût pour le problème d'identification</td>
  </tr>
  <tr>
    <td>[M], [C], [K]</td>
    <td>matrices de masse, d'amortissement et de raideur du système physique</td>
  </tr>
  <tr>
    <td>[M], [C], [K]</td>
    <td>matrices de masse, d'amortissement et de raideur du système projeté sur les p modes naturels du système linéaire sous-jacent</td>
  </tr>
  <tr>
    <td>p(A)</td>
    <td>densité de probabilité de l'énergie modale A</td>
  </tr>
  <tr>
    <td>q<sub>i</sub></td>
    <td>i-ème composante modale du SLPA</td>
  </tr>
  <tr>
    <td>[S]</td>
    <td>matrice de répartition géométrique de la sollicitation pour le système physique</td>
  </tr>
  <tr>
    <td>[S]</td>
    <td>matrice de répartition géométrique de la sollicitation pour le système projeté</td>
  </tr>
  <tr>
    <td>[S<sub>q</sub>]</td>
    <td>matrice de DSP de q vecteurs des q<sub>i</sub> composantes modales du SLPA</td>
  </tr>
  <tr>
    <td>[S<sub>U</sub>]</td>
    <td>matrice de DSP du déplacement U du système physique</td>
  </tr>
  <tr>
    <td>[S<sub>W</sub>]</td>
    <td>matrice de DSP de la sollicitation W</td>
  </tr>
  <tr>
    <td>[S<sub>X</sub>]</td>
    <td>matrice de DSP du déplacement X du système projeté</td>
  </tr>
  <tr>
    <td>U</td>
    <td>vecteur déplacement du système physique</td>
  </tr>
  <tr>
    <td>W</td>
    <td>sollicitation du système physique</td>
  </tr>
  <tr>
    <td>X</td>
    <td>vecteur déplacement du système projeté suivant [Φ]<sub>p</sub></td>
  </tr>
  <tr>
    <td>&delta;<sub>i</sub></td>
    <td>facteur de participation du mode i à l'énergie modale totale A<sup>2</sup></td>
  </tr>
  <tr>
    <td>&theta;</td>
    <td>vecteur des paramètres du modèle non-linéaire</td>
  </tr>
  <tr>
    <td>&xi;</td>
    <td>constante d'amortissement visqueux du système physique</td>
  </tr>
  <tr>
    <td>&phi;<sub>i</sub></td>
    <td>phase du i-ème mode du système projeté</td>
  </tr>
  <tr>
    <td>[Φ]<sub>p</sub></td>
    <td>matrice modale des p premiers modes naturels du système linéaire sous-jacent au système physique</td>
  </tr>
  <tr>
    <td>&Psi;<sub>eq,i</sub></td>
    <td>forme modale du i-ème mode du SLE obtenu par LSG</td>
  </tr>
  <tr>
    <td>&Omega;<sub>i</sub></td>
    <td>pulsation du i-ème mode (non-linéaire) du système projeté</td>
  </tr>
  <tr>
    <td>&omega;<sub>eq,i</sub></td>
    <td>pulsation du i-ème mode du SLE obtenu par LSG</td>
  </tr>
  <tr>
    <td>[&omega;<sub>i</sub><sup>2</sup>]<sub>p</sub></td>
    <td>matrice diagonale des carrés des p premières pulsations naturelles du système linéaire sous-jacent au système physique</td>
  </tr>
  <tr>
    <td>[.]<sup>T</sup></td>
    <td>matrice transposée</td>
  </tr>
  <tr>
    <td>E</td>
    <td>espérance mathématique</td>
  </tr>
</table>

variables aléatoires dont la loi de probabilité est celle des « amplitudes modales » de la réponse stationnaire du système non-linéaire. Cette loi de probabilité, qui peut, sous certaines hypothèses, être calculée, dépend des modes couplés non-linéaires, de l'amortissement et de la DSP de de la sollicitation. Pour un système à n DDL, le calcul d'une telle approximation de la matrice de DSP nécessite, pour chaque pulsation, l'évaluation numérique d'une intégrale n-uple et la résolution d'un problème aux valeurs propres non-linéaire (calcul des modes couplés non-linéaires) en chaque point d'intégration.

Pour réduire le coût de calcul une version simplifiée de la méthode a été proposée dans [8, 9]. Les fréquences modales sont définies comme fonctions d'une seule variable aléatoire, la norme $L^2$ du vecteur amplitude, tandis que les formes modales sont données, en un seul calcul, par Linéarisation Stochastique Gaussienne (LSG). La détermination de la DSP approchée n'exige alors que le calcul d'une intégrale simple. Les approximations de la DSP ainsi calculées qui restent toujours très sensibles aux composantes non-linéaires, peuvent être avantageusement employées dans une procédure d'identification à

l'aide d'une fonction coût tenant compte de la différence entre l'approximation et la DSP mesurée.

Cet article est consacré à la description de la procédure d'identification. Les étapes d'approximation de la DSP sont rappelées section 2. La section 3 est plus particulièrement dédiée à la construction du SLPA. Le problème d'identification est abordé à la section 4 et une mise en oeuvre à partir de données recueillies sur une structure réelle constituée d'une poutre sur supports non-linéaires est présentée section 5. Signalons que dans la suite de cet exposé, nous utiliserons indifféremment le terme fréquence ou pulsation pour désigner les pulsations.

## 2 Problème direct : approximation de la DSP

Soit un système mécanique à n DDL, décrit par un système différentiel du second ordre de la forme

$$[M] \ddot{U}(t)+2\xi [C] \dot{U}(t)+[K] U(t)+F(U(t)) = \sqrt{4\xi} [S] W(t) \quad (1)$$
$[\mathcal{M}]$, $[\mathcal{C}]$ et $[\mathcal{K}]$ sont des matrices $n \times n$ réelles, symétriques et définies positives, $\xi$ est une constante vérifiant $0 < \xi \ll 1$. La fonction $\mathcal{F}(U)$ qui modélise la force de rappel non-linéaire, est supposée monotone, croissante et impaire. $[S]$ est une matrice $n \times m$ réelle et $W(t)$ un processus stochastique à valeurs dans $\mathbb{R}^m$, stationnaire, large bande, de valeur moyenne nulle (la matrice de DSP est $[S_W(\omega)]$).

La construction d'une approximation de la matrice de DSP $[S_U(\omega)]$ de la réponse stationnaire $U(t)$ utilise les étapes suivantes :

**étape 1.** Réduction du système (1) par projection sur les $p$ ($p < n$) premiers modes du système linéaire sous jacent ($\mathcal{F}(U) = 0$) d'où l'équation

$$[M] \ddot{X}(t) + 2\xi [C] \dot{X}(t) + F(X(t)) = \sqrt{4\xi}[S]W(t) \quad (2)$$

$U = [\Phi]_p X, [M] = [\Phi]_p^T [\mathcal{M}] [\Phi]_p, [C] = [\Phi]_p^T [C] [\Phi]_p, [S] = [\Phi]_p^T [S]$ et $F(X) = [K]X + [\Phi]_p^T \mathcal{F}([\Phi]_p X)$ avec $[K] = [\Phi]_p^T [\mathcal{K}] [\Phi]_p$. $[\Phi]_p$ représente la matrice modale constituée des $p$ premiers modes naturels c'est-à-dire que $[\Phi]_p$ est une matrice $n \times p$ satisfaisant le problème aux valeurs propres classique $[\mathcal{K}] [\Phi]_p = [\mathcal{M}] [\Phi]_p [\omega_i^2]_p$ où $[\omega_i^2]_p$ est la matrice diagonale des carrés des $p$ premières fréquences naturelles;

**étape 2.** Approximation de la matrice de DSP $[S_X(\omega)]$ à l'aide d'un SLPA $[7, 8]$;

**étape 3.** Approximation de la matrice de DSP $[S_U(\omega)]$

$$[S_U(\omega)] = [\Phi]_p [S_X(\omega)][\Phi]_p^T \quad (3)$$

Dans ce qui suit seule l'étape 2 est détaillée.

# 3 Construction du SLPA et approximation de la DSP

## 3.1 Résultats préliminaires. LSG et amplitudes de Rayleigh

Le Système Linéaire Équivalent (SLE) issu de la LSG [10] appliquée à (2) est de la forme

$$[M] \ddot{X}(t) + 2\xi [C] \dot{X}(t) + [K_{eq}]X(t) = \sqrt{4\xi}[S]W(t) \quad (4)$$

où la matrice de raideur $[K_{eq}]$ minimise la moyenne quadratique $\mathbb{E}(\mathcal{E}^T \mathcal{E})$ de l'écart $\mathcal{E} = F(X) - [K_{eq}]X$. La loi de probabilité de la réponse non-linéaire n'étant pas connue, la loi de probabilité (qui est une gaussienne) de la réponse du système linéaire (4) est utilisée pour évaluer l'espérance mathématique $\mathbb{E}(\cdot)$ et un procédé itératif conduit à la matrice de raideur équivalente $[K_{eq}]$.

Au couple déplacement-vitesse $(X, \dot{X})$, solution de (4), est associé le couple amplitude-phase $(a, \varphi)$ par la transformation

$$\left\{ \begin{array}{l} X(t) = \sum_{i=1}^{p} \Psi_{eq,i} a_i(t) \cos(\omega_{eq,i} t + \varphi_i(t)) \\ \dot{X}(t) = - \sum_{i=1}^{p} \Psi_{eq,i} \omega_{eq,i} a_i(t) \sin(\omega_{eq,i} t + \varphi_i(t)) \end{array} \right. \quad (5)$$

où $a_i$ (respectivement $\varphi_i$) représente la $i$-ème composante de l'amplitude (respectivement de la phase) et $(\omega_{eq,i}^2, \Psi_{eq,i})$ désignent les $p$ modes propres du SLE solutions du problème aux valeurs propres

$$[M][\Psi_{eq}][\omega_{eq}^2] = [K_{eq}][\Psi_{eq}] , \quad [\Psi_{eq}]^T [M][\Psi_{eq}] = [I] \quad (6)$$

où $[\Psi_{eq}]$ désigne la matrice des vecteurs propres.

Le principe de moyenne stochastique [11] permet de décrire l'évolution des composantes modales $a_i$, à l'aide de l'Équation Différentielle Stochastique de Itô (EDSI) suivante (pour simplifier l'écriture la dépendance temporelle des variables sera toujours omise),

$$\begin{aligned} da_i = -\xi \left( c_i a_i - \frac{1}{a_i} g_i^T [S_W^R(\omega_{eq,i})] g_i \right) dt \\ + (2\xi g_i^T [S_W^R(\omega_{eq,i})] g_i)^{\frac{1}{2}} dw_i \end{aligned} \quad (7)$$

où $c_i = \Psi_{eq,i}^T [C] \Psi_{eq,i}$, $g_i^T = \omega_{eq,i}^{-1} \Psi_{eq,i}^T [S]$ et $[S_W^R]$ désigne la partie réelle de $[S_W]$. On peut vérifier à partir de (7) que les amplitudes $a_i$ sont distribuées suivant une loi de Rayleigh et stochastiquement indépendantes.

## 3.2 Décomposition modale avec fréquences fonction d'une variable vectorielle

Suivant [8], au vecteur déplacement-vitesse, solution de l'équation (2), est associé un couple amplitude-phase par la transformation

$$X(t) = \sum_{i=1}^{p} \Psi_{eq,i} a_i(t) \cos \phi_i(t) \quad (8)$$

$$\dot{X}(t) = - \sum_{i=1}^{p} \Psi_{eq,i} \Omega_i(a(t)) a_i(t) \sin \phi_i(t) \quad (9)$$

$$\phi_i(t) = \int_0^t \Omega_i(a(s))ds + \varphi_i(t), \quad (10)$$

où $a_i(t)$ (respectivement $\varphi_i(t)$) désigne la $i$-ème composante de $a(t)$ (respectivement $\varphi(t)$), $\phi_i(t)$ désigne l'angle associé et les vecteurs $\Psi_{eq,i}$ satisfont (6). Ce développement est une généralisation de la décomposition modale (5), les fonctions $a \mapsto \Omega_i(a)$ ($\mathbb{R}_+^p \mapsto \mathbb{R}_+$) sont appelées fréquences modales. En ne conservant que les coefficients des termes en $\cos \phi_j$ pour $j = 1, \dots, p$ dans le développement de Fourier de la fonction $2\pi$-multipériodique

$$\phi = (\phi_1, \dots, \phi_p) \mapsto F \left( \sum_{j=1}^{p} \Psi_{eq,j} a_j \cos \phi_j \right) \quad (11)$$
la méthode du premier harmonique basée sur l'équilibrage
harmonique suivant les termes cos φ<sub>j</sub> pour j = 1, ..., p,
appliquée à l'équation (2), fournit les équations

$$
[M] \Psi_{eq,i} \Omega_i^2 = \frac{2}{a_i (2\pi)^p} \int_{[0,2\pi]^p} F \left( \sum_{j=1}^{p} \Psi_{eq,j} a_j \cos \phi_j \right) \\
\times \cos \phi_i \, d\phi_1 \cdots d\phi_p, \quad i = 1, \dots, p \quad (12)
$$

qui définissent, pour $a = (a_1, a_2, \dots, a_p)$ fixé, les
fréquences $\Omega_i$.

En intégrant par parties relativement aux $\phi_i$, les
fréquences $\Omega_i$ s'expriment sous la forme

$$
\Omega_i^2(a) = \frac{2}{(2\pi)^p} \int_{[0,2\pi]^p} \Psi_{eq,i}^T \left[ \frac{\partial}{\partial X} F \left( \sum_{j=1}^p \Psi_{eq,j} a_j \cos \phi_j \right) \right] \times \Psi_{eq,i} \sin^2 \phi_i \, d\phi_1 \cdots d\phi_p \quad (13)
$$

Chaque fréquence est une fonction du vecteur ampli-
tude a. Pour i = 1, ..., p, les couples (Ωᵢ(a), Ψ<sub>eq,i</sub>) se-
ront, par abus de langage, appelés modes couplés non-
linéaires. En effet, dans la littérature, la notion de modes
non-linéaires fait référence à des propriétés intrinsèques
du système, ce qui n'est pas le cas ici. Les Ψ<sub>eq,j</sub>, et par
suite les Ωᵢ, dépendent de l'excitation.

La méthodologie et le principe de moyenne stochas-
tique décrits dans [7], appliqués à (2), conduisent aux
2p-EDSI couplées qui caractérisent la dynamique des
composantes de a et ϕ

$$
\begin{equation}
\begin{split}
d a_i = -\xi \left( \right. & \left( \frac{c_i}{D_i} a_i - g_i^T \frac{\partial [S_W^R(\Omega_i(a))]g_i}{\partial a_i} \right. \\
& \left. - \frac{1}{a_i} g_i^T [S_W^R(\Omega_i(a))]g_i \right) \mathrm{d}t \\
& \left. + (2\xi g_i^T [S_W^R(\Omega_i(a))]g_i)^{\frac{1}{2}} \mathrm{d}w_i \right)
\end{split}
\tag{14}
\end{equation}
$$

$$
d\varphi_i = -\xi \frac{\eta_i}{a_i} dt + \frac{1}{a_i} (2\xi g_i^T [S_W^R(\Omega_i(a))] g_i)^{\frac{1}{2}} d\bar{w}_i \quad (15)
$$

où $[S_W^R(\cdot)]$ désigne la partie réelle de $[S_W(\cdot)]$ (sup-
posée définie positive), $\bar{D}_i = 1 + (2\Omega_i)^{-1}a_i \frac{\partial \Omega_i}{\partial a_i}$,
$g_i^T = (\Omega_i \bar{D}_i)^{-1} \Psi_{eq,i}^T [S]$,
$\eta_i = \frac{2}{\Omega_i \bar{D}_i} \frac{\partial}{\partial a_i} (\Omega_i \bar{D}_i \times$
$\int_0^\infty g_i^T [S_W^R(\nu)] g_i \sin(\Omega_i \nu) d\nu)$,
$\{w_i(t)\}$ et $\{\bar{w}_j(t)\}$, pour $i = 1, \dots, p$ et $j = 1, \dots, p$, représentent 2p processus de
Wiener unité, scalaires et indépendants. Par comparai-
son avec l'équation (7), les $g_i$ sont maintenant fonction
de a par l'intermédiaire des fréquences $\Omega_i$ et les $c_i$ restent
inchangés.

L'équation (14) permet d'exprimer la densité de pro-
babilité de a et, suivant [12], de construire un SLPA qui
 conduit à une approximation de la matrice de DSP. Cette
 approche fait appel au calcul d'intégrales de dimension p.
 Afin de simplifier la mise en œuvre, le vecteur amplitude a
 est remplacé par une variable scalaire qui tient compte de
 la répartition énergétique de chaque composante modale.

3.3 Décomposition modale avec fréquences fonction
d'une variable scalaire

Soit le processus « énergie modale » défini par [8,9]

$$
A^2(t) = \sum_{j=1}^{p} a_j^2(t) \qquad (16)
$$

L'évolution de A, déduite des équations (14) et (15), est
régie par l'EDSI

$$
dA = -\xi \sum_{i=1}^{p} \left( \frac{a_i}{A} f_i + \left( \frac{g_i^T [S_W^R(\Omega_i(a))] g_i}{A} \right) \left( 1 - \frac{a_i^2}{A^2} \right) \right) dt \\
+ \sum_{i=1}^{p} \frac{a_i}{A} \left( 2\xi g_i^T [S_W^R(\Omega_i(a))] g_i \right)^{\frac{1}{2}} dw_i \quad (17)
$$

avec la notation $f_i = \frac{c_i}{D_i} a_i - g_i^T \frac{\partial [S_W^R(\Omega_i(a))]g_i}{\partial a_i}$ —
$\frac{1}{a_i} g_i^T [S_W^R(\Omega_i(a))]g_i$. Cette EDSI dépend des amplitudes
modales et ne suffit pas à caractériser le processus $\{A(t)\}$.

De façon à découpler les équations (14) et (17),
l'influence des composantes $a_i(t)$ est prise en compte
dans (17) par une projection de la forme $a_i(t) = \delta_i A(t)$.
Chaque coefficient est déterminé en minimisant, par rap-
port à $\delta_i^2$, l'expression

$$
\mathbb{E} [a_i^2(t) - \delta_i^2 A^2(t)]^2 \quad (18)
$$

On vérifie simplement que $\sum_{i=1}^p \delta_i^2 = 1$. Le facteur de participation $\delta_i^2$ rend compte de la part de l'énergie de la composante $i$ à l'énergie totale.

Le calcul des $\delta_i$ nécessite de connaître la loi de probabilité invariante de (14). On se restreint ici à l'approximation donnée par la LSG. Les $a_i$ étant alors stochastiquement indépendantes et distribuées suivant une loi de Rayleigh, les $\delta_i$ s'expriment sous la forme

$$
\delta_i^2 = \frac{\mathbb{E} a_i^4 + \mathbb{E} a_i^2 \sum_{j \neq i} \mathbb{E} a_j^2}{\sum_i \mathbb{E} a_i^4 + \sum_{k \neq j} \mathbb{E} a_k^2 \mathbb{E} a_j^2} \quad (19)
$$

et se calculent simplement à l'aide des moments donnés
par le SLE (4).

De l'EDSI (17), découplée en substituant $\delta_i A$ à $a_i$, on
déduit la densité de probabilité invariante de $\{A(t)\}$ (sous
l'hypothèse de la convergence des intégrales)

$$
p(A) = C_N \left[ \frac{A}{\sqrt{\alpha(A)}} \exp \left( - \int_{0}^{A} \frac{\beta(a)}{\alpha(a)} a \, da \right) \right. \\
\left. \times \exp \left( \int_{0}^{A} \frac{\gamma(a)}{a \alpha(a)} da \right) \right] \quad (20)
$$

avec $\alpha(A) = \sum_i \delta_i^2 g_i^T(A) [S_W^R(\Omega_i(A))] g_i(A)$, $\beta(A) = \sum_i \delta_i^2 \frac{c_i}{D_i(A)}$, $\gamma(A) = \sum_i 2(1-\delta_i^2) g_i^T(A) [S_W^R(\Omega_i(A))] g_i(A)$, $C_N$ constante de normalisation. Dans ces expressions, la dépendance de $g_i$ vis-à-vis des $p$ coefficients $\delta_j$ est omise.
En substituant $\delta_i A$ à $a_i$, la décomposition modale (8) devient

$$X(t) = \sum_{i=1}^{p} \Psi_{eq,i} \delta_i A(t) \cos \phi_i(t) \quad (21)$$

et les fréquences déduites de (13), avec $a_i = \delta_i A$, s'écrivent

$$\Omega_i^2(A) = \frac{2}{(2\pi)^p} \int_{[0,2\pi]^p} \Psi_{eq,i}^T \left[ \frac{\partial}{\partial X} F(X(t)) \right] \times \Psi_{eq,i} \sin^2 \phi_i \, d\phi_1 \cdots d\phi_p \quad (22)$$

On dispose maintenant d'une famille de modes $(\Omega_i(A), \Psi_{eq,i})$ dont les fréquences dépendent d'une variable scalaire. Le caractère non intrinsèque des modes est accentué par la dépendance des $\Omega_i$ vis-à-vis des $\delta_j$.

### 3.4 Approximation de la DSP

Par analogie avec l'analyse modale linéaire et suivant la méthodologie proposée dans [12], à chaque mode non-linéaire est associée une composante modale dont la dynamique est décrite par le SLPA

$$\ddot{q}_i(t) + 2\xi c_i \dot{q}_i(t) + \Omega_i^2(A) q_i(t) = \sqrt{4\xi} B_i(A) \Psi_{eq,i}^T [S] W(t) \quad (23)$$

où $\Omega_i$, $\Psi_{eq,i}$ et $c_i$ sont définis en (22), (6) et (7). $\mathcal{A}$ est une variable aléatoire scalaire (indépendante de $\{W(t)\}$) de densité de probabilité $p(A)$.

La fonction $B_i$ est introduite pour assurer la conservation du moment d'ordre 2 entre la réponse du système linéaire (23) pour $\mathcal{A} = A$ et la i-ème composante de (21). Ce qui donne

$$B_i(A) = A \sqrt{\frac{\delta_i^2 c_i \Omega_i^2(A)}{2s_i^2(\Omega_i(A))}} \quad (24)$$

avec $s_i^2(\omega) = \Psi_{eq,i}^T[S][S_W(\omega)][S^T]\Psi_{eq,i}$.

Finalement une approximation de la matrice de DSP de la réponse non-linéaire est donnée par

$$[S_X(\omega)] \approx \int_{\mathbb{R}_+} [\Psi_{eq}]^T [S_q(\omega; A)][\Psi_{eq}] p(A) dA \quad (25)$$

où $[S_q(\omega; A)]$ est la matrice de DSP de la réponse de (23) pour $\mathcal{A} = A$. La matrice de DSP de l'excitation $[S_W(\omega)]$ est nécessaire dans le calcul des $\Psi_{eq,i}$, $\delta_i$, de $p(A)$ et $[S_q(\omega)]$.

### 3.5 Procédure numérique

L'approximation (25) repose sur des développements analytiques techniques mais sa mise en œuvre numérique,

faisant appel à des procédures classiques, a la structure simple suivante :

1. calcul de la matrice $[K_{eq}]$ du système LSG (4);

2. résolution du problème aux valeurs propres (6);

3. calcul des facteurs de participation $\delta_i^2$ (19);

4. calcul de la matrice de DSP (25) par formule de quadrature nécessitant l'évaluation

(a) des fréquences $\Omega_i(A)$ (22);

(b) de la probabilité $p(A)$ (20), en chaque point d'intégration.

Si $\mathcal{F}(U)$ est une fonction polynomiale, l'expression (22) s'obtient sous forme explicite en fonction des coefficients du polynôme simplifiant ainsi l'étape 4(a).

## 4 Problème inverse : identification des non-linéarités

Les termes linéaires de l'équation (1) étant supposés connus, un modèle paramétrique, $\mathcal{F}(U, \theta)$, est choisi pour décrire la non-linéarité, $\theta$ fait référence au vecteur des paramètres.

L'identification de $\theta$ est réalisée par calage de la DSP approchée de la réponse $U$ sur la mesure correspondante. La fonction de coût retenue s'écrit

$$\begin{align*} J(\theta) &= \sum_{j=1}^{n} J_j(\theta) \\ &= \sum_{j=1}^{n} \sum_{i=N_1}^{N_2} (\log_{10}(S_{U_j}^{\theta}(\omega_i)) - \log_{10}(S_{U_j}^{m}(\omega_i)))^2 \quad (26) \end{align*}$$

où $S_{U_j}^{\theta}(\omega_i)$ représente l'approximation théorique de la i-ème composante de $U(t)$, donnée par (3) et (25), qui dépend de $\theta$ par les modes non-linéaires, alors que $S_{U_j}^m(\omega_i)$ désigne la DSP mesurée associée. La bande de fréquences sur laquelle est réalisée l'identification, suggérée par la réponse expérimentale, fixe les valeurs de $N_1$ et $N_2$. Le problème de miminum non-convexe est résolu par un algorithme classique.

## 5 Application : conditions aux limites non-linéaires

### 5.1 Dispositif expérimental et modèle

Le banc d'essais a été élaboré à l'École Centrale de Lyon dans le cadre du COST Structural Dynamics. Le montage expérimental est schématisé figure 1. La structure vibrante comprend une poutre en acier, encastrée en $z = 0$, et une non-linéarité localisée à l'extrémité $z = L = 0,633$ m, réalisée à l'aide d'une lame d'acier excitée en grands déplacements. La poutre est de section carrée $A = 0,196 \times 10^{-3}$ m$^2$.

Le signal d'entrée est un processus bruit blanc à bande limitée (incluant les deux premières fréquences
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Schéma du dispositif experimental. Détail des modèles des supports.</figcaption></figure>

de résonance), de moyenne nulle et d'intensité ajustable. Cette force d'excitation est délivrée par le vibreur *e*, mesurée par le transducteur *d* au point $z_e = 0,063$ m, les réponses en accélération sont mesurées par les accéléromètres *a*, *b* et *c* aux points $z_1 = 0,220$ m, $z_2 = 0,420$ m et $z_3 = 0,605$ m. Des mesures sont disponibles pour trois niveaux du signal d'entrée, faible (1,8 N), moyen (14,2 N) et fort (17,4 N). Le terme de niveau fait référence à l'écart-type de la DSP de l'entrée.

La poutre est modélisée suivant la théorie d'Euler-Bernoulli en supposant connus les paramètres géométriques et physiques du matériau (module de Young : $E = 2 \times 10^{11}$ N·m⁻², masse volumique : $\rho = 7900$ kg·m⁻³).

À l'extrémité $z = 0$ les conditions aux limites retenues sont

$$ -EI \frac{\partial^2 u(0)}{\partial z^2} = k_r^0 \frac{\partial u(0)}{\partial z} \quad (27) $$

$$ EI \frac{\partial^3 u(0)}{\partial z^3} = 0 \quad (28) $$

$k_r^0$ désignant le coefficient de raideur de rotation (voir 1, Fig. 1).

L'élément de jonction entre la poutre et la lame est représenté par une masse ponctuelle *m*, la non-linéarité localisée (voir 2–3, Fig. 1) est modélisée par les conditions aux limites suivantes

$$ -EI \frac{\partial^2 u(L)}{\partial z^2} = k_r^L \frac{\partial u(L)}{\partial z} + \lambda_r^L \left( \frac{\partial u(L)}{\partial z} \right)^3 \quad (29) $$

$$ EI \frac{\partial^3 u(L)}{\partial z^3} = k_d^L u(L) + \lambda_d^L (u(L))^3 \quad (30) $$

où les raideurs de rotation $k_r^L$ et de translation $k_d^L$ caramérisent la partie linéaire, $\lambda_r^L$ et $\lambda_d^L$ la partie non-linéaire du modèle. Cette forme cubique simplifie le calcul des $\Omega_i(A, \delta)$.

Une méthode d'éléments finis fournit les termes non dissipatifs de (1), l'équation (2) est alors obtenue par projection sur les deux premiers modes ($p = 2$), seuls disponibles expérimentalement. De façon classique une matrice diagonale d'amortissement $[C]$ est introduite.

Les coefficients de raideur linéaire $k_r^0$, $k_d^L$, $k_r^L$ et les taux d'amortissement $\xi_1$, $\xi_2$ sont d'abord identifiés, ces valeurs étant connues on recherche les paramètres de raideur non-linéaire $\lambda_r^L$ et $\lambda_d^L$.

<table><caption>**Tableau 1.** Paramètres linéaires identifiés.</caption><thead><tr><th>$k_r^0$ (N·m)</th><th>$k_d^L$ (N·m⁻¹)</th><th>$k_r^L$ (N·m)</th><th>$\xi_1$</th><th>$\xi_2$</th></tr></thead><tbody><tr><td>7192</td><td>3444</td><td>234</td><td>0,0020</td><td>0,0019</td></tr></tbody></table>

<table><caption>**Tableau 2.** Fréquences naturelles identifiées.</caption><thead><tr><th>$f_1^\theta$ (Hz)</th><th>$f_1^m$ (Hz)</th><th>$f_2^\theta$ (Hz)</th><th>$f_2^m$ (Hz)</th></tr></thead><tbody><tr><td>24,907</td><td>25,000</td><td>134,063</td><td>134,100</td></tr></tbody></table>

## 5.2 Identification des paramètres linéaires

L'identification de $k_r^0$, $k_d^L$, $k_r^L$, $\xi_1$ et $\xi_2$ se fait avec les données relatives à l'entrée faible niveau (1,8 N). La fonction coût tient compte des accélérations en $z_1$, $z_2$ et $z_3$ soit

$$ J(\theta) = J_1(\theta) + J_2(\theta) + J_3(\theta) \quad (31) $$

$J_j(\theta)$ est donné par (26), $N_1$ et $N_2$ correspondent à une bande de 20 Hz autour des fréquences de résonance. Les résultats sont donnés tableau 1. Les valeurs des deux premières fréquences de résonance, $f_i^\theta = \frac{\omega_i^\theta}{2\pi}$, calculées pour ces paramètres, sont données tableau 2 où $f_i^m$ désigne les valeurs mesurées correspondantes.

## 5.3 Identification des paramètres non-linéaires

L'identification des paramètres non-linéaires $\lambda_d^L$ et $\lambda_r^L$ se fait à fort niveau d'entrée, avec une fonction coût et une plage de fréquences identiques au cas linéaire. Les valeurs obtenues sont données tableau 3 (méthode 1).

À titre de comparaison, les résultats obtenus en [13], à partir de la méthode globale (méthode 2), sont également indiqués.

Pour illustrer l'insuffisance de l'approximation LSG à décrire la DSP de la réponse à une excitation de fort niveau (17,4 N), la DSP obtenue par LSG est comparée à la DSP mesurée (Fig. 2a). Le modèle LSG reproduit assez bien la première fréquence de résonance, suretime son amplitude et ignore l'étalement spectral. Par contre l'identification à l'aide du LSPA conduit à une DSP en bon accord avec la mesure (Fig. 2b). Pour le mode 1, la translation de la fréquence de résonance et l'étalement spectral autour de celle-ci sont bien reproduits. La non-linéarité a peu d'effet sur la seconde résonance. La résonance secondaire (sous-harmonique et/ou effet de couplage) qui apparaît sur la réponse expérimentale n'est pas prise en compte par le modèle et ne peut donc pas être retrouvée.

Pour valider les résultats, l'approximation de la réponse à un niveau moyen d'entrée (14,2 N) calculée avec les paramètres identifiés à fort niveau (17,4 N), se trouve en bon accord avec la réponse expérimentale (Fig. 3a). Inversement, la réponse théorique à fort niveau d'entrée calculée avec les paramètres identifiés à niveau moyen, ne reproduit pas entièrement l'étalement spectral (Fig. 3b).
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Approximations de la DSP de la réponse en $z_3$ (accéléromètre c) à fort niveau de sollicitation, fournies (a) par le modèle LSG et (b) par le modèle SLPA identifiés à fort niveau de sollicitation. Ligne continue : modèle ; ligne pointillée : mesure.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Approximations de la DSP de la réponse en $z_3$ (accéléromètre c) (a) à moyen niveau de sollicitation, le modèle SLPA étant identifié à fort niveau et (b) à fort niveau de sollicitation le modèle SLPA étant identifié à moyen niveau. Ligne continue : modèle ; ligne pointillée : mesure.</figcaption></figure>

<table><caption>**Tableau 3.** Paramètres non-linéaires identifiés.</caption><thead><tr><th>méthode</th><th>niveau (N)</th><th>δ<sub>1</sub></th><th>δ<sub>2</sub></th><th>λ<sup>L</sup><sub>d</sub> (N·m<sup>-3</sup>)</th><th>λ<sup>L</sup><sub>r</sub> (N·m)</th></tr></thead><tbody><tr><td>1</td><td>17,4</td><td>0,9523</td><td>0,3053</td><td>2,73 × 10<sup>9</sup></td><td>-7,97 × 10<sup>6</sup></td></tr><tr><td>2</td><td>17,4</td><td></td><td></td><td>2,49 × 10<sup>9</sup></td><td>-5,11 × 10<sup>6</sup></td></tr><tr><td>1</td><td>14,2</td><td>0,9580</td><td>0,2869</td><td>2,18 × 10<sup>9</sup></td><td>-12,55 × 10<sup>6</sup></td></tr><tr><td>2</td><td>14,2</td><td></td><td></td><td>1,78 × 10<sup>9</sup></td><td>-7,61 × 10<sup>6</sup></td></tr></tbody></table>

L'information modale donnée par le SLPA peut être utilisée pour compléter les résultats. Le comportement, en fonction de A, des fréquences $Ω_1$ et $Ω_2$, calculées à partir de l'équation (22), est visualisé sur les figures 4a et 4b respectivement, pour les niveaux moyen et fort. On constate que le premier mode est durcissant, le second molissant. La faible différence entre les niveaux moyen et fort se traduit par un faible écart entre les fréquences correspondantes. Les formes modales, non représentées, sont très peu affectées par la non-linéarité dans le modèle LSG.

La figure 5a représente, pour une excitation fort niveau, le mode 1 calculé à l'aide du modèle identifié à fort

niveau (ligne continue) et à moyen niveau (ligne pointillée). La figure 5b représente, dans les mêmes conditions, la densité de probabilité. L'étalement en fréquences insuffisant, constaté sur la DSP (Fig. 3b), est dû à l'écart entre les modes, les densités de probabilité ayant pratiquement le même support.

Pour une excitation moyen niveau, le calcul du mode 1 à partir du modèle identifié à fort niveau ou à moyen niveau, donne des fréquences très proches, indicernables sur la figure 6a. Il en est de même pour la densité de probabilité, représentée figure 6b. Cette situation est en accord avec l'étalement en fréquences de la DSP bien reproduit (Fig. 3a).
<figure><img src="image_5.png" /><figcaption>Fig. 4. Évolutions (a) de la pulsation modale 1 et (b) de la pulsation modale 2 fournies par les modèles SLPA identifiés à fort niveau (ligne continue) et à moyen niveau (ligne pointillée).</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 5. Évolutions à fort niveau de sollicitation (a) de la pulsation modale 1 et (b) de la densité de probabilité fournies par les modèles SLPA identifiés à fort niveau (ligne continue) et à moyen niveau (ligne pointillée).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 6. Évolutions à moyen niveau de sollicitation (a) de la pulsation modale 1 et (b) de la densité de probabilité fournies par les modèles SLPA identifiés à fort niveau (ligne continue) et à moyen niveau (ligne pointillée).</figcaption></figure>
Références

[1] C. Soize, O. Le Fur, Modal identification of weakly non-linear multidimensional dynamical systems using a stochastic linearization method with random coefficients, Mechanical Systems and Signal Processing, 11 (1997) 37–49

[2] J.S. Bendat, P.A. Palo, R.N. Coppolino, A general identification technique for nonlinear differential equation of motion, Probabilistic Engineering Mechanics, 7 (1992) 43–61

[3] A.F. Vakakis, L.I. Manevitch, Y.V. Mikhlin, V.N. Pilipchuk, Zevin A.A., Normal Modes and Localization in Nonlinear Systems, John Wiley & Sons, New York, 1996

[4] R. Bouc, The power spectral density of response for a strongly nonlinear random oscillator, Journal of Sound and Vibration, 175 (1994) 317–331

[5] C. Soize, Stochastic linearization method with random parameters for SDOF nonlinear dynamical systems : prediction and identification procedures, Probabilistic Engineering Mechanics, 10 (1995) 143–152

[6] R. Bouc, M. Defilippi, Multimodal nonlinear spectral response of a beam with impact under random input, Probabilistic Engineering Mechanics, 12 (1997) 163–170

[7] S. Bellizzi, R. Bouc, Analysis of multi-degree of freedom strongly nonlinear systems with random input. Part I : Nonlinear modes and stochastic averaging, Probabilistic Engineering Mechanics, 14 (1999) 229–244

[8] R. Bouc, M. Defilippi, A simplified approach of solving multidimensional nonlinear random vibration with application to nonlinear plates. Dans P.D. Spanos éditeur, Computational Stochastic Mechanics, Third International Conference on Computational Stochastic Mechanics, Greece june 1998, A.A. Balkama, Rotterdam, Brookfield, 1999

[9] R. Bouc, Multidimensional Nonlinear Random Vibration : A Linearization Method using Nonlinear Coupled Modes, in Bouc R., Soize C. (Eds.), Progress in stochastic structural dynamics, Journée Nationale Dynamique Stochastique des Structures, Châtillon, France, juin 1999, Publication du LMA-CNRS Marseille, 1999

[10] J.B Roberts, P.D. Spanos, Random Vibration and Statistical Linearization, John Wiley & Sons, Chichester, 1990

[11] J.B. Roberts, P.D. Spanos, Stochastic averaging : an approximate method of solving random vibration problems, International Journal of Non-Linear Mechanics, 21 (1986) 111–134

[12] S. Bellizzi, R. Bouc, Analysis of multi-degree of freedom strongly nonlinear systems with random input. Part II : equivalent linear system with random matrices and power spectral density matrix, Probabilistic Engineering Mechanics 14, (1999) 245–256

[13] S. Bellizzi, M. Defilippi, Nonlinear mechanical systems identification using linear systems with random parameters, International Conference on Structural Systems Identification, Université de Kassel, Allemagne, septembre 2001