r4.05.02 Approche stochastique pour l’analyse sismique

Résumé :

Ce document présente une méthode de calcul probabiliste pour déterminer la réponse d’une structure soumise à une excitation aléatoire de type sismique à partir des interspectres de l’excitation aux points d’appui de la structure. La réponse est elle-même exprimée sous forme d’interspectres.

Principe de la démarche

Position du problème considéré et principe général

On se place dans le cas d’une structure multi-appuyée , c’est-à-dire que la structure possède \(m\) degré de liberté-appuis, chacun étant soumis à sa propre excitation (non nécessairement égale partout). On suppose que la structure est représentée par un modèle éléments finis comportant \(n\) degrés de liberté. On cherche la réponse en un nombre fini (et faible) de \(l\) degrés de liberté.

On suppose que la grandeur excitation est de type mouvement imposé et se traduit par une famille d’accélérogrammes \({\mathrm{g}}_{j}(\mathrm{t})\) pour chacun des degrés de liberté-appuis \(j\) , \(j=1,m\) .

Le mouvement absolu de la structure est décomposé classiquement en mouvement d’entraînement et mouvement relatif.

Le calcul de la réponse en interspectres de puissance est réalisé par recombinaison modale .

A la suite de ce calcul modal, un calcul de réponse dynamique aléatoire se décompose en trois parties:

• définition de l’interspectre de puissance excitateur,

• calcul de l’interspectre de puissance réponse.

Ces deux premières parties font l’objet de la commande DYNA_ALEA_MODAL [U4.53.22].

La restitution de l’interspectre de puissance réponse sur base physique est réalisée avec la commande REST_SPEC_PHYS [U4.63.22].

• calcul de paramètres statistiques à partir de l’interspectre de puissance résultat.

Cette dernière étape est traitée par la commande POST_DYNA_ALEA [R7.10.01] [U4.84.04].

Décomposition du mouvement

Les décompositions et projections suivantes sont détaillées dans la documentation de référence relative à la résolution par calcul transitoire d’un calcul sismique [R4.05.01]. Nous n’en retenons ici que les grandes lignes.

Soit \({X}_{a}\) le vecteur déplacement absolu (de dimension \(n\) ) de tous les degrés de liberté de la structure.

La réponse totale dite absolue \({X}_{a}\) de la structure s’exprime comme la somme d’une contribution relative \({X}_{r}\) et de la contribution d’entrainement \({X}_{e}\) due aux déplacements d’ancrage (soumis à des accélérations représentées par un accélérogramme \({g}_{j}(t)\) en chacun des degrés de liberté-appuis \(j\) , \(j=1,m\) ).

Soient \(M,K\) et \(C\) les matrices de masse, de rigidité et d’amortissement du problème, limitées aux degrés de liberté non appuyés.

L’équation du mouvement s’écrit alors dans le repère lié au mouvement relatif:

