r5.07.01 Calcul de modes non-linéaires avec l’opérateur MODE_NON_LINE#

Résumé:

L’opérateur MODE_NON_LINE permet de calculer les modes non-linéaires d’un système linéaire conservatif autonome doté de non-linéarités de choc localisées. On décrit ici les méthodes numériques utilisées pour le calcul des modes non-linéaires.

Définitions des modes non-linéaires#

La première définition a été donnée par Rosenberg dans les années 1960, mais c’est au début des années 1990 que la définition la plus solide théoriquement des modes non-linéaires a été donnée en utilisant la théorie des systèmes dynamiques. Ainsi un mode non-linéaire est défini comme «une variété invariante de dimension 2 de l’espace des phases, tangente à un point d’équilibre stable du système linéarisé» (cf. [Bib2]).

Pour les systèmes conservatifs, on peut définir un mode non-linéaire comme une «famille de solutions périodiques». Cette définition est beaucoup plus attractive car donnant accès à des outils numériques performants comme les méthodes de continuation. C’est pourquoi c’est cette définition qui a été choisi ici, par conséquent MODE_NON_LINE permet le calcul de famille de solutions périodiques.

Modélisation#

On suppose qu’en l’absence de contacts, la structure étudiée a un comportement linéaire et que les équations régissant son équilibre dynamique ont été discrétisées par différences finies ou par éléments finis. On obtient un système discret d’équations différentielles du second ordre à \(n\) degrés de libertés.

De façon générale, ces équations prennent la forme suivante:

\(\mathrm{M}.\ddot{\mathrm{U}}(t)+\mathrm{K}.\mathrm{U}(t)+{\mathrm{F}}^{\mathrm{nl}}(\mathrm{U}(t))=0,\)

  • \(\mathrm{M}\) est la matrice de masse du système,

  • \(K\) est la matrice de rigidité élastique du système,

  • \({\mathrm{F}}^{\mathrm{nl}}\) est le vecteur des forces non-linéaires de chocs,

On se place dans le cadre de non-linéarités localisées, c’est-à-dire que les forces s’appliquent sur un nombre relativement faible de degrés de liberté. On notera \({I}_{\mathit{nl}}\) l’ensemble des indices de ces degrés de liberté, et \({I}_{l}\) son complémentaire, c’est-à-dire \({I}_{n}=\lbrace 1,\dots ,n\rbrace ={I}_{\mathit{nl}}\cup {I}_{l}\) . La taille du vecteur \({I}_{\mathit{nl}}\) est égale à \({n}_{\mathit{nl}}\) .

Pour pouvoir utiliser l’algorithme EHMAN qui combine la méthode d’équilibrage harmonique (EH) et la méthode asymptotique numérique (MAN), il est nécessaire d’adopter un certain formalisme. Celui-ci consiste à réécrire le problème sous la forme d’un système d’ordre un, possédant des non-linéarités au plus quadratiques, sois un système d’inconnue \(S\) de la forme

(3070)#\[A(\dot{S})=C+L(S)+Q(S,S),\]

où \(C\) est une constante, \(L\) est une forme linéaire et \(Q\) une forme quadratique. On regroupe en général les termes du membre de droite sous l’appellation \(R(S)\) , le problème à résoudre se mets alors sous forme condensée

(3071)#\[A(\dot{S})=R(S)\]

Pour cela, on introduit deux vecteurs d’inconnues supplémentaires \({\mathrm{F}}^{\mathrm{nl}}(t)\) et \(\mathrm{Z}(t)\) , ainsi que des équations supplémentaires pour définir ces vecteurs. Le vecteur \({\mathrm{F}}^{\mathrm{nl}}(t)\) décrit la force de contact et le vecteur \(\mathrm{Z}(t)\) est appelé vecteur des variables auxiliaires. Le vecteur \(\mathrm{Z}(t)\) est de taille \({n}_{z}\) . Cette réécriture relativement abstraite, détaillée par le biais d’exemples, conduit à résoudre un système de la forme

(3072)#\[\begin{split}\lbrace \begin{array}{c}\dot{\mathrm{U}}(t)=\mathrm{V}(t),\\ \mathrm{M}.\dot{\mathrm{V}}(t)=-\mathrm{K}\mathrm{U}(t)-\left(\begin{array}{c}0\\ {\mathrm{F}}^{\mathrm{nl}}(\mathrm{U}(t))\end{array}\right),\\ 0=\mathrm{G}({\mathrm{U}}^{\mathrm{nl}}(t),{\mathrm{F}}^{\mathrm{nl}}(t),\mathrm{Z}(t)),\end{array}\end{split}\]

