r4.02.01 Éléments finis en acoustique#

Résumé:

Ce document décrit en acoustique linéaire stationnaire à basse fréquence en milieu borné les équations utilisées, les formulations variationnelles qui en découlent ainsi que la traduction correspondante en éléments finis isoparamétriques à l’aide d’une formulation usuelle à une inconnue \(p\) (pression acoustique). Ce document correspond à PHENOMENE=”ACOUSTIQUE”.

Équations et conditions aux limites du problème#

Équations et conditions aux limites#

L’équation à résoudre en domaine fréquentiel est l’équation de Helmholtz [bib2]:

(898)#\[(\Delta +{k}^{2})p=s\]

\(p\) est une grandeur complexe désignant la pression acoustique et \(s\) , également complexe, représente les termes sources du problème. Le paramètre \(k\) désigne le nombre d’onde du problème traité; il peut être complexe ou réel, suivant que la propagation s’effectue ou non dans un domaine comportant une certaine atténuation modélisée par le biais d’une célérité complexe, par exemple si le matériau est «poreux» [bib6]:

(899)#\[k=\frac{\omega}{c}\]

avec \(c\) désignant la célérité du son, qui peut être complexe dans le cas d’une propagation en milieu «poreux» ou viscoacoustique et \(\omega\) est un réel dans tous les cas, qui désigne la pulsation telle que:

(900)#\[\omega =2\pi f\]

\(f\) est la fréquence de travail du problème harmonique. Nous représentons sur la figure le domaine confiné quelconque où s’applique l’équation de Helmholtz () et les conditions aux frontières. \(\Omega\) est un ouvert borné de \({ℝ}^{3}\) de frontière \(\partial \Omega\) régulière, partitionnée en \(\partial {\Omega}_{v}\) et \(\partial {\Omega}_{z}\) :

(901)#\[\partial \Omega =\partial {\Omega}_{v}\cup \partial {\Omega}_{z}\]
../../../../_images/100010C00000263E000012B587364503D60EC012.svg

Figure 2.1-1: Configuration du problème acoustique en domaine borné

L’équation () est à résoudre dans un domaine borné \(\Omega\) . Les conditions aux limites à prendre en compte sur la frontière \(\partial \Omega\) du domaine \(\Omega\) s’expriment sous leur forme la plus générale comme:

(902)#\[\alpha p+\beta \frac{\partial p}{\partial n}=\gamma\]

Classiquement, \(\partial /\partial n\) désigne l’opérateur de dérivée normale. \(\alpha\) , \(\beta\) et \(\gamma\) sont des opérateurs complexes, qui peuvent être des scalaires, ou des opérateurs intégraux suivant que la frontière d’application de la condition à la limite est à réaction locale ou à réaction non locale (cas de l’interaction fluide-structure).

Les développements actuellement réalisés dans Code_Aster concernent uniquement des conditions à la limite à réaction locale, pour lesquelles \(\alpha\) , \(\beta\) et \(\gamma\) sont des scalaires; les cas spécifiables sont les suivants:

  • Le cas \(\alpha =0\) , \(\beta \ne 0\) et \(\gamma \ne 0\) désigne une frontière du domaine à vitesse vibratoire imposée. En effet, il existe une relation reliant le gradient de pression acoustique complexe à la vitesse vibratoire normale particulaire complexe. C’est-à-dire:

(903)#\[\frac{\partial p}{\partial n}=-j\omega {\rho}_{0}{V}_{n}\]

\({\rho}_{0}\) désigne la masse volumique du fluide considéré, et on impose \({V}_{n}\) , la vitesse vibratoire normale à la paroi (\({V}_{n}=v\cdot n\)\(n\) désigne le vecteur unitaire de la normale extérieure à la frontière \(\partial \Omega\) ). Une paroi étanche vérifie ainsi la condition naturelle: \(\frac{\partial p}{\partial n}=0\) .

  • Le cas \(\alpha \ne 0\) , \(\beta \ne 0\) et \(\gamma =0\) concerne une frontière à impédance acoustique \(Z\) imposée. L’impédance acoustique \(Z\) est définie comme le rapport de la pression à la vitesse vibratoire particulaire normale au voisinage de la paroi à impédance imposée:

  • Le cas \(\alpha \ne 0\) , \(\beta =0\) et \(\gamma \ne 0\) représente le cas où l’on impose la pression acoustique \(p\) à une frontière (le plus souvent \(\gamma =0\) , correspondant à \(p=0\) ).

(904)#\[Z=\frac{p}{{V}_{n}}\]
  • Le cas \(\alpha \ne 0\) , \(\beta =0\) et \(\gamma \ne 0\) représente le cas où l’on impose la pression acoustique \(p\) à une frontière (le plus souvent \(\gamma =0\) , correspondant à \(p=0\) ).