\(\begin{array}{c}{{\begin{array}{}\text{M}\ddot{X}\end{array}}_{r}(t)+\text{C {}\dot{X}}_{r}(t)+{\text{K X}}_{r}(t)\text{= -}\text{M}\ddot{{X}_{e}}(t)+{F}_{\text{ext}}\end{array}\)

\({F}_{\text{ext}}\) : vecteur des forces extéŽrieures

En général les forces extérieures sont nulles lors d’un calcul de réponses sismiques.

Décomposition sur la base modale

Le calcul de réponse en interspectres de puissance est réalisé par recombinaison modale et fait appel, en mouvement imposé, à une base modale qui comprend à la fois des modes dynamiques et des modes statiques.

Soit \(\Phi =\left\lbrace {\phi }_{i,i=1,n}\right\rbrace\) la matrice \((n,n)\) des modes dynamiques calculés pour le système conservatif associé, en maintenant les m appuis bloqués.

Soit \(\Psi =\left\lbrace {\psi}_{j,j=1,m}\right\rbrace\) la matrice \((n,m)\) des modes statiques. Le mode \({\Psi}_{j}\) correspond à la déformée de la structure sous un déplacement unitaire imposé au degré de liberté-appui \(j\) , les autres degrés de liberté‑appuis étant bloqués.

Le déplacement imposé des ancrages \({X}_{s}(t)\) est relié à \({X}_{e}(t)\) par la relation : \({X}_{e}(t)=\Psi {X}_{s}(t)\) .

Les composantes de l’accélération des points d’ancrage \({\ddot{X}}_{s}(t)\) sont les accélérogrammes \({g}_{j}(t)\) , \(j=1,m\) .

On peut donc écrire \(\ddot{{X}_{e}}(t)={\Psi \ddot{X}}_{s}(t)=\sum_{j\text{=1}}^{m}{\psi}_{j}{g}_{j}(t)\text{}\) .

On effectue le changement de variable \({X}_{r}(t)=\Phi .q(t)\) , \(q(t)\) est le vecteur des coordonnées généralisées. En prémultipliant l’équation du mouvement par \({}^{T}\text{}\Phi\) , on obtient - en l’absence de forces extérieures autres que l’excitation sismique - l’équation projetée sur la base des modes dynamiques:

\({}^{T}\text{}\Phi M\Phi \ddot{q}(t)+{}^{T}\text{}\Phi C\Phi \dot{q}(t)+{}^{T}\text{}\Phi K\Phi q(t)=-{}^{T}\text{}\Phi M\Phi {\ddot{X}}_{s}(t)\)

On suppose que la matrice d’amortissement est une combinaison linéaire des matrices de masse et de rigidité (hypothèse d’amortissement de Rayleigh constant sur la structure ou hypothèse de Basile permettant un amortissement diagonal). La base \(\Phi\) qui orthogonalise les matrices \(M\) et \(K\) , orthogonalise donc aussi la matrice \(C\) .

Compte-tenu de cette hypothèse, l’équation précédente se décompose en \(n\) équations scalaires découplées sous la forme:

pour \(i=1,n\)

où l’on a noté :

\({\mu}_{i}={}^{T}\text{}{\Phi}_{i}M{\Phi}_{i}\) la masse modale

\({k}_{i}={}^{T}\text{}{\Phi}_{i}K{\Phi}_{i}\) la rigidité modale

\({\omega}_{i}=\sqrt{\frac{{k}_{i}}{{\mu}_{i}}}\) la pulsation modale

\({\xi}_{i}=\frac{{}^{T}\text{}{\Phi}_{i}.C.{\Phi}_{i}}{2{\mu}_{i}{\omega}_{i}}\) l’amortissement modal réduit »

\({p}_{ij}=\frac{{}^{T}\text{}{\Phi}_{i}.M.{\Psi}_{j}}{{\mu}_{i}}\) le facteur de participation modale de l’appui \(j\) sur le mode dynamique \(i\) .

La solution \({q}_{i}(t)\) de cette équation correspond à la réponse du mode dynamique \(i\) à l’ensemble de l’excitation sismique.

On peut encore décomposer le problème en introduisant l’inconnue \({d}_{ij}(t)\) solution de l’équation différentielle: \({\ddot{\mathrm{d}}}_{ij}+2{{\xi}_{i}{\omega}_{i}\dot{\mathrm{d}}}_{ij}+{\omega}_{i}^{2}{\mathrm{d}}_{ij}={g}_{j}(t)\) , cette dernière équation correspond à la réponse du mode dynamique \(i\) à l’accélération

. Le déplacement relatif sur la base physique s’exprime alors:

\({X}_{r}(t)=-\sum_{i=1}^{n}\sum_{j=1}^{m}{p}_{ij}{d}_{ij}(t){\varphi}_{i}\)

L’information sur la position du point d’appui est contenue dans le facteur de participation modale.

Réponse harmonique

On a donc décomposé la réponse totale de la structure en une contribution relative et une contribution différentielle due aux déplacements des ancrages telles que :

\(\begin{array}{}{X}_{a}(t)={X}_{r}(t)+{X}_{e}(t)\\ \text{avec}\\ \lbrace \begin{array}{}{\ddot{X}}_{e}(t)={\Psi \ddot{X}}_{s}(t)=\sum_{j=1}^{m}{\psi}_{j}{g}_{j}(t)\\ {X}_{r}(t)\text{= -}\sum_{i=1}^{n}\sum_{j=1}^{m}{{p}_{ij}d}_{ij}(t){\varphi}_{i}\text{}\text{où}{d}_{ij(t)}\text{est solution de}\text{}{\ddot{d}}_{ij}+2{\xi}_{i}{\omega}_{i}{\dot{d}}_{ij}+{\omega}_{i}^{2}{d}_{ij}={g}_{j}(t)\end{array}\end{array}\)

La résolution de cette dernière équation différentielle par la méthode de la transformation de Fourier fait intervenir les fonctions de transfert modales \({h}_{i}(\omega )\) telles que : . \({h}_{i}(\omega )=\frac{1}{{\omega}_{i}^{2}-{\omega}^{2}+2i{\xi}_{i}{\omega}_{i}\omega }\)

On obtient donc : \({d}_{ij}(\omega )={h}_{i}(\omega ).{g}_{j}(\omega )\text{et}{\ddot{d}}_{ij}(\omega )\text{=-}{\omega}^{2}{h}_{i}(\omega ).{g}_{j}(\omega )\)

La réponse harmonique totale de la structure se déduit des formules précédentes par recombinaison modale.

\(\begin{array}{}{\ddot{X}}_{a}(\omega )={\ddot{X}}_{r}(\omega )+{\ddot{X}}_{e}(\omega )\\ {\ddot{X}}_{a}(\omega )={\omega}^{2}\sum_{i=1}^{n}\sum_{j=1}^{m}{p}_{ij}{h}_{j}(\omega ){\varphi}_{j}(\omega ){j}_{i}+\sum_{j\text{=1}}^{m}{\Psi}_{j}{g}_{j}(\omega )\end{array}\)

On fait alors apparaitre la matrice complexe \((n,m)\) , dite matrice de transfert \(H(\omega )\) suivante:

\(H(\omega )={\omega}^{2}P.h(\omega )\Phi +\Psi\)

où \(P\) est la matrice des facteurs de participation, \(h(\omega )\) le vecteur des fonctions de transfert modales \({h}_{i}(\omega )\) .

La réponse totale de la structure vaut \({\ddot{X}}_{a}(\omega )=H(\omega )\ddot{E}(\omega )\) où \(\ddot{E}(\omega )\) est le vecteur de m lignes constitué des transformées de Fourier des accélérations

aux \(m\) degrés de liberté-appuis.

On voit que cette expression détermine la réponse en accélération. Ceci impose ensuite d’intégrer deux fois la réponse pour obtenir le déplacement, ce problème est présenté dans [bib4]. Un des intérêts supplémentaires de la méthode que nous proposons ici est de s’abstraire de cette difficulté.

La réponse dynamique aléatoire

Rappel sur les densités spectrales de puissance [bib2]

Définitions

Soit un signal probabiliste défini par sa densité de probabilité \({p}_{x}({x}_{1},{t}_{1},\dots ,{x}_{n},{t}_{n})\) . Cette densité de probabilité permet de calculer les fonctions moments du signal.

Moment d’ordre 1 ou espérance du signal:

\({\mu}_{X}(t)=E\left[X(t)\right]=\underset{-\infty }{\overset{+\infty }{\int}}X{p}_{x}(x,t)\text{dx}\)

Moments d’ordre 2 ou intercorrélation de deux signaux :

\({\rho}_{\text{XY}}({t}_{1},{t}_{2})=E\left[X({t}_{1})\overline{Y({t}_{2})}\right]=\underset{-\infty }{\overset{+\infty }{\int}}x\stackrel{ˉ}{y}p(x,{t}_{1};y,{t}_{2})\text{dxdy}\)

Lorsque le signal est stationnaire, l’intercorrélation ne dépend que de \(\tau ={t}_{2}-{t}_{1}\) .

On l’écrit \({R}_{\text{XY}}(t)=E\left[X(t)\overline{Y(t-\tau )}\right]\)

Densité spectrale de puissance et interspectre

On définit \({S}_{\text{XY}}(\omega )\) l’interspectre de puissance ou densité interspectrale de puissance entre deux signaux probabilistes stationnaires par la transformée de Fourier de la fonction d’intercorrélation, ce qu’on écrit :

\({S}_{\text{XY}}(\omega )=\frac{1}{2\pi }\underset{-\infty }{\overset{+\infty }{\int}}{R}_{\text{XY}}(\tau ){e}^{-i\omega \tau }d\tau\)

La formule inverse s’écrit : \({R}_{\text{XY}}(t)=\underset{-\infty }{\overset{+\infty }{\int}}{S}_{\text{XY}}(\omega ){e}^{i\omega \tau }d\omega\)

\({S}_{\text{XY}}(\omega )\) est généralement complexe et vérifie la relation de symétrie : \({S}_{\text{YX}}(\omega )=\overline{{S}_{\text{XY}}(\omega )}\) .

Lorsque \(X=Y\) , \({S}_{\text{XX}}(\omega )\) s’appelle autospectre de puissance ou densité spectrale de puissance(DSP) . Cette fonction a la propriété d’être réelle et toujours positive.

Relations entre la DSP et les autres caractéristiques du signal

Remarque:

La plupart du temps, le signal est défini sur un temps limité, sa transformée de Fourier n’existe pas, on définit alors une transformée de Fourier estimée sur une période de longueur \(T\) par:

\({\stackrel{ˆ}{X}}_{T}(\omega )=\frac{1}{2\pi }\underset{-T/2}{\overset{T/2}{\int}}X(t){e}^{-i\omega \tau }d\tau\) .

On a alors les relations suivantes avec cette transformée de Fourier estimée :

\(\begin{array}{}{S}_{\text{XY}}(\omega )=\underset{T\to +\infty }{\lim}\frac{2\pi }{T}E\left[{\stackrel{ˆ}{X}}_{T}(w)\overline{{\stackrel{ˆ}{Y}}_{T}(\omega )}\right]\\ {S}_{\text{XX}}(\omega )=\underset{T\to +\infty }{\lim}\frac{2\pi }{T}E\left[{\stackrel{ˆ}{X}}_{T}(w)\overline{{\stackrel{ˆ}{X}}_{T}(\omega )}\right]\end{array}\)

Lien entre l’autospectre de puissance et la puissance du signal :

La puissance d’un signal est égale à sa variance. Pour un signal centré, la variance vaut : \({\sigma}_{X}^{2}={R}_{\text{XX}}(0)\) .

On a donc : \({\sigma}_{X}^{2}={R}_{\text{XX}}(0)=\underset{-\infty }{\overset{+\infty }{\int}}{S}_{\text{XX}}(\omega )d\omega\) .

Les équations de mouvement

La réponse totale de la structure est déterminée par la relation : \({\ddot{X}}_{a}(\omega )=H(\omega )\ddot{E}(\omega )\) ,

où \(\ddot{E}(\omega )\) est le vecteur de \(m\) lignes constitué des excitations représentées par les transformées de Fourier des accélérogrammes \({g}_{j}(t)\) aux \(m\) degrés de liberté-appuis,

\(H(\omega )\) est la matrice de transfert définie par \(H(\omega )={\omega}^{2}\text{ph}(\omega )\Phi +\Psi\)

où \(p\) est la matrice des facteurs de participation,

\(H(\omega )\) le vecteur des fonctions de transfert modales \({h}_{i}(\omega )\)

\(\Phi\) base des modes dynamiques

\(\Psi\) base des modes statiques

elle comporte \(n\) lignes (= nombre de degrés de liberté libres de la structures) et \(m\) colonnes.

Matrice «interspectrale-excitation»

NB :

Cette appellation « matrice interspectrale-excitation » est abusive : elle signifie « matrice de densité interspectrale de puissance de l’excitation ».

On suppose que l’excitation sismique peut être considérée comme un signal stationnaire - compte-tenu des rapports entre les temps représentatifs - et centré. Ceci permet d’utiliser un certain nombre de résultat de l’analyse probabiliste. On s’intéresse alors à la réponse stationnaire du système à une excitation stationnaire.

On note \({S}_{\ddot{E}\ddot{E}}(\omega )\) la matrice des interspectres de puissance correspondant à l’excitation. Sa donnée est explicitée dans le chapitre 4.

Pour mémoire nous rappelons ici qu’elle est calculée à partir de transformées de Fourier des accélérations. C’est une matrice (\(m\times m\) ). Le terme ij correspond à l’interspectre entre les signaux \({\ddot{E}}^{i}\) et \({\ddot{E}}^{j}\) soit encore entre les transformées de Fourier des accélérogrammes \({g}_{i}\) et \({g}_{j}\) .

Réponse dynamique aléatoire

On a vu que l’interspectre de puissance entre deux signaux probabilistes est la transformée de Fourier de la fonction d’intercorrélation des deux signaux. On l’applique à la réponse totale de la structure:

\({S}_{{{\ddot{X}}_{a}\ddot{X}}_{a}}(\omega )=\frac{1}{2\pi }\underset{-\infty }{\overset{+\infty }{\int}}{R}_{{{\ddot{X}}_{a}\ddot{X}}_{a}}(\tau ){e}^{-i\omega \tau }d\tau =\frac{1}{2\pi }\underset{-\infty }{\overset{+\infty }{\int}}E\left[{\ddot{X}}_{a}(t){}^{T}\text{}\overline{{\ddot{X}}_{a}(t-\tau )}\right]{e}^{-i\omega \tau }d\tau\)

On travaille alors dans le domaine temporel pour exprimer la fonction d’intercorrélation de la réponse totale \({R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(t,t')\) .

On note \(h(t)\) la réponse impulsionnelle du système : \(h(t)={\text{TF}}^{-1}[H(\omega )]\)

et \(\ddot{e}(t)\) la transformée de Fourier inverse de la DSP excitatrice : \(\ddot{e}(t)={\text{TF}}^{-1}[\ddot{E}(\omega )]\)

Par transformée de Fourier inverse de la relation : \({\ddot{X}}_{a}(\omega )=H(\omega )\ddot{E}(\omega )\)

on a \({\ddot{\mathrm{X}}}_{\mathrm{a}}(t)=\mathrm{h}\times \ddot{\mathrm{e}}(t)=\underset{R}{\int}\mathrm{h}(u)\ddot{\mathrm{e}}(t-u)\text{du}\)

\(\begin{array}{}{R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(t,t')=E\left[{\ddot{X}}_{a}(t){}^{T}\text{}\overline{{\ddot{X}}_{a}(t')}\right]\\ {R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(t,t')=E\left[\underset{R}{\int}h(u)\ddot{e}(t-u)\text{du}{}^{T}\text{}\overline{\underset{R}{\int}h(v)\ddot{e}(t'-v)\text{dv}}\right]\\ {R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(t,t')=E\left[\underset{R}{\int}\underset{R}{\int}h(u)\ddot{e}(t-u){}^{T}\text{}\overline{\ddot{e}(t'-v)}{}^{T}\text{}\overline{h(v)}\text{dv}\text{du}\right]\end{array}\)

On suppose dans cette analyse le système déterministe, on peut donc sortir la réponse impulsionnelle du calcul de l’espérance mathématique. Il vient :

\({R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(t,t')=\underset{R}{\int}\underset{R}{\int}h(u)E[\ddot{e}(t-u){}^{T}\text{}\overline{\ddot{e}(t'-v)}]{}^{T}\text{}\overline{h(v)}\text{dv}\text{du}\)