où \(\mathrm{G}\) est un opérateur quadratique, qui définit le contact de façon implicite. \(\mathrm{G}\) peut-être décomposé de la manière suivante

(3073)#\[\mathrm{G}(\mathrm{S}(t))={\mathrm{G}}_{\mathrm{c}}+{\mathrm{G}}_{l}(\mathrm{S}(t))+{\mathrm{G}}_{\mathrm{q}}(\mathrm{S}(t),\mathrm{S}(t)).\]

\({\mathrm{G}}_{\mathrm{c}}\) est une application constante, \({\mathrm{G}}_{l}\) linéaire et \({\mathrm{G}}_{\mathrm{q}}\) bilinéaire.

Contact unilatéral#

../../../../_images/1000000000000251000001C246248D4DB86753D8.png ../../../../_images/100000000000035A0000029DD35F6E0747EBC59F.png

La non-linéarité de contact unilatéral, sous sa forme non-régulière, s’écrit :

(3074)#\[\begin{split}{\mathrm{f}}^{\mathrm{nl}}(t)=\lbrace \begin{array}{c}\alpha ({u}^{\mathit{nl}}(t)-g)\mathit{si}{u}^{\mathit{nl}}(t)⩾g,\\ 0\mathit{si}{u}^{\mathit{nl}}(t)⩽g,\end{array}\end{split}\]

où \(\alpha\) est la raideur de contact, \(g\) le jeu entre le nœud de la structure et le nœud de contact. La régularisation proposée s’écrit sous la forme implicite

(3075)#\[{\mathrm{f}}^{\mathrm{nl}}(t)({\mathrm{f}}^{\mathrm{nl}}(t)-\alpha ({u}^{\mathit{nl}}(t)-g))-\alpha \eta =0,\]

où \(\eta\) représente le paramètre de régularisation.

Contact bilatéral#

../../../../_images/1000000000000388000001C63E14CC67E35C27B0.png ../../../../_images/1000000000000351000002A464ADF2963A80BE90.png

La non-linéarité de contact bilatéral, sous sa forme non-régulière, s’écrit

(3076)#\[\begin{split}{\mathrm{f}}^{\mathrm{nl}}(t)=\lbrace \begin{array}{c}\alpha ({u}^{\mathit{nl}}(t)-g)\mathit{si}{u}^{\mathit{nl}}(t)⩾g,\\ 0\mathit{si}\mid {u}^{\mathit{nl}}(t)\mid ⩽g,\\ \alpha ({u}^{\mathit{nl}}(t)+g)\mathit{si}{u}^{\mathit{nl}}(t)⩽-g,\end{array}\end{split}\]

où \(\alpha\) est la raideur de contact, \(g\) le jeu entre le nœud de la structure et le nœud de contact. La régularisation proposée s’écrit sous la forme implicite suivante

(3077)#\[{f}^{\mathit{nl}}(t)({f}^{\mathit{nl}}(t)-\alpha ({u}^{\mathit{nl}}(t)-g))({f}^{\mathit{nl}}(t)-\alpha ({u}^{\mathit{nl}}(t)+g))-{\alpha}^{2}\eta {u}^{\mathit{nl}}(t)=0,\]

où \(\eta\) représente le paramètre de régularisation. On peut le réécrire sous une forme quadratique,

(3078)#\[\begin{split}\lbrace \begin{array}{c}{\alpha}^{2}{g}^{2}{f}^{\mathit{nl}}(t)-\eta {\alpha}^{2}{u}^{\mathit{nl}}(t)-{f}^{\mathit{nl}}(t)z(t)=0,\\ z(t)-{({f}^{\mathit{nl}}(t)-\alpha {u}^{\mathit{nl}}(t))}^{2}=0.\end{array}\end{split}\]

Contact annulaire#

../../../../_images/10000000000002430000024A75452B565D1CA9A2.png ../../../../_images/10000000000003700000030292CC0E0A415B890B.png

On introduit \(z(t)\) , la distance radiale par rapport à la position d’équilibre \(({e}_{x},{e}_{y})\) du nœud de choc, dont la position courante est \(({u}_{x}^{\mathit{nl}},{u}_{y}^{\mathit{nl}})\) . \(z(t)\) s’écrit