Formulation en pression#

Expression mathématique du problème#

On utilise la procédure standard visant à poser le problème aux éléments finis isoparamétriques classiques. On suppose la solution du problème suffisamment régulière, \(p\in {H}^{2}\left(\Omega \right)\) . On multiplie l’équation de Helmholtz () sans terme source \(s\) par une fonction test \(\varphi ` . On intègre sur :math:\)Omega ` et on utilise la formule de Green. D’après (), la frontière \(\partial \Omega\) du domaine \(\Omega\) , se subdivise en deux zones, une zone à vitesse vibratoire imposée, \(\partial {\Omega}_{v}\) et une zone à impédance acoustique imposée, \(\partial {\Omega}_{Z}\) . L’équation obtenue peut être réécrite sous la forme:

(905)#\[{\int}_{\Omega}\mathrm{grad}(p)\cdot \mathrm{grad}(\varphi )\mathit{dV}-{\int}_{\Omega}\frac{{\omega}^{2}}{{c}^{2}}p\varphi \mathit{dV}+j{\int}_{\partial {\Omega}_{Z}}\frac{{\rho}_{0}\omega }{Z}p\varphi \mathit{dS}+j{\int}_{\partial {\Omega}_{v}}{\rho}_{0}\omega {V}_{n}\varphi \mathit{dS}=0\]

\(\mathit{dV}\) représente un élément de volume différentiel dans \(\Omega\) et \(\mathit{dS}\) représente un élément de surface sur \(\partial \Omega\) . La vitesse vibratoire particulaire est ensuite déterminée par:

(906)#\[\mathrm{v}=\frac{j}{{\rho}_{0}\omega }\mathrm{grad}(p)\]

Discrétisation par éléments finis#

Dans le cas des éléments finis isoparamétriques classiques, les intégrales élémentaires sont \({\mathrm{K}}^{e}\) , \({\mathrm{M}}^{e}\) , \({\mathrm{C}}^{e}\) et \({\mathrm{S}}^{e}\) suivant la décomposition indiquée par () (\({\mathrm{K}}^{e}\) est la matrice de «rigidité», \({\mathrm{M}}^{e}\) la matrice de «masse», \({C}^{e}\) la matrice d’amortissement et \({\mathrm{S}}^{e}\) le vecteur source). Deux d’entre elles proviennent d’intégrales volumiques, les deux autres sont le résultat d’intégrales respectivement sur une surface vibrante et sur une surface à impédance imposée. On supposera que les coordonnées globales d’un élément peuvent s’écrire grâce à la donnée de \(m\) fonctions de forme élémentaires \({H}_{i}\) :

(907)#\[{\mathit{OM}}^{e}={\sum}_{i=1}^{m}{N}_{i}{\mathit{OM}}_{i}^{e}\]

On se donne en outre, des fonctions de base \({N}_{i}\) , pour décrire la pression acoustique élémentaire. La pression à l’intérieur d’un élément pourra s’écrire:

(908)#\[{p}^{e}(x,y,z)={\sum}_{i=1}^{m}{N}_{i}{p}_{i}^{e}\]

\({p}_{i}^{e}\) est la pression acoustique au nœud \(i\) de l’élément \(e\) . Dans le cas des éléments finis isoparamétriques, les fonctions de base \({N}_{i}\) sont égales aux fonctions de forme \({H}_{i}\) . Sur chaque élément du domaine, le problème aux éléments finis en pression acoustique s’écrit:

(909)#\[\left({\mathrm{K}}^{e}-{\omega}^{2}{\mathrm{M}}^{e}+j\omega {\mathrm{C}}^{e}\right){p}^{e}=-j\omega {\mathrm{S}}^{e}\]

\({p}^{e}\) est le vecteur-colonne des valeurs nodales de la pression acoustique sur l’élément.

La matrice de « rigidité »#

La matrice de «rigidité» \({\mathrm{K}}^{e}\) correspond au calcul du premier terme de (), soit:

(910)#\[{\int}_{{\Omega}^{e}}\left(\mathrm{grad}(p)\cdot \mathrm{grad}(\varphi )\right)\mathit{dV}\]

Elle correspond à l’énergie cinétique du fluide acoustique. Elle admet comme terme général:

(911)#\[{K}_{ij}^{e}={\int}_{{\Omega}^{e}}\left(\nabla {N}_{i}\nabla {N}_{j}\right).\mathit{dV}\]

Remarque:

à cause du choix de la pression comme variable principale et de la dérivation seconde en temps, le premier terme de () correspond à l’énergie cinétique *.* Pour une raison pratique on le nomme malgré tout matrice de «rigidité».

La matrice de « masse »#

La matrice de «masse» \({\mathrm{M}}^{e}\) correspond au calcul du deuxième terme de (), soit:

(912)#\[{M}_{ij}^{e}={\int}_{{\Omega}^{e}}\frac{1}{{c}^{2}}p\varphi \mathit{dV}\]

Elle correspond à l’énergie de compressibilité du fluide acoustique. Elle admet comme terme général:

(913)#\[{M}_{ij}^{e}={\int}_{{\Omega}^{e}}\frac{1}{{c}^{2}}{N}_{i}{N}_{j}.\mathit{dV}\]

Remarque1:

Acause du choix de la pression comme variable principale et de la dérivation seconde en temps, le deuxièmeterme de () correspond à l’énergie élastique. Pour une raison pratique on le nomme malgré tout matrice de «masse».

Remarque2:

Dans le cas d’une propagation en milieu viscoacoustique hystérétique, la célérité du son \(c\) qui intervient dans () est un nombre complexe.

La matrice d’amortissement#

La matrice d’amortissement \({\mathrm{C}}^{e}\) correspond au calcul du troisième terme de (), soit:

(914)#\[{C}_{ij}^{e}={\int}_{\partial {\Omega}_{Z}^{e}}\frac{{\rho}_{0}}{Z}p\varphi \mathit{dS}\]

Elle admet comme terme général, sur les frontière concernées:

(915)#\[{C}_{ij}^{e}={\int}_{\partial {\Omega}_{Z}^{e}}\frac{{\rho}_{0}}{Z}{N}_{i}{N}_{j}\mathit{dS}\]

Le vecteur source#

Le vecteur source \({\mathrm{S}}^{e}\) correspond au calcul du dernier terme de (), soit:

(916)#\[{S}_{i}^{e}={\int}_{\partial {\Omega}_{v}}{\rho}_{0}{V}_{n}\varphi \mathit{dS}=0\]

Il admet comme terme général:

(917)#\[ \begin{align}\begin{aligned}{S}_{i}^{e}={\int}_{\partial {\Omega}_{v}^{e}}{\rho}_{0}{V}_{n}{N}_{i}\mathit{dS}\\{S}_{i}^{e}={\int}_{\partial {\Omega}_{v}^{e}}{\rho}_{0}{V}_{n}{N}_{i}\mathit{dS}\end{aligned}\end{align} \]

Commandes spécifiques à la modélisation acoustique#

Lors d’une étude par modélisation en éléments finis acoustiques avec Code_Aster on utilise des commandes générales et des commandes qui sont propres à l’acoustique, ou dont les mots-clés et options sont particulières à cette discipline; nous en présentons ci-dessous la liste.

Les valeurs nodales de la pression acoustique sur l’élément fini (degrés de liberté) sont notées: PRES.

Définition des caractéristiques des milieux de propagation#

Il est nécessaire de donner la masse volumique (valeur réelle) et la célérité de propagation (valeur complexe) ; on utilise pour cela la commande DEFI_MATERIAU `[U4.43.01] <https://www.code-aster.org/V2/doc/default/fr/man_u/u4/u4.43.01.pdf>`_avec les mots-clés suivants:

mot-clé facteur:

FLUIDE

mots-clés:

RHO

(masse volumique \({\rho}_{0}\) )

CELE_C

(célérité complexe \({c}_{0}\) )

Exemple:

air = DEFI_MATERIAU(FLUIDE=_F(RHO= 1.3, CELE_C = (“RI”,343.0,0.,))

Dans ce cas: \({c}_{0}=343.0\) est réel (la partie imaginaire – qui représenterait une dissipation d’énergie matérielle – est nulle).

Conditions aux limites#

On doit affecter des valeurs de vitesse vibratoire normale par face (ou arête en bidimensionnel) aux mailles définissant les frontières sources, et aussi des valeurs d’impédance acoustique par face (arête en bidimensionnel) aux mailles définissant les frontières à impédance imposée. On utilise la commande spécifique à l’acoustique AFFE_CHAR_ACOU [U4.44.04]

Calcul des matrices élémentaires#

Les différentes matrices élémentaires (rigidité, masse et amortissement) sont calculées par des options spécifiques. On emploie la commandeCALC_MATR_ELEM [U4.61.01] avec le mot-clé OPTION pour lequel on spécifie les valeurs d’affectation possibles :

Mots-clés:

OPTION

“RIGI_ACOU”

Équation ()

“MASS_ACOU”

Équation ()

“AMOR_ACOU”

Équation ()

Ces matrices sont de type complexe.

Remarque:

Les matrices assemblées peuvent être obtenuesdirectement avec la macro-commande ASSEMBLAGEet les mêmes options.

Calcul du vecteur source élémentaire#

Le vecteur élémentaire est calculé par une option spécifique; il faut obligatoirement indiquer le chargement. On emploie la commande CALC_VECT_ELEM [U4.61.02] avec le mot-clé OPTION

Calcul de la solution#

Après assemblage des matrices et vecteur élémentaires la solution harmonique peut se calculer directement avec la commande DYNA_VIBRA [U4.53.03].

On peut également calculer avec les matrices assemblées les modes acoustiques parfaits (sans viscosité, donc avec une célérité réelle) sur le domaine borné modélisé. Il faut donc utiliser les parties réelles des matrices assemblées en entrée de la commande de calcul modal [U4.52.02]. Le concept du résultat est alors du type ou “MODE_ACOU’.

Post-traitements#

À partir du résultat de la résolution, des commandes de post-traitement permettent d’obtenir les champs nodaux de grandeurs acoustiques suivantes:

  • Niveau \({L}_{p}\) de pression acoustique \(p\) en \(\mathit{dB}\) : \({L}_{p}=20{\log}_{10}\left[\frac{\mid p\mid }{2,0\times {10}^{-5}}\right]\)

  • partie réelle de la pression acoustique

  • partie imaginaire de la pression acoustique

  • intensité acoustique active \(I=\frac{1}{2}\text{Re}\left[p{v}^{\text{*}}\right]\)

  • intensité acoustique réactive \(J=\frac{1}{2}\text{Im}\left[p{v}^{\text{*}}\right]\)

Ces champs sont calculés par utilisation de la commande de post-traitement CALC_CHAMP (le concept du résultat est du type “ACOU_HARMO” ou “MODE_ACOU”) avec les mots-clés RESULTAT et OPTION.

« Mot-clé : », « RESULTAT », « « , «  » « Mot-clé : », « ACOUSTIQUE », « “PRAC_ELNO” », « (niveau de pression en décibels) » « « , « « , « « , « (partie réelle de la pression) » « « , « « , « « , « (partie imaginaire de la pression) » « « , « « , « “INTE_ELNO” », « (intensité active) » « « , « « , « « , « (intensité réactive) »

Conclusion#

Une modélisation a donc été ajoutée, permettant de faire des calculs d’acoustique intérieure en basses fréquences pour des géométries complexes par des éléments finis acoustiques classiques. La formulation a été validée par comparaison à une solution analytique; des cas tests sont présentés dans le manuel de validation V7 sous la codification AHLV100.

Bibliographie#

    1. BOUIZI: “Mise en œuvre d’un code de calcul d’éléments finis en vue du traitement de l’équation de Helmholtz en espace clos” - Travail de fin d’études, E.C.L. 1986.

    1. BOUIZI: “Analyse spectrale de l’équation de Helmholtz.” - Rapport de DEA, Ecole Centrale de Lyon, 1986.

    1. BOUIZI: “Éléments finis mixtes en acoustique linéaire stationnaire : Développement du code AIRMEF” - Département Acoustique. DER - EDF. HE‑24/88.02. 1988.

    1. BOUIZI, M. COURTADE, D. JEANDEL, E. LUZZATO, A. MIGNOT, C. SURRY.: “Conditions de compatibilité de Brezzi-Babuska pour des méthodes d’éléments finis mixtes conformes en Mécanique et Acoustique” AUM, Actes du 8ème Congrès Français de Mécanique, Nantes, 1987.

    1. BOUIZI, M. COURTADE, D. JEANDEL, E. LUZZATO, C. SURRY: “Traitement de l’équation de Helmholtz par un code d’éléments finis mixtes en espace clos” GAMI, Colloque Vibrations Chocs, 1988, ECL, 1988.

    1. BOUIZI : “Résolution des équations de l’Acoustique linéaire par une méthode d’éléments finis mixtes”. Thèse présentée devant l’École Centrale de Lyon -Spécialité : Mécanique -. Soutenue le 02/03/89.

    1. HABASQUE: “Validation expérimentale du code de calcul d’acoustique interne, Basse Fréquence”. Rapport de stage de DEA. ECL 1986 (+ Projet de fin d’études).

  1. A.ADOBES: “Etude numérique et expérimentale des champs d’ondes stationnaires” Rapport DER / EDF - HE-2287.22

    1. STIFKENS, A.ADOBES: “Bilan de l’intégration des éléments finis classiques dans Aster” - Rapport DER / EDF - HP-64/91.149

    1. STIFKENS: “Intégration des éléments finis acoustiques mixtes dans Aster” Rapport DER / EDF – HP-61/92.081.