L’excitation est supposée un processus stationnaire, l’intercorrélation ne dépend donc que de l’écart de temps \(\tau =t-t'\) :

\({R}_{\ddot{E}\ddot{E}}(t-t'-u+v)=E[\ddot{e}(t-u){}^{T}\text{}\overline{\ddot{e}(t'-v)}]={R}_{\ddot{E}\ddot{E}}(\theta )\text{pour}\theta =t-t'-u+v=\tau -u+v\)

d’où \({R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(t,t')=\underset{R}{\int}\underset{R}{\int}h(u){R}_{\ddot{E}\ddot{E}}(\theta ){}^{T}\text{}\overline{h(v)}\text{dv}\text{du}={R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(\tau )\) ce qui justifie a posteriori l’approche.

On reporte maintenant cette expression dans l’expression de la densité spectrale de puissance de la réponse:

\({S}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(\omega )=\frac{1}{2\pi }\underset{-\infty }{\overset{+\infty }{\int}}{R}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}(\tau ){e}^{-i\omega \tau }d\tau =\frac{1}{2\pi }\underset{-\infty }{\overset{+\infty }{\int}}\underset{R}{\int}\underset{R}{\int}h(u){R}_{\ddot{E}\ddot{E}}(\tau -u+v){}^{T}\text{}\overline{h(v)}{e}^{-i\omega \tau }\text{dv}\text{du}\mathrm{dt}\)

En répartissant les variables muettes d’intégration on fait apparaître les transformées de Fourier respectives de \(h(u),{R}_{\ddot{E}\ddot{E}}(t-u+v),{}^{T}\text{}\overline{h(v)}\) , il vient finalement:

\({S}_{{\ddot{X}}_{a}{\ddot{X}}_{a}}=H(\omega ).{S}_{\ddot{E}\ddot{E}}(\omega ).\overline{{}^{T}\text{}H(\omega )}\)

avec \(H(\omega )={\omega}^{2}\text{P.h}(\omega )\Phi +\Psi\)

Compte-tenu des relations entre les transformées de Fourier du déplacement, de la vitesse et de l’accélération, on a de plus:

\(\begin{array}{}{S}_{{\dot{X}}_{a}{\dot{X}}_{a}}=\frac{-1}{{\omega}^{2}}H(\omega ){S}_{\ddot{E}\ddot{E}}(\omega )\overline{{}^{T}\text{}H(\omega )}\\ {S}_{{X}_{a}{X}_{a}}=\frac{1}{{\omega}^{4}}H(w){S}_{\ddot{E}\ddot{E}}(\omega )\overline{{}^{T}\text{}H(\omega )}\end{array}\)

Ces relations permettent d’exprimer la réponse de la structure par la DSP du déplacement ou de la vitesse.

Remarques :

• Selon l’expression donnée à \(H(\omega )\) , on exprime la DSP du déplacement (respectivement de la vitesse ou de l’accélération) total, relatif ou différentiel :

mouvement absolu : \(H(\omega )={\omega}^{2}\text{P. h}(\omega )\Phi +\Psi\)

mouvement relatif : \(H(\omega )={\omega}^{2}\text{P. h}(\omega )\Phi \text{}\)

mouvement différentiel (ie d’entrainement) : \(H(\omega )=\Psi\)

• Il est d’usage, lors d’un calcul avec Code_Aster, de restreindre la matrice de la fonction de transfert aux lignes des \(l\) degrés de liberté d’observation. Ceci permet d’alléger d’autant les calculs dès que \(l\) est petit devant \(n\) .

Application dans Code_Aster

L’ensemble de l’approche spectrale pour le calcul sismique est traité dans la commande DYNA_ALEA_MODAL [U4.53.22]. Les données sont regroupées sous trois mots clé facteurs et un mot clé simple.

La base modale est constituée des modes dynamiques calculés par la commande CALC_MODES [U4.52.02] stockés dans un concept de type mode_meca récupéré par le mot-clé facteur BASE_MODALE, d’une part ; des modes statiques calculés par la commande MODE_STATIQUE [U4.52.14] stockés dans un concept de type mode_stat récupéré par le mot‑clé simple MODE_STAT, d’autre part. Le mot clé facteur BASE_MODALE possède également les arguments qui permettent de déterminer la bande de fréquence ou les modes retenus pour le calcul et les amortissements correspondants.

Les données correspondant à l’excitation sont rassemblées sous le mot-clé facteur EXCIT (cf.paragraphe [§4]) : on y précise le type d’excitation au sens de la GRANDEUR : excitation en déplacement ou en effort, le ou les nœuds NOEUD et composantes NOM_CMP excités, le nom des interspectres ou autospectres INTE_SPEC, fonctions complexes préalablement lues ou calculées, respectivement par les opérateurs LIRE_INTE_SPEC [U4.36.01] ou CALC_INTE_SPEC [U4.36.03] et stockés dans une table d’interspectre de concept tabl_intsp qui s’appliquent en chaque degré de liberté excité.

Sous le mot-clé facteur REPONSE se trouvent les données liées au choix de la discrétisation.

La commande DYNA_ALEA_MODAL fournit la réponse sous forme de densité spectrale de puissance sur base modale. Pour obtenir la restitution des DSP sur base physique, on utilisera REST_SPEC_PHYS [U4.63.22] qui permet de préciser le type de grandeur de la réponse (déplacement ou effort), aux “points d’observation” (nœud-composante) du résultat. En présence d’une réponse de type déplacement, on précisera ici aussi si la réponse correspond au déplacement absolu, relatif ou différentiel.

REST_SPEC_PHYS fournit une table d’interspectres qui contient selon la demande de l’utilisateur, la matrice interspectrale en déplacement \({S}_{\text{XX}}\) , en vitesse \({S}_{\dot{X}\dot{X}}\) , ou en accélération \({S}_{\ddot{X}\ddot{X}}\) pour une expression dans le repère absolu (indice \(a\) ), le repère relatif (indice \(r\) ) ou d’entrainement (indice \(e\) ).

Chaque « combinaison » précédente nécessite un appel spécifique à la commande REST_SPEC_PHYS.

Définition de la matrice interspectrale de puissance excitatrice

L’excitation sismique est par nature, nous l’avons dit, aléatoire. Aussi elle peut être connue non pas par son expression temporelle mais sous forme fréquentielle par une densité spectrale de puissance dite aussi interspectre.

Lorsqu’il y a plusieurs appuis, ils peuvent être excités par des excitations identiques ou différentes, ce dernier cas est celui du multi-appuis.

Pour \(m\) appuis, on définit la matrice de densité interspectrale de puissance d’ordre \(m\) , ou par abus de langage l’interspectre d’ordre \(m\) , qui est une matrice (\(m\times m\) ) de fonctions complexes dépendant de la fréquence.

Les termes diagonaux représentent les densités “auto-” spectrales de puissances -ouautospectres- aux points d’excitation, les termes extra-diagonaux correspondent aux densités interspectrales entre les excitations en deux points d’appui distincts (chaque ligne ou colonne de la matrice représente en fait un point d’appui en maillage physique ou un mode en calcul modal). Par définition de ces termes, il s’en déduit que les matrices de densité interspectrales de puissance manipulées sont hermitiennes. (Voir [bib2] ou documentation de Référence associée à la commande POST_DYNA_ALEA [R7.10.01])

Nous présentons ci-après les différentes commandes de Code_Aster qui permettent d’obtenir une matrice de densité interspectrale de puissance.

Lecture sur un fichier

La façon la plus élémentaire de définir une matrice de densité interspectrale de puissance est de donner, “à la main”, les valeurs aux différents pas de fréquence.

On utilise alors l’opérateur LIRE_INTE_SPEC [U4.36.01].

LIRE_INTE_SPEC lit dans un fichier “interspectre excitation”. Le format du fichier dans lequel est consignée la matrice interspectrale est simple : on décrit successivement la fonction de chaque terme de la matrice interspectrale ; pour chaque fonction, on donne une ligne par fréquence en indiquant la fréquence, les parties réelle et imaginaire du nombre complexe; ou la fréquence, le module et la phase du nombre complexe (mot clé FORMAT).

Exemple de fichier interspectre excitation (pour une matrice réduite à un terme) :

INTERSPECTRE DIM = 1 FONCTION_C I = 1 J = 1 NB_POIN = 4 VALEUR = 2.9999 0. 0. 3. 1. 0. 13. 1. 0. 13.0001 0. 0. FINSF FIN

Obtention d’un interspectre à partir de fonctions du temps

On peut déduire la matrice de densité interspectrale de puissance à partir de fonctions du temps. On utilise alors l’opérateur CALC_INTE_SPEC [U4.36.03] dans Code_Aster [bib3].

A partir d’une liste de \(N\) fonctions du temps, cet opérateur permet de calculer l’interspectre de puissance \(N\times N\) qui leur correspond.

Pour chaque terme de la matrice interspectrale (\(N\times N\) ) on utilise la démarche suivante [bib3].

Pour calculer l’interspectre de deux signaux on utilise la relation de Wiener-Khintchine [bib7] qui permet d’établir une formule de calcul de la densité spectrale de puissance par la transformée de Fourier d’échantillons finis des signaux \(x(t)\) et \(y(t)\) .

Il vient alors:

\({S}_{xy}(f)=\underset{T\to \infty }{\lim}\frac{1}{T}\mathrm{E}[{\mathrm{X}}_{k}(f,T).{\mathrm{Y}}_{k}\times (f,T)]\)

où \(\lbrace \begin{array}{c}{X}_{k}(f,T)=\text{TF}[{x}_{k}](f)=\underset{0}{\overset{T}{\int}}{x}_{k}(t){e}^{-\mathrm{i2}\pi f}\text{dt}\\ {Y}_{k}(f,T)=\text{TF}[{y}_{k}](f)=\underset{0}{\overset{T}{\int}}{y}_{k}(t){e}^{-\mathrm{i2}\pi f}\text{dt}\end{array}\)

sont les transformées de Fourier discrètes de \(«x»\) et \(«y»\) .

Lorsque l’on s’intéresse à des signaux issus de mesures, on ne dispose la plupart du temps que de signaux connus de façon discrète, de même un résultat de calcul transitoire est un signal discret.

Une approximation de l’interspectre des signaux discrets \(x[n]\) et \(y[n]\) définis sur \(L\) points espacés de \(\Delta t\) , découpés en \(p\) blocs de \(q\) points est obtenue par la relation:

\(\begin{array}{}{\stackrel{ˆ}{S}}_{xy}\left[k\right]=\frac{1}{\mathrm{pq}\Delta t}\sum_{i=1}^{p}{X}^{(i)}[k]{Y}^{(i)}\ast [k]\\ {X}^{(i)}[k]=\Delta t\sum_{n=0}^{q}{x}^{(i)}[n]{e}^{-\mathrm{2i}\pi \text{kn}/q}\\ {Y}^{(i)}[k]=\Delta t\sum_{n=0}^{q}{y}^{(i)}[n]{e}^{-\mathrm{2i}\pi \text{kn}/q}\end{array}\)

Les différents blocs peuvent ou non se recouvrir. Les valeurs \(p\) et \(q\) sont au choix de l’utilisateur.

Cette méthode est celle du périodogramme de WELCH [bib8].

Le calcul se fait sur une fenêtre qui se déplace sur le domaine de définition des fonctions. L’utilisateur spécifie dans la commande la longueur de la fenêtre d’analyse, le décalage entre deux fenêtres de calcul successives et le nombre de points par fenêtre.

Excitations prédéfinies ou reconstituées à partir de fonctions complexes existantes

On peut souhaiter définir une matrice de densité interspectrale de puissance de différentes façons:

• par un bruit blanc : les valeurs sont constantes

• selon la formule analytique de KANAI-TAJIMI utile en calcul sismique (bruit blanc filtré),

• ou en reprenant des fonctions complexes existantes.

On utilise alors l’opérateur DEFI_INTE_SPEC [U4.36.02].

Fonctions complexes existantes

Il suffit sous le mot-clé facteur PAR_FONCTION de donner le nom de la fonction pour chaque paire d’indice NUME_ORDRE_I, NUME_ORDRE_J, correspondant à la matrice triangulaire supérieure (en raison de son hermiticité).

Bruit blanc

Un bruit blanc se caractérise par une valeur constante sur tout le domaine de définition considéré. Sous le mot-clé facteur CONSTANT, on donne cette valeur (VALE_R ou VALE_C) sur la bande de fréquence [FREQ_MIN, FREQ_MAX] pour chaque paire d’indice INDI_I, INDI_J, correspondant à la matrice triangulaire supérieure (en raison de son hermiticité). Pour définir parfaitement la fonction, on précise l’interpolation et les prolongements.

Bruit blanc filtré par KANAI-TAJIMI [bib9]

Pour une structure appuyée sur le sol, il est courant de prendre comme excitation la densité spectrale de puissance de Kanaï-Tajimi. Cette densité spectrale représente le filtrage d’un bruit blanc par le sol, considéré comme un système à un degré de liberté. Les paramètres de la formule permettent de jouer sur la fréquence centrale et la largeur de bande du spectre.

Le spectre \(G(\omega )\) s’exprime par la relation suivante:

\(\begin{array}{}G(\omega )=\frac{{\omega}_{g}^{4}+4{\xi}_{g}^{2}{\omega}_{g}^{2}{\omega}^{2}}{{({\omega}_{g}^{2}-{\omega}^{2})}^{2}+4{\xi}_{g}^{2}{\omega}_{g}^{2}{\omega}^{2}}{G}_{0}\\ \begin{array}{cc}{\omega}_{g}=2\pi f& \text{pulsation propre}\\ {\xi}_{g}& \mathrm{amortissement}\mathrm{total}\\ {G}_{0}& \mathrm{niveau}\mathrm{du}\mathrm{bruit}\mathrm{blanc}\mathrm{avant}\mathrm{filtrage}\end{array}\end{array}\)

L’utilisateur doit spécifier la fréquence propre \({f}_{g}\) du filtre, l’amortissement modal \({\xi}_{g}\) et le niveau du bruit blanc \({G}_{0}\) (= VALE_R) avant le filtrage ; ainsi que comme pour toute fonction : l’interpolation, les profils extérieurs et le domaine de définition (bande de fréquence).

Par défaut un sol courant est bien représenté par les valeurs \({f}_{g}=2.5\mathrm{Hz}\) et \({\xi}_{g}=0.6\) .

Exemple d’utilisation pour un bruit blanc filtré par KANAI_TAJIMI :

Interex =        DEFI_INTE_SPEC (

DIMENSION : 1 KANAI_TAJIMI : ( NUME_ORDRE_I : 1 indices du terme de la matrice de densité NUME_ORDRE_J : 1 interspectrale de puissance FREQ_MOY : 2.5 fréquence propre AMOR : 0.6 amortissement modal VALE_R : 1 niveau du bruit blanc INTERPOL : ‘LIN’ interpolation linéaire PROL_GAUCHE : ‘CONSTANT’ prolongation PROL_DROIT : ‘ CONSTANT’ FREQ_MIN : 0. domaine de définition FREQ_MAX : 200. PAS : 1. ) ) ;

Autres types d’excitation

Les calculs des paragraphes précédents ont été effectués dans le cadre de l’hypothèse d’une excitation en mouvement imposé sur un degré de liberté. Moyennant quelques modifications il est possible d’utiliser la même approche pour une excitation en effort [§4.4.1] ou par des sources fluides [§ 4.4.2], celle-ci étant exprimée dans un élément fini [§4.4.3] ou sur une fonction de forme de la structure [§4.4.4].

Dans la suite de ce paragraphe, on suppose l’excitation aléatoire connue et fournie par l’utilisateur sous le forme d’une DSP, densité spectrale de puissance.

Cas de l’excitation en forces imposées

Sous le mot-clé EXCIT on a GRANDEUR = EFFO.

Quand l’excitation aux appuis est de type force imposée, l’équation générale du mouvement est:

\(M\ddot{X}(t)+C\dot{X}(t)+\text{KX}(t)=\sum_{j=1}^{m}{F}_{j}\)

La réponse de la structure est alors calculée sur une base de modes dynamiques \(\Phi =\left\lbrace {\phi }_{i,i=1,n}\right\rbrace\) , ces modes étant calculés en supposant les appuis excitateurs libres. On ne distingue pas, dans ce cas de mouvement absolu, relatif et différentiel et on n’utilise pas de modes statiques.

On définit le facteur de participation modale sous la forme : \({P}_{ij}=\frac{{}^{T}\text{}{\phi }_{i}{F}_{j}}{{\mu}_{i}}\)

Les réponses transitoires, harmoniques et aléatoires ont les mêmes expressions que les réponses du mouvement relatif de l’excitation multi-appui dans le cas général [§3]. (Ce qui correspond à l’absence de modes statiques). La force excitatrice est représentée en chaque degré de liberté-appui par sa DSP sous forme d’un terme équivalent à \({S}_{\ddot{E}\ddot{E}}(\omega )\) .

Excitation par des sources fluides

Les sources fluides apparaissent, par exemple, dans l’étude d’un réseau de tuyauteries. Elles correspondent à des organes actifs ou à des branchements de tuyauteries secondaires. Elles sont généralement des sources de pression ou des sources de débit. Ces différents types de source sont présentés ci-après en fonction de leur mise en forme mathématique et de ce que fait Code_Aster dans chaque configuration.

Ces sources fluides ne sont pas directement des excitations sismiques mais peuvent être induites par un séisme. La résolution du problème mécanique fait appel aux même méthodes, du fait de leur caractère aléatoire, ce qui justifie leur présentation ici.

La modélisation du réseau de tuyauterie est supposée réalisée à l’aide de poutre vibro-acoustique de Code_Aster .

La réponse à des sources fluides se calcule dans le cadre de la réponse à des forces imposées (cf.[§4.4.1]), dans ce cadre on s’intéresse à des réponses de grandeur de type “déplacement” (GRANDEUR = DEPL_R sous le mot-clé REPONSE).

Les sources de pression et de force , pour des raisons de modélisation des sources fluides sont représentées par des dipôles [bib5], il est donc nécessaire de donner deux points d’application .

Source de débit-volume : GRANDEUR = SOUR_DEBI_VOLU sous le mot-clé EXCIT

Un débit volume s’exprime en \({m}^{3}/s\) , sa densité spectrale de puissance en \({({m}^{3}/s)}^{2}/\mathrm{Hz}\) .

Une source de débit-volume est considérée, dans la formulation \(P-\phi\) des éléments de tuyauterie avec fluide, comme un effort imposé sur le degré de liberté \(\phi\) du nœud d’application de la source [R4.02.02].

L’utilisateur fournit la DSP de débit volume \({S}_{\text{vv}}(\omega )\) , la DSP \({S}_{\text{vv}}^{'}(\omega )\) appliquée en effort sur le degré de liberté \(\phi\) est : \({S}_{\text{vv}}^{'}(\omega )={(\rho \omega )}^{2}{S}_{\text{vv}}(\omega )\)

où \(\rho\) est la masse volumique du fluide.

Source de débit-masse : GRANDEUR = SOUR_DEBI_MASS sous le mot-clé EXCIT

Un débit-masse s’exprime en \(\mathit{kg}/s\) , sa densité spectrale de puissance en \({(\mathrm{kg}/s)}^{2}/\mathrm{Hz}\) . Le débit‑masse est le produit du débit-volume par la masse volumique du fluide.

L’utilisateur fournit la DSP de débit-masse \({S}_{\text{mm}}(\omega )\) , la DSP \({S'}_{\text{mm}}(\omega )\) appliqué en effort sur le degré de liberté \(\phi\) est : \({S'}_{\text{mm}}\omega ={\omega}^{2}{S}_{\text{mm}}\omega\)

Source de pression : GRANDEUR = SOUR_PRESS sous le mot-clé EXCIT

Une source de pression est appliquée dans Aster en un dipôle \({P}_{1}{P}_{2}\) .

Pour une source de pression dont la DSP est \({S}_{\text{PP}}(\omega )\) , exprimée en \({\mathrm{Pa}}^{2}/\mathrm{Hz}\) , Aster construit une matrice de densité interspectrale de puissance \({S'}_{\text{pp}}(\omega )\) qui est appliquée en force imposée sur le degré de liberté \(\phi\) des points \({P}_{1}\) et \({P}_{2}\) .

\({S'}_{\text{PP}}(\omega )={S}_{\text{PP}}(\omega )(\begin{array}{cc}{(\frac{S}{\text{dx}})}^{2}& -{(\frac{S}{\text{dx}})}^{2}\\ -{(\frac{S}{\text{dx}})}^{2}& {(\frac{S}{\text{dx}})}^{2}\end{array})\)

où \(S\) est la section fluide, \(\text{dx}\) la distance entre les deux points \({P}_{1}\) et \({P}_{2}\) .

Source de force : GRANDEUR = SOUR_FORCE sous le mot-clé EXCIT

La force correspond simplement au produit de la pression par la section fluide du tube : \(F=\mathrm{PS}\) . Elle est donc aussi appliquée sur un dipôle \({P}_{1}{P}_{2}\) .

Pour une source de force dont la DSP est \({S}_{\text{FF}}(\omega )\) , exprimée en \({N}^{2}/\mathrm{Hz}\) , Aster applique en force imposée sur le degré de liberté \(f\) des points \({P}_{1}\) et \({P}_{2}\) , (distants de \(\mathrm{dx}\) ), la matrice de densité interspectrale de puissance \({S'}_{\text{FF}}(\omega )\) telle que:

\({S'}_{\text{FF}}(\omega )={S}_{\text{FF}}(\omega )(\begin{array}{cc}{(\frac{1}{\text{dx}})}^{2}& -{(\frac{1}{\text{dx}})}^{2}\\ -{(\frac{1}{\text{dx}})}^{2}& {(\frac{1}{\text{dx}})}^{2}\end{array})\)

Excitation répartie sur une fonction de forme

Si la densité spectrale de puissance de l’excitation \(E(\omega )\) correspond à un effort imposé sur une fonction de forme \({f}_{i}\) , \(E(\omega )\) donne la dépendance fréquentielle du niveau de l’excitation.

La pondération spatiale de l’effort est représentée dans Code_Aster par un champ aux nœuds qui ne dépend pas de la fréquence : mot clé CHAM_NO sous le mot-clé facteur EXCIT. Ce champ aux nœuds est un “vecteur assemblé”. Du point de vue théorique le formalisme de calcul est le même que précédemment (excitation en force imposée [§4.4.1]), pour un vecteur de force en second membre égal à \({f}_{i}\) .

Applications

Ces différents types d’excitation sont repris dans les tests de validation, et sont présentés par des exemples dans le rapport [bib6]. En particulier les excitations de type fluide sont dans le test : tuyau soumis à des excitations fluides aléatoires [V2.02.105] (SDLL105). Les excitations sur des fonctions de forme sont testées dans le cas test : poutre soumise à une excitation aléatoire répartie [V2.02.106] (SDLL106).

Bibliographie

1. 16. LABBE et H. NOE : “Stochastic approach for the seismic design of nuclear power plant equipements”. Nuclear Engineering and Design 129 (1991) 367‑379.

2. 1. DUMOND Rapport EDF DER HP62/95.021B : Post traitement d’un calcul de mécanique vibratoire sous excitations aléatoires dans Code_Aster . Note de référence de la commande POST_DYNA_ALEA.

3. 7. JACQUART Rapport EDF DER HP61/93.073 : Générations de signaux aléatoires de densité spectrale donnée : Note de principe et cahier des charges de l’intégration à ASTER.

4. Fe. WAECKEL Rapport EDF DER HP62/95.017B : Méthode pour le calcul par superposition modale de la réponse sismique d’une structure multi-supportée.

5. 16. THOMAS : Prise en compte des sources acoustiques dans les modèles de tuyauteries en mécanique. Bulletin de la DER - série A. Nucléaire Hydraulique Thermique n° 2 1991 pp19-36.

6. 3. DUVAL Rapport EDF DER HP-61/92.148 : Réponse dynamique sous excitations aléatoires dans Code_Aster : principes théoriques et exemples d’utilisation.

7. BENDAT et PIERSOL : Engineering applications of correlation and spectral analysis. John Wiley and Son 1980.

8. MARPLE : Digital spectral analysus with applications. Prentice Hall 1987.

9. 8. TAJIMI A statistical method of determining the maximum response of a building structure during an earthquake. Proc 2nd world Conf. Earthquake Eng. Tokyo and Kyoto, Japan (1960) pp 751-797.

Description des versions du document

Version Aster	Auteur(s) Organisme(s)	Description des modifications
4	A.Dumond EDF-R&D/AMA	Texte initial
14	F.Voldoire EDF-R&D/ERMES	Petites corrections