(3079)#\[z(t)=\sqrt{{\left({u}_{x}^{\mathit{nl}}(t)-{e}_{x}\right)}^{2}+{\left({u}_{y}^{\mathit{nl}}(t)-{e}_{y}\right)}^{2}}.\]

La non-linéarité de contact annulaire, sous sa forme non-régulière, s’écrit

(3080)#\[\begin{split}{\mathrm{f}}^{\mathrm{nl}}(t)=\lbrace \begin{array}{c}\alpha (z(t)-g)\mathit{si}z(t)>g,\\ 0\mathit{si}z(t)⩽g,\end{array}\end{split}\]

où \(g\) désigne le rayon de l’anneau (autrement dit le jeu) et \(\alpha\) (\(\alpha >0\) ) désigne le coefficient de raideur associé. En introduisant les composantes \({f}_{x}^{\mathit{nl}}\) et \({f}_{y}^{\mathit{nl}}\) de la force non-linéaire dans les directions \(x\) et \(y\) , la régularisation proposée s’écrit sous la forme implicite suivante

(3081)#\[\begin{split}\lbrace \begin{array}{c}{\mathrm{f}}^{\mathrm{nl}}(t)\left({\mathrm{f}}^{\mathrm{nl}}(t)-\alpha (z(t)-g)\right)-\alpha \eta =0,\\ {z}^{2}(t)-{({u}_{x}^{\mathit{nl}}(t)-{e}_{x})}^{2}-{({u}_{y}^{\mathit{nl}}(t)-{e}_{y})}^{2}=0,\\ {f}_{x}^{\mathit{nl}}(t)z(t)–{\mathrm{f}}^{\mathrm{nl}}(t){u}_{x}^{\mathit{nl}}(t)=0,\\ {f}_{y}^{\mathit{nl}}(t)z(t)–{\mathrm{f}}^{\mathrm{nl}}(t){u}_{y}^{\mathit{nl}}(t)=0,\end{array}\end{split}\]

où \(\eta\) représente le paramètre de régularisation.

Procédure de calcul des modes non-linéaires#

On cherche les solutions périodiques du système

(3082)#\[\begin{split}\lbrace \begin{array}{c}\mathrm{M}.\ddot{\mathrm{U}}(t)+\mathrm{K}\mathrm{U}(t)+\left(\begin{array}{c}0\\ {\mathrm{F}}^{\mathrm{nl}}(\mathrm{U}(t))\end{array}\right)=0,\\ \mathrm{G}({\mathrm{U}}^{\mathrm{nl}}(t),{\mathrm{F}}^{\mathrm{nl}}(t),\mathrm{Z}(t))=0,\end{array}\end{split}\]

en écrivant le vecteur \(\mathrm{U}(t)\) à l’aide d’une série de Fourier que l’on tronque à l’ordre \({H}_{l}\) (en supposant que l’influence des harmoniques supérieures est négligeable), soit

(3083)#\[\mathrm{U}(t)={\mathrm{U}}_{0}+\sum_{k=1}^{{H}_{l}}{\mathrm{U}}_{\mathrm{Ck}}\cos(k\omega t)+{\mathrm{U}}_{\mathrm{Sk}}\sin(k\omega t).\]

On effectue le même développement pour les vecteurs d’inconnues \({\mathrm{F}}^{\mathrm{nl}}(t)\) et \(\mathrm{Z}(t)\) , mais avec un ordre de troncature différent, qu’on appelle \({H}_{\mathit{nl}}\) . Il vient alors

(3084)#\[\begin{split}\lbrace \begin{array}{c}{\mathrm{F}}^{\mathrm{nl}}(t)={\mathrm{F}}_{0}^{\mathrm{nl}}+\sum_{k=1}^{{H}_{\mathit{nl}}}{\mathrm{F}}_{\mathrm{Ck}}^{\mathrm{nl}}\cos(k\omega t)+{\mathrm{F}}_{\mathrm{Sk}}^{\mathrm{nl}}\sin(k\omega t),\\ \mathrm{Z}(t)={\mathrm{Z}}_{0}+\sum_{k=1}^{{H}_{\mathit{nl}}}{\mathrm{Z}}_{\mathrm{Ck}}\cos(k\omega t)+{\mathrm{Z}}_{\mathrm{Sk}}\sin(k\omega t).\end{array}\end{split}\]

