# Une approche ondulatoire adaptative pour l'acoustique dans le domaine des moyennes fréquences

Benjamin Sourcis, Pierre Ladevèze, Hervé Riou, Béatrice Faverjon

► To cite this version:

Benjamin Sourcis, Pierre Ladevèze, Hervé Riou, Béatrice Faverjon. Une approche ondulatoire adaptative pour l'acoustique dans le domaine des moyennes fréquences. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01413786

HAL Id: hal-01413786

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

Submitted on 11 Dec 2016

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 approche ondulatoire adaptative pour l'acoustique dans le domaine des moyennes fréquences

B. Sourcis¹, P. Ladevèze¹,², H. Riou¹, B. Faverjon¹

¹ LMT-Cachan (ENS de Cachan/CNRS UMR8535/Université Paris 6/Universud Paris)
61, avenue du Président Wilson, F-94235 CACHAN CEDEX, France
{sourcis, ladeveze, riou, faverjon} @lmt.ens-cachan.fr

² Chaire de la Fondation EADS 'Techniques Avancées en Calcul des Structures'

**Résumé** — Une approche multiéchelle appelée « Théorie Variationnelle des rayons Complexes » (TVRC) a été récemment proposée pour le calcul de structures faiblement amorties dans le domaine des moyennes fréquences. Cette méthode est une approche ondulatoire, car elle utilise un ensemble d'ondes planes pour approximer la solution. L'objectif de ce travail est d'étendre la méthode au problèmes d'acoustique en 3-D. Mais en comparaison avec les problèmes 2-D, cette dimension supplémentaire engendre un coût de calcul élevé. Pour pallier ce problème nous proposons une procédure adaptative de calcul basée sur un indicateur d'erreur. Le problème peut être ainsi discrétisé efficacement en terme de directions de propagation des ondes, car seules les ondes ayant une participation énergétique significative dans la solution sont prises en compte.

**Mots clés** — acoustique, procédure adaptative, TVRC.

## 1 Introduction

Les méthodes des éléments finis (FEM) et des éléments de frontières (BEM) sont classiquement utilisées pour résoudre des problèmes d'acoustique. Toutefois, l'utilisation de fonctions polynômiales conduit à la construction de gros modèles numériques et, en pratique, limite l'application de ces techniques prédictives au domaine des basses fréquences. En outre, l'effet de pollution [1] nuit à la robustesse de la FEM quand le nombre d'onde augmente. En effet, pour les grands nombres d'onde, le degré de raffinement du maillage nécessaire pour maintenir l'erreur commise par le cacul raisonnablement petite requiert un grand nombre de noeuds, ce qui rend prohibitif le coût des calculs dans le domaine des moyennes fréquences.

Pour pallier les problèmes liés à la méthode classique des éléments finis, de nombreuses techniques ont été proposées, comme par exemple les méthodes de bases réduites [12], ou encore la Galerkin FEM [5]. Par rapport à la méthode classique, ces techniques permettent de réduire le coût des calculs de façon drastique. Toutefois, la gamme de fréquences traitées par ces méthodes reste limitée aux basses fréquences.

D'autres approches approximent la solution du problème de vibrations (ou d'acoustique) par des fonctions oscillantes qui sont solutions exactes de l'équation d'équilibre. Ces méthodes incluent, par exemple, une utilisation particulière de la méthode de partition de l'unité (PUM) [13], la méthode variationnelle ultra-faible (UWVF) [2], la Discontinuous Galerkin Method (DGM) [4], la Wave Boundary Element Method (WBEM) [9], ou encore la Wave Based Method (WBM) [3]. La Théorie Variationnelle des Rayons Complexes (TVRC), succintement présentée dans ce papier,
appartient à cette catégorie de techniques numériques, que nous qualifieront d'approches « ondu-
latoires ». La TVRC a été initialement proposée dans [6] afin de trouver une solution approchée
aux problèmes de vibrations dans le domaine moyennes fréquences, qui s'avèrent trop onéreux si
l'on utilise la FEM ou la BEM. Elle est basée sur formulation variationnelle du problème qui a été
développée dans le but de permettre aux approximations faites dans chaque élément d'être a priori
indépendantes les unes des autres. Ainsi, dans chaque sous-domaines, tout type de fonctions de
forme peut être utilisé, ce qui donne une grande flexibilité et robustesse à la méthode. Ces dernières
années, la méthode a été développée avec succès pour les vibrations structurales [7, 10, 8]. Une
première utilisation de la TVRC pour résoudre des problèmes de Helmholtz en 2-D a récemment
été proposée dans [11]. Les travaux en cours sont axés vers une extension de la méthode au 3-D
qui engendre coûts de calcul élevés. Par conséquent, nous proposons une stratégie adaptative de
calcul basée sur un indicateur d'erreur permettant de discrétiser de façon adpatée le portrait d'am-
plitudes. Le principal inconvénient d'une telle représentation est la quasi-non-unicité de la solution
(dans un sens numérique). En effet, différentes distributions d'amplitudes conduisent à des solu-
tions numériquement équivalentes. Pour pallier ce problème, nous proposons une technique de
régularisation.

2 Présentation de la TVRC

2.1 Problème de référence en acoustique linéaire

<figure><img src="image_2.png" /><figcaption>Figure 1 – Problème de référence.</figcaption></figure>

Par souci de simplicité, nous présentons ici le problème d'une cavité acoustique fermée Ω séparée en 2 sous-cavités Ω₁ et Ω₂ par une interface notée Γ (voir figure 1). Sur le bord des ces sous-cavités sont appliquées trois types de conditions aux limites (Dirichlet, Neuman et Robin). Le problème de référence consiste à trouver les champs de pression acoustique (p₁, p₂) ∈ H¹(Ω₁)× H¹(Ω₂) tels que :

$$
\Delta p_e + k^2 p_e = 0 \quad \text{dans } \Omega_e, \quad e \in \{1, 2\} \tag{1}
$$

$$
\left|
\begin{array}{ll}
p_e = p_{de} & \text{sur } \partial_p \Omega_e \\
\mathcal{L}_v(p_e) = v_{de} & \text{sur } \partial_v \Omega_e, \quad e \in \{1, 2\} \\
p_e - Z_e \mathcal{L}_v(p_e) = h_{de} & \text{sur } \partial_z \Omega_e
\end{array}
\right|
$$

$$
\begin{gather*}
p_1 - p_2 = 0 \quad \text{sur } \Gamma \\
\mathcal{L}_v(p_1) + \mathcal{L}_v(p_2) = 0 \quad \text{sur } \Gamma
\end{gather*}
\tag{2}
$$