Nota Bene: Le choix d’ordre différent découle directement des grandeurs représentées. Pour un système choquant peu, au voisinage du mode linéaire, la réponse \(\mathrm{U}(t)\) sera proche d’un sinus pur. On choisira donc un nombre d’harmonique faible. On augmentera \({H}_{l}\) lorsqu’on cherchera à obtenir des solutions associées à des niveaux d’énergie élevés. En revanche, même pour des forces de choc faible, le contenu spectral peut être important, en particulier si le choc est raide. Il faudra donc retenir un nombre important d’harmonique, et on aura donc en général \({H}_{\mathit{nl}}\gg {H}_{l}\) *.*

En pratique, on pourra commencer à réaliser un calcul en choisissant des ordres relativement faible, pour avoir une idée du comportement principal. En revanche, avec un nombre d’harmonique faible, on ne sera pas en mesure de capter correctement les bifurcations. On augmentera donc graduellement les valeurs de \({H}_{l}\) et \({H}_{\mathit{nl}}\) pour faire apparaître des comportements plus complexes.

Après avoir effectué un équilibrage harmonique, on obtient un système algébrique sous-déterminé dont les inconnues sont les coefficients de Fourier de \(\mathrm{U}(t)\) , \({\mathrm{F}}^{\mathrm{nl}}(t)\) et \(\mathrm{Z}(t)\) et l’inconnue supplémentaire, la pulsation propre \(\omega\) . On regroupe ces inconnues dans un vecteur unique afin de mettre le système sous la forme (). On peut alors utiliser la MAN tel qu’indiqué dans la référence [Bib3]. Pour cela, on développe en série entière le vecteur d’inconnu \(S\) en fonction d’un paramètre de chemin \(a\)

(3085)#\[S(a)={S}_{0}+\sum_{k=1}^{{N}_{\mathit{MAN}}}{a}^{k}{S}_{k},\]

où \({S}_{0}\) correspond au vecteur d’initialisation de l’algorithme, et \({S}_{k}\) pour \(k=1,\dots ,{N}_{\mathit{MAN}}\) sont les coefficients de la série entière qui restent à déterminer. Ensuite, on développe en série de Taylor la fonction \(R\) , au voisinage du vecteur d’initialisation \({S}_{0}\)

(3086)#\[0=R(S(a))=R({S}_{0})+\frac{\mathit{dR}}{\mathit{dS}}({S}_{0})\left(\sum_{k=1}^{{N}_{\mathit{MAN}}}{a}^{k}{S}_{k}\right)+\frac{{d}^{2}R}{{\mathit{dS}}^{2}}({S}_{0})\left(\sum_{k=1}^{{N}_{\mathit{MAN}}}{a}^{k}{S}_{k}\right)\left(\sum_{l=1}^{{N}_{\mathit{MAN}}}{a}^{l}{S}_{l}\right)\]

Cette relation est vraie quelque soit \(a\) , on passe alors à la résolution d’une suite de \({N}_{\mathit{MAN}}\) systèmes linéaires possédant la même matrice tangente, dépendant les uns des autres de manière récursive. La résolution de ces systèmes permet d’obtenir les vecteurs \({S}_{k}\) . A noter qu’on peut calculer analytique-ment cette matrice car la fonction \(R\) est quadratique, soit

(3087)#\[\frac{\mathit{dR}}{\mathit{dS}}({S}_{0}){e}_{i}=L({e}_{i})+Q({e}_{i},{S}_{0})+Q({S}_{0,}{e}_{i}),\]

où \({e}_{i}\) est le vecteur de la base canonique. On obtient ainsi une branche de solutions \(S(a)\) avec \(a\in [0,{a}_{max}]\) , où \({a}_{max}\) représente le domaine de validité de la série entière. Dès lors, on peut obtenir les branches de solutions périodiques qui forment les modes non-linéaires du modèle.

../../../../_images/10000000000004A80000030CFB22F3678E0FF2EA.png

La méthode permet également de localiser les bifurcations simples par la détection de série géométrique dans la représentation en série entière.

Bibliographie#

[Bib1] E.H. MOUSSI : “Analyses de structures vibrantes dotées de non-linéarités localisées à jeu à l’aide des modes non-linéaires” - Thèse de doctorat de l’université d’Aix-Marseille, 2013.

[Bib2] G. KERSCHEN, M. PEETERS, J.C. GOLINVAL et A.F. VAKAKIS : “Nonlinear normal modes, part I : A useful framework for the structural dynamicist” - Mechanical Systems and Signal Processing, 2009.

[Bib3] B. COCHELIN, C. VERGEZ : “A high order purely frequency-based harmonic balance formulation for continuation of periodic solutions” - Journal of Sound and Vibration, 2009.