où k = ω/c₀ désigne le nombre d'onde (avec c₀ la vitesse du son et ω la pulsation), L_v est l'opé-
rateur vitesse défini par L_v(□_e) = i/(ρ₀ω)∇(□_e)·n_e, (avec n_e la normale sortante au domaine Ω_e,
ρ₀ la masse volumique du fluide et i = √-1 l'unité imaginaire pure), p_de la pression imposée
sur ∂_pΩ_e, v_de la vitesse imposée sur ∂_vΩ_e, h_de l'excitation imposée sur ∂_ZΩ_e et Z_e une fonction
d'impédance complexe dont la partie réelle introduit de l'amortissement dans le modèle. (1) est
l'équation de Helmholtz homogène.
## 2.2 Reformulation et discrétisation du problème

La TVRC repose sur une formulation faible du problème de référence qui permet aux approximations faites sur chaque sous-cavité d'être a priori indépendantes. Ceci donne à la TVCR une grande souplesse puisque n'importe quel type de fonctions de forme peut être utilisé. Soit $S_{0,ad}^{(e)}$ l'espace des fonctions $p_e$ satisfaisant (1) $(S_{0,ad}^{(e)} = \{p_e \in H^1(\Omega_e) \mid \Delta p_e + k^2 p_e = 0 \text{ dans } \Omega_e\})$. Le problème aux dérivées partielles (1)-(2) est équivalent au problème variationnel suivant : trouver $(p_1, p_2) \in S_{ad}^{(1)} \times S_{ad}^{(2)}$ tel que

$$ \mathcal{R} \left\{ \sum_{e=1,2} \int_{\partial_p \Omega_e} (p_e - p_{de}) \cdot \overline{\mathcal{L}_v[\delta p_e]} ds + \sum_{e=1,2} \int_{\partial_v \Omega_e} (\overline{\mathcal{L}_v[p_e] - v_{de}}) \cdot \delta p_e ds \right. \\ + \sum_{e=1,2} \int_{\partial_z \Omega_e} \left( \left( (1-Z_e \mathcal{L}_v)[p_e] - h_{de} \right) \cdot \overline{\mathcal{L}_v(\delta p_e)} + \left( \overline{(\mathcal{L}_v - 1/Z_e)[p_e] + h_{de}/Z_e} \right) \cdot \delta p_e \right) ds \\ \left. + \frac{1}{2} \int_{\Gamma} \left( (p_1-p_2) \cdot \overline{\mathcal{L}_v[\delta p_1 - \delta p_2]} + \overline{\mathcal{L}_v[p_1+p_2]} \cdot (\delta p_1 + \delta p_2) \right) ds \right\} = 0, \\ \forall (\delta p_1, \delta p_2) \in S_{0,ad}^{(1)} \times S_{0,ad}^{(2)} \quad (3) $$

où $\overline{\square}$ et $\mathcal{R}\{\square\}$ désignent respectivement le conjugué et la partie réelle de la quantité complexe $\square$. La preuve de cette équivalence est une extension directe de la démonstration figurant dans [6]. Notons que l'existence et l'unicité de la solution sont assurées si et seulement si $\mathcal{R}(Z_e) > 0$, car c'est cette quantité qui introduit de l'amortissement dans le modèle. Notons également que la continuité faible des champs au niveau de l'interface $\Gamma$ est directement prise en compte dans la formulation (3). Aucune technique de pénalisation ou de multiplicateurs de Lagrange ne sont requises pour assurer cette continuité, ce qui montre la souplesse et la robustesse de la TVRC.

La TVRC cherche une solution du problème sous la forme :

$$ p_e(\mathbf{x}) = \int_{C_e} A_e(\mathbf{k}_e) e^{i\mathbf{k}_e \cdot \mathbf{x}} dC_e \quad (4) $$

où $\mathbf{k}_e$ est le vecteur d'onde tel que $e^{i\mathbf{k}_e \cdot \mathbf{x}}$ vérifie exactement (1), et $C_e$ désigne la courbe décrite par le vecteur d'onde $\mathbf{k}_e$ quand ce dernier parcoure toutes les directions du plan (en 2-D). En prenant $\mathbf{k}_e = k[\cos\theta, \sin\theta], \theta \in [0, 2\pi[C_e$ est alors un cercle), on a l'ensemble des rayons propagatifs parcourant toutes les directions du plan (voir figure 2(a)). L'expression (4) présente deux échelles : une échelle lente constituée par la distribution des amplitudes $A_e(\mathbf{k}_e)$ qui varie a priori lentement quand le vecteur d'onde $\mathbf{k}_e$ varie, ce qui en d'autres termes signifie que cette distribution oscille lentement, et une échelle rapide constituée par le terme $e^{i\mathbf{k}_e \cdot \mathbf{x}}$ caractérisant le phénomène vibratoire et qui varie rapidement avec $\mathbf{k}_e$ et la variable spatiale $\mathbf{x}$. Cette échelle est prise en compte de façon analytique. L'approximation (4) est discrétisée de la façon suivante :

$$ p_e^h(\mathbf{x}) = \int_{C_e} A_e^h(\mathbf{k}_e) e^{i\mathbf{k}_e \cdot \mathbf{x}} dC_e = \sum_{l=1}^{L_e} a_l^{(e)} \int_{C_e^l} e^{i\mathbf{k}_e \cdot \mathbf{x}} dC_e \quad (5) $$

où $A_e^h$ représente l'approximation de $A_e$, $C_e^l$ une discrétisation angulaire du cercle $C_e$ et $a_l^{(e)}$ la valeur de $A_e^h$ supposée constante sur le secteur angulaire $C_e^l$. Une telle discrétisation est représentée sur la figure 2(b) : le cercle $C_e$ est approximé par huit secteurs angulaires $C_e^l$, et dans chacun de ces secteurs la distribution d'amplitude est approximée par la valeur constante $a_l^{(e)}$. En injectant (5) dans la formulation variationnelle (3), on aboutit à la construction d'un système linéaire de $(L_1 + L_2)$ équations algébriques dont la solution est le vecteur des $(L_1 + L_2)$ amplitudes généralisées $\mathbf{a} = [\mathbf{a}_1^T \mathbf{a}_2^T]^T$ :
<figure><img src="image_3.png" /><figcaption>Figure 2 – Exemple de fonction de forme et discrétisation du portrait d'amplitudes.</figcaption></figure>

$$ \begin{bmatrix} \mathbf{K}_{11} & \mathbf{K}_{12} \\ \mathbf{K}_{21} & \mathbf{K}_{22} \end{bmatrix} \cdot \begin{bmatrix} \mathbf{a}_1 \\ \mathbf{a}_2 \end{bmatrix} = \begin{bmatrix} f_1 \\ f_2 \end{bmatrix} \qquad (6) $$

En résolvant (6) et en injectant le résultat dans (5), on obtient une solution approchée du problème d'acoustique.

## 2.3 Régularisation des amplitudes

L'approximation (5) possède un propriété paradoxale. En effet, la solution $p_e^h(\mathbf{x})$ peut être représentée par des portraits d'amplitudes $A_e^h(\theta)$ très différents, ce qui représente un inconvénient majeur de la méthode. En conséquence de cette propriété, la matrice TVRC est très mal conditionnée. C'est aujourd'hui un aspect bien connu des approches ondulatoires présentées en introduction. Bien que la convergence de la méthode en fonction du nombre de degrés de liberté (ddls) soit exponentielle [11], elle s'accompagne également d'une rapide inflation du nombre de conditionnement. Pour obtenir une solution significative en terme de distribution d'amplitudes, une technique de régularisation est utilisée. Au lieu de résoudre (6), nous cherchons la solution régularisée $\mathbf{a}_{\lambda_p}$ du problème suivant

$$ \begin{bmatrix} \mathbf{K}_{11} + \lambda_p \mathbf{H}_1 & \mathbf{K}_{12} \\ \mathbf{K}_{21} & \mathbf{K}_{22} + \lambda_p \mathbf{H}_2 \end{bmatrix} \cdot \begin{bmatrix} \mathbf{a}_{1,\lambda_p} \\ \mathbf{a}_{2,\lambda_p} \end{bmatrix} = \begin{bmatrix} f_1 \\ f_2 \end{bmatrix} \qquad (7) $$

où l'opérateur $\mathbf{H}_e$ peut être vu comme la norme de la dérivée seconde des amplitudes généralisées $\mathbf{a}_e$ par rapport à la variable angulaire $\theta$, et $\lambda_p$ est un paramètre régularisation qui contrôle le poids du terme de régularisation. Soit $\mathbf{b}$ un vecteur dans $\mathbb{C}^{L_e}$, $\mathbf{H}_e$ est construit comme suit

$$ \delta \mathbf{b}^* \cdot \mathbf{H}_e \cdot \mathbf{b} = C(\text{mes}(\Omega_e)k) \left( \delta \left| \frac{2b_1 - b_2 - b_n}{\Delta\theta^2} \right|^2 + \sum_{i=2}^{L_e-1} \delta \left| \frac{2b_i - b_{i+1} - b_{i-1}}{\Delta\theta^2} \right|^2 + \delta \left| \frac{2b_{L_e} - b_1 - b_{L_e-1}}{\Delta\theta^2} \right|^2 \right) \qquad (8) $$

où $\Delta\theta$ désigne l'écart angulaire de chaque faisceau. La matrice représentant $\mathbf{H}_e$ est creuse et symétrique et son coût de calcul est négligeable. Pour illustrer l'efficacité de cette technique, la figure 3 montre des solutions TVRC obtenues avec et sans régularisation. Le choix du paramètre $\lambda_p$ ne sera pas discuté ici.

## 2.4 Vers une stratégie adaptative de calcul

Une fois le comportement des amplitudes maitrisés, il devient possible de mettre en place une démarche adaptative de calcul. Le choix retenu ici s'appuie sur un indicateur d'erreur défini à
<figure><img src="image_4.png" /><figcaption>Figure 3 – Comparaison entre des solutions du système (6) (a), et des solutions régularisées du système (7) (b).</figcaption></figure>

l'échelle de la sous-cavité. La stratégie se déroule en trois étapes notées E :

E1. une première analyse globale utilisant une discrétisation grossière des portraits d'amplitudes est effectuée.

E2. des calculs utilisant une discrétisation très fine sont effectués à l'échelle de la sous-cavité, en applicant sur les bords de cette dernière la trace des champs calculés à l'étape précédente. A ce niveau est défini un indicateur d'erreur à la solution. Les énergies potentielles ($e_p = 1/(4\rho_0c_0) \int_{\Omega_e} p_e^h p_e^h d\Omega_e$) associées à trois solutions sont calculées en applicant sur le bord de $\Omega_e$ la pression, la vitesse normale et les deux (la TVRC permet en effet de traiter les problèmes dits « mal posés » en ajoutant un terme à (3)), et l'indicateur est défini par :

$$I_{\Omega_e} = \frac{2e_p^{p,v} - (e_p^v + e_p^h)}{2e_p^{p,v} + (e_p^v + e_p^h)} \quad (9)$$

Cet indicateur permet de vérifier la compatibilité des champs de pression et de vitesse calculés à l'étape E1. En fonction de sa valeur, on revient éventuellement à l'étape précédente, ou alors on calcule une discrétisation adaptée qui est un intermédiaire entre deux grilles pré-définies.

E3. une analyse globale utilisant les grilles calculées à l'étape précédente est effectuée.

Cette procédure peut être éventuellement rendue itérative pour atteindre le degré de précision voulu. La figure 4 montre le taux de convergence de la procédure adaptative comparé à celui utilisant une discrétisation régulière des portraits d'amplitudes.

## 3 Conclusion

Ce papier présente le premier développement d'une stratégie adaptative de calcul dans le cadre des approches ondulatoires. Cette stratégie sera plus amplement présentée et illsutrée avec des exemples 2-D.
<figure><img src="image_5.png" /></figure>

Figure 4 – Comparaison des taux convergences entre une discrétisation régulière et en utilisant la procédure adaptative de calcul.

## 3.1 Références bibliographiques

### Références

[1] I. Babuška, S.A. Sauter. Is the pollution effect of the FEM avoidable for the Helmholtz equation considering high wave numbers? *SIAM Review*, 42(3), 451–484, 1997.

[2] O. Cessenat, B. Despres. Application of an Ultra Weak Variational Formulation of Elliptic PDEs to the Two-Dimensional Helmholtz Problem. *SIAM Journal on Numerical Analysis*, 35(1), 255–299, 1998.

[3] W. Desmet, P. Sas, D. Vandepitte. An indirect Trefftz method for the steady-state dynamic analysis of coupled vibro-acoustic systems. *Computer Assisted Mechanics and Engineering Sciences*, 8, 271–288, 2001.

[4] C. Farhat, I. Harari, U. Hetmaniuk. A discontinuous Galerkin method with plane waves and Lagrange multipliers for the solution of Helmholtz problems in the mid-frequency regime. *Computer Methods in Applied Mechanics and Engineering*, 192, 1389–1419, 2003.

[5] I. Harari, T. J. R. Hughes. Galerkin/least-squares finite element methods for the reduced wave equation with non-reflecting boundary conditions in unbounded domains. *Computer Methods in Applied Mechanics and Engineering*, 98(3), 411–454, 1992.

[6] P. Ladevèze. A new computational approach for structure vibrations in the medium frequency range. *Compte Rendus de l'Académie des Sciences de Paris*, 322(II b), 849–856, 1996.

[7] P. Ladevèze, P. Rouch, H. Riou, X. Boniheust. Analysis of Medium-Frequency Vibrations in a Frequency Range. *Journal of Computational Acoustics*, 11(2), 255–284, 2003.

[8] P. Ladevèze, M. Chevreuil. A new computational method for transient dynamics including the low- and the medium-frequency ranges. *International Journal for Numerical Methods in Engineering*, 64(4), 503–527, 2005.

[9] E. Perrey-Debain, J. Trevelyan, P. Bettess. Wave boundary elements : a theoretical overview presenting applications in scattering of short waves. *Engineering Analysis with Boundary Elements*, 28(2), 131–141, 2004.

[10] H. Riou, P. Ladevèze, P. Rouch. Extension of the Variational Theory of Complex Rays to shells for medium-frequency vibrations. *Computer Methods in Applied Mechanics and Engineering*, 272(1-2), 341–360, 2004.

[11] H. Riou, P. Ladevèze, B. Sourcis. The multiscale VTCR approach applied to acoustics problems. *Journal of Computational Acoustics*, 16(4), 487–505, 2008.

[12] C. Soize. Reduced models in the medium frequency range for the general dissipative structural dynamic systems. *European Journal of Mechanics, A/Solids*, 17, 657–685, 1998.

[13] T. Strouboulis, I. Babuška, R. Hidajat. The generalized finite element method for Helmholtz equation : Theory, computation, and open problems. *Computer Methods in Applied Mechanics and Engineering*, 195, 4711–4731, 2006.