Computational approach for the acoustic modelling of large aftertreatment devices with multimodal incident sound fields

Abstract

The influence of multimodal incident sound fields on the acoustic behaviour of large aftertreatment devices incorporating a monolith is modelled and analysed in detail. The analytical mode matching method is applied to the compatibility conditions of the three-dimensional acoustic fields at the device geometric discontinuities, leading to the computation of the complex wave amplitudes in all the subdomains involved and the corresponding device transmission loss. To have a realistic model, three-dimensional propagation must be considered in the inlet/outlet ducts and chambers, while one-dimensional wave propagation has to be assumed along the small capillaries of the aftertreatment device monolith (such catalytic converters and particulate filters); therefore, the monolith can be replaced by a plane wave four-pole transfer matrix from an acoustical point of view. On the other hand, for large aftertreatment device inlet ducts such as those found in heavy-duty and off-road engines, the usual models with plane incident wave excitation are not accurate since the onset of higher order incident modes in the inlet duct is expected for the frequency range of interest. Therefore, a variation of the acoustic attenuation performance is likely to occur depending on these modes, similar to the results previously found in the case of large dissipative silencers. Results are presented for three different multimodal incident sound field hypotheses: equal modal amplitude, equal modal power and equal modal energy density. A relevant influence on the sound attenuation is found for the test problems considered in the current investigation.

Introduction

Large silencers are commonly used to attenuate noise in many ventilation applications, gas turbine systems, as well as large internal combustion engines (trucks, freight trains, power generation, etc.), among others. For their acoustic modelling and design, the well-established computational methods considered for the relatively small automotive exhaust silencers^1–9 can be extended. As these models for small devices usually assume that plane waves are present in the inlet/outlet ducts even if higher order acoustic modes are propagating inside the main body of the silencer, some problems may arise when modelling and computing the acoustic attenuation performance of large devices. If numerical methods such as the finite element method (FEM) are applied,^5–8 a large increase in the number of degrees of freedom is required to achieve sufficiently accurate predictions at higher frequencies. When analytical approaches^2–4,10 are taken into account, the increase in the number of propagating higher order modes for large configurations significantly affects the silencer performance. ¹¹ This influence is associated not only with the modes inside the central chamber, but also with the modes present in the incident sound field within the inlet port.^11–16 It is worth noting that, for large inlet ducts, the excitation is not restricted to a plane wave only and this incident propagation hypothesis no longer holds; thus, a model for multimodal excitation is necessary.

The aforementioned considerations can be extended to a number of aftertreatment devices (ATD) commonly used for the control of pollutant emissions that also have a relevant contribution to noise attenuation, such catalytic converters and particulate filters. One of the main features of these ATD is that they incorporate a central monolith reactor, usually made of ceramic or metallic material with a large number of small parallel channels with thin walls running axially through the core. ¹⁷ While the ATD acoustic attenuation performance can be assessed accurately for sizes typical of small automotive exhaust systems assuming incident plane wave excitation,^18–23 further research is required for large ATD since the ducts and chambers involved are big enough to allow higher order modes to propagate not only inside the main body of the device but also along the inlet/outlet ports. Therefore, it cannot be assumed that the ATD incident sound pressure field consists solely of a plane acoustic wave.

Following this introduction, Section ‘Overview of the acoustic problem’ presents the acoustic problem under consideration. To have a realistic model, three-dimensional (3D) propagation must be considered in the inlet/outlet ducts and chambers for large devices, while 1D wave propagation has to be assumed along the small capillaries of the central monolith reactor of the ATD; therefore, these monoliths can be replaced by plane wave four-pole transfer matrices from an acoustical point of view.^20–23 The mode matching method is then applied in subsection ‘Mode matching method’ to the compatibility conditions of the 3D acoustic fields at the ATD geometric discontinuities, leading to the computation of the complex wave amplitudes in all the subdomains involved and the corresponding transmission loss (TL). On the other hand, for large ATD inlet ducts such as those found in heavy-duty and off-road engines, the usual models with incident plane wave excitation are not accurate since the onset of higher order incident modes in the inlet duct is expected for the frequency range of interest. Therefore, a TL variation is likely to occur depending on these modes, similar to the results found in large dissipative silencers. ¹¹ Therefore, some details are given in subsection ‘Multimodal excitation’ for three different multimodal incident sound field hypotheses:^1,12 equal modal amplitude (EMA), equal modal power (EMP) and equal modal energy density (EMED). To conclude this section, the last part is devoted to a brief description of the acoustic models associated with sound propagation in the monolith of a catalytic converter as well as a particulate filter (see subsection ‘Acoustic modelling of the ATD monolith’). Results are presented in Section ‘Results and discussion’ for a number of configurations. First, a comparison is provided regarding the results associated with an ATD incorporating a catalytic converter monolith. These are computed by the proposed mode matching method and also through a FE formulation; in addition, experimental measurements are carried out in a laboratory set-up. As it will be shown, both modelling techniques exhibit an excellent agreement, with undistinguishable TL curves, while benchmarking against experimental attenuation is also very good. Then, further results are shown for two large axisymmetric configurations of ATD incorporating, on the one hand, a catalytic converter monolith, and on the other, a particulate filter monolith. TL computations are carried out with increasing number of higher order modes in the incident acoustic field; in addition, the influence of the mean flow Mach number is assessed in the case of the particulate filter for illustration purposes. As it will be shown, a relevant influence of the multimodal incident sound field on the sound attenuation is found for the test problems considered in the current investigation. To conclude, some final remarks are given in Section ‘Conclusions’.

Overview of the acoustic problem

Figure 1 shows a circular concentric ATD including a monolith. As shown in recent studies,^20–23 the sound attenuation of an exhaust device incorporating a monolith can be properly predicted if the monolith is replaced by a plane wave four-pole matrix providing a relationship between the acoustic fields at both sides of the capillary region. Therefore, the presence of higher order modes in the cylindrical inlet/outlet regions is combined with 1D wave propagation within the capillary ducts of the central monolith.

Figure 1.

Scheme of a circular concentric ATD incorporating a monolith. The latter is replaced by a transfer matrix to model its acoustic behaviour.

In all the rigid ducts involved (A, B, D and E), the acoustic pressure and velocity fields can be written in terms of a series expansion. For example, the solution of the wave equation in region A is given by^1–3

$$P_A\left(r,z\right)=\sum _{n=0}^{\infty }\left(A_n^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{A,n}z}+A_n^{-}{\mathrm{e}}^{\mathrm{j}{k}_{A,n}z}\right){\mathrm{J}}_0\left({\alpha }_{n}r/R_1\right)$$

(1)

$$U_A(r,z)=\frac{1}{{\rho }_0\omega }\sum _{n=0}^{\infty }{k}_{A,n}(A_n^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{A,n}z}-A_n^{-}{\mathrm{e}}^{\mathrm{j}{k}_{A,n}z}){\mathrm{J}}_0\left({\alpha }_{n}r/R_1\right)$$

(2)

where P_A and U_A are the acoustic pressure and axial velocity fields, $A_n^{+}$ and $A_n^{-}$ are the complex wave amplitudes (unknowns of the problem) associated with incident and reflected waves, ${\mathrm{J}}_0\left({\alpha }_{n}r/R_1\right)$ is the Bessel function of the first kind and order zero defining the transversal pressure mode (also called transversal modal function or eigenfunction along this work) corresponding to the circular cross-section under analysis (α_n are the roots satisfying the rigid wall boundary condition^1,2), $k_{A,n}$ is the axial wavenumber, ω is the angular frequency, ρ₀ is the density of the air and j is the imaginary unit. The axial wavenumber is given by

$$k_{A,n}=\pm \sqrt{k_0^2-{\left({\alpha }_n/R_1\right)}^2}$$

(3)

k ₀ = ω/c₀ being the wavenumber in air, with c₀ denoting the speed of sound. The signs associated with the axial wavenumbers are chosen to guarantee exponential attenuation of evanescent modes.^2–4

Mode matching method

The complete acoustic field of the system requires the computation of the wave amplitudes $A_n^{\pm }$ , $B_n^{\pm }$ , $D_n^{\pm }$ and $E_n^{\pm }$ in all the ducts. First, the compatibility conditions of the acoustic fields are taken into account, being applied at the interfaces between duct A and central chamber B (inlet expansion), duct D and central chamber E (outlet contraction), ²⁴ as well as the coupling of acoustic pressures and axial velocities at both sides of the monolith through the aforementioned plane wave four-pole transfer matrix. These can be expressed as detailed next. ²³

At the inlet expansion, the continuity of acoustic pressure and axial velocity are written as

$$P_A\left(r,z=0\right)=P_B\left(r,z=0\right)r\in S_A$$

(4)

$$U_A\left(r,z=0\right)=U_B\left(r,z=0\right)r\in S_A$$

(5)

while the rigid wall boundary condition at the left end plate is given by

$$U_B(r,z=0)=0r\in S_B-S_A$$

(6)

The acoustic coupling associated with the monolith at the interfaces S_B ≡ S_D leads to

$$\begin{array}{cc}\multicolumn{1}{c}{P_B(r,z=L_B)=}& T_{11}^m{P}_D(r,z\prime =0)+T_{12}^m{U}_D(r,z\prime =0)\hfill \\ r\in S_B\equiv S_D\hfill \end{array}$$

(7)

$$\begin{array}{cc}\multicolumn{1}{c}{U_B(r,z=L_B)=}& T_{21}^m{P}_D(r,z\prime =0)+T_{22}^m{U}_D(r,z\prime =0)\hfill \\ r\in S_B\equiv S_D\hfill \end{array}$$

(8)

or, in matrix form

$$\begin{array}{ccc}\multicolumn{1}{c}{\left\{\begin{array}{c}{P}_B(r,z=L_B)\hfill \\ \multicolumn{1}{c}{U_B(r,z=L_B)}\end{array}\right\}=}& \left[\begin{array}{cc}{T}_{11}^m\hfill & T_{12}^m\hfill \\ \multicolumn{1}{c}{T_{21}^m}& T_{22}^m\hfill \end{array}\right]\left\{\begin{array}{c}{P}_D(r,z\prime =0)\hfill \\ \multicolumn{1}{c}{U_D(r,z\prime =0)}\end{array}\right\}\hfill & r\hfill \\ r\in S_B\equiv S_D\hfill \end{array}$$

(9)

where $T_{11}^m$ , $T_{12}^m$ , $T_{21}^m$ and $T_{22}^m$ are the acoustic four-poles^1,2 associated with the monolith. It is worth noting that the characteristics of the monolith, including its length L_C (see Figure 1 for details), are implicitly considered in the four-poles. These will be explicitly detailed in subsection ‘Acoustic modelling of the ATD monolith’.

Finally, similar equations to those considered for the expansion can be written for the contraction (ducts D and E). The continuity of acoustic pressure and axial velocity gives

$$P_D\left(r,z\prime =L_D\right)=P_E\left(r,z\prime =L_D\right)r\in S_E$$

(10)

$$U_D\left(r,z\prime =L_D\right)=U_E\left(r,z\prime =L_D\right)r\in S_E$$

(11)

and the zero axial velocity at the right end rigid plate gives

$$U_D(r,z\prime =L_D)=0r\in S_D-S_E$$

(12)

To compute the wave amplitudes $A_n^{\pm }$ , $B_n^{\pm }$ , $D_n^{\pm }$ and $E_n^{\pm }$ in all the ducts, the mode matching method can be applied to equations (4)–(12). Basically, a weighted integration procedure is applied whose weighting functions are the transversal modal functions of the ducts.^{3–5,10,23,24} Hereinafter, N_m is the number of higher order modes considered to expand the acoustic pressure and velocity fields (see equations (1) and (2)) for regions B and D, associated with the central chambers, while N_a and N_e higher order modes are taken into account for the inlet and outlet ducts A and E, respectively. The algebraic system of equations required to compute the unknown wave amplitudes is generated next.^23,25

Starting with the expansion, equation (4) expressing continuity of pressure is multiplied by the transversal modal function ${\mathrm{J}}_0\left({\alpha }_{s}r/R_1\right)$ , with s = 0, 1, 2, …, N_a, and integrated over S_A. This leads to

$$\begin{array}{c}\multicolumn{1}{c}{{\int }_{S_A}{P}_A(r,z=0){\mathrm{J}}_0\left({\alpha }_{s}r/R_1\right)\mathrm{dS}}\\ \multicolumn{1}{c}{={\int }_{S_A}{P}_B(r,z=0){\mathrm{J}}_0\left({\alpha }_{s}r/R_1\right)\mathrm{dS}}\end{array}$$

(13)

The orthogonality properties of the eigenfunctions^1–4,23 lead, for s = 0, 1, 2, …, N_a, to

$$\begin{array}{c}\multicolumn{1}{c}{(A_s^{+}+A_s^{-}){\int }_{S_A}{\mathrm{J}}_0{\left({\alpha }_{s}r/R_1\right)}^2\mathrm{dS}=\sum _{n=0}^{N_m}(B_n^{+}+B_n^{-})}\\ {\int }_{S_A}{\mathrm{J}}_0\left({\alpha }_{n}r/R_2\right){\mathrm{J}}_0\left({\alpha }_{s}r/R_1\right)\mathrm{dS}\hfill \end{array}$$

(14)

The compatibility of the axial velocity at the expansion is given by equations (5) and (6); these are multiplied by the transversal modal function ${\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)$ , with s = 0, 1, 2, …, N_m. The former product is integrated over S_A and the latter over S_B– S_A. After summation yields

$$\begin{array}{c}\multicolumn{1}{c}{{\int }_{S_A}{U}_A(r,z=0){\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)\mathrm{dS}}\\ \multicolumn{1}{c}{={\int }_{S_B}{U}_B(r,z=0){\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)\mathrm{dS}}\end{array}$$

(15)

Similar to the previous case, for s = 0, 1, 2, …, N_m, the orthogonality properties of the transversal modal functions lead to

$$\begin{array}{c}\sum _{n=0}^{N_a}{k}_{A,n}(A_n^{+}-A_n^{-}){\int }_{S_A}{\mathrm{J}}_0\left({\alpha }_{n}r/R_1\right){\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)\mathrm{dS}\hfill \\ =k_{B,s}(B_s^{+}-B_s^{-}){\int }_{S_B}{\mathrm{J}}_0{\left({\alpha }_{s}r/R_2\right)}^2\mathrm{dS}\hfill \end{array}$$

(16)

The corresponding integrals in equations (14) and (16) can be computed analytically through the expressions detailed in the works,^3,4,23 leading to a considerable reduction of the computational effort and speeding up transmission loss calculations.

After considering the expansion, the application of the mode matching method focuses on the monolith. Thus, equations (7) and (8) are multiplied by the transversal modal function ${\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)$ , with s = 0, 1, 2, …, N_m. Integrating over S_B ≡ S_D, taking advantage of the orthogonality relations and removing common factors yield the following equations for s = 0, 1, 2, …, N_m

$$\begin{array}{c}{B}_s^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{B,s}{L}_B}+B_s^{-}{\mathrm{e}}^{\mathrm{j}{k}_{B,s}{L}_B}=T_{11}^m(D_s^{+}+D_s^{-})\hfill \\ +T_{12}^m\frac{k_{D,s}}{{\rho }_0\omega }(D_s^{+}-D_s^{-})\hfill \end{array}$$

(17)

$$\begin{array}{c}\frac{k_{B,s}}{{\rho }_0\omega }(B_s^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{B,s}{L}_B}-B_s^{-}{\mathrm{e}}^{\mathrm{j}{k}_{B,s}{L}_B})=T_{21}^m(D_s^{+}+D_s^{-})\hfill \\ +T_{22}^m\frac{k_{D,s}}{{\rho }_0\omega }(D_s^{+}-D_s^{-})\hfill \end{array}$$

(18)

As indicated in Denia et al., ²³ it is worth noting here that very simple algebraic expressions have been obtained, where neither integrations nor modal summations appear, relating directly wave amplitudes with equal modal number. These equations do not depend on the transversal cross section, provided that its geometry is axially uniform.

Finally, the equations and integrals associated with the outlet contraction between chamber D and duct E are obtained. Similar to the approach presented for the expansion, equations (10)–(12) associated with the contraction provide a final set of equations as follows. The pressure condition (10) is multiplied by the eigenfunction ${\mathrm{J}}_0\left({\alpha }_{s}r/R_3\right)$ , with s = 0, 1, 2, …, N_e, and integrated over S_E, leading to

$$\begin{array}{c}\sum _{n=0}^{N_m}(D_n^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{D,n}{L}_D}+D_n^{-}{\mathrm{e}}^{\mathrm{j}{k}_{D,n}{L}_D}){\int }_{S_E}{\mathrm{J}}_0\left({\alpha }_{n}r/R_2\right){\mathrm{J}}_0\left({\alpha }_{s}r/R_3\right)\mathrm{dS}\hfill \\ =(E_s^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{E,s}{L}_D}+E_s^{-}{\mathrm{e}}^{\mathrm{j}{k}_{E,s}{L}_D}){\int }_{S_E}{\mathrm{J}}_0{\left({\alpha }_{s}r/R_3\right)}^2\mathrm{dS}\hfill \end{array}$$

(19)

Note that the orthogonality properties have already been applied in the right hand side of equation (19). In addition, the velocity conditions (11) and (12) provide, after multiplication by ${\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)$ , integration and summation, the following expression

$$\begin{array}{c}{k}_{D,s}(D_s^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{D,s}{L}_D}-D_s^{-}{\mathrm{e}}^{\mathrm{j}{k}_{D,s}{L}_D}){\int }_{S_D}{\mathrm{J}}_0{\left({\alpha }_{s}r/R_2\right)}^2\mathrm{dS}\hfill \\ =\sum _{n=0}^{N_e}{k}_{E,n}(E_n^{+}{\mathrm{e}}^{-\mathrm{j}{k}_{E,n}{L}_D}-E_n^{-}{\mathrm{e}}^{\mathrm{j}{k}_{E,n}{L}_D}){\int }_{S_E}{\mathrm{J}}_0\left({\alpha }_{n}r/R_3\right){\mathrm{J}}_0\left({\alpha }_{s}r/R_2\right)\mathrm{dS}\hfill \end{array}$$

(20)

with s = 0, 1, 2, …, N_m, and where, again, orthogonality has been considered. Also, the integrals in equations (19) and (20) can be computed analytically to reduce the computational expenditure. The reader is referred to the works^3,4,23 for further details.

Multimodal excitation

Transmission loss can be evaluated to determine the acoustic attenuation of the aftertreatment device. ² First, the modal amplitudes $A_n^{\pm }$ , $B_n^{\pm }$ , $D_n^{\pm }$ and $E_n^{\pm }$ have to be computed through equations (13)–(20), which provide an algebraic system with (N_a + 1) + 4 (N_m + 1) + (N_e + 1) equations. The additional conditions for TL computation are associated with an anechoic termination, so $E_n^{-}=0$ $\forall n$ , and with the excitation by means of the incident wave amplitudes $A_n^{+}$ . Regarding the latter, as indicated previously, plane wave propagation in the inlet/outlet ducts of large ATD no longer holds. Thus, a procedure to account for multimodal excitation is described here. This can be summarized in five main steps as follows: step 1 is associated with the computation of the number of propagating modes in each duct for a given excitation frequency. For the inlet, this includes the plane wave mode (always present) and the non-planar excitation modes up to of N_a, obtained through the condition

$$k_{A,n}=\sqrt{k_0^2-{\left({\alpha }_n/R_1\right)}^2}\ge 0\mathrm{for}n=0,1,...,N_a$$

(21)

Note that N_a can vary during an ATD acoustic computation, increasing as the excitation frequency exceeds the cut-on frequencies of the successive transverse higher order modes. ² Step 2 requires the computation of the excitation wave amplitudes $A_n^{+}$ of the inlet duct.^1,11,12 If only a plane wave (PW) is included in the analysis, this yields

$$A_0^{+}=1;A_n^{+}=0forn>0$$

(22)

Alternative mode mixtures include equal modal amplitude (EMA), equal modal power (EMP) and equal modal energy density (EMED), the corresponding expression being as follows ¹ :

$$A_n^{+}=1\mathrm{for}n\le N_a(\mathrm{EMA})$$

(23)

$$A_n^{+}=\sqrt{\left({\delta }_n^2/N_a\right)\left(1/\mathrm{Re}\left(k_{A,n}/k_0\right)\right)}\mathrm{for}n\le N_a(\mathrm{EMP})$$

(24)

$$A_n^{+}=\sqrt{\left({\delta }_n^2/\sum _j^{N_a}\mathrm{Re}\left(k_{A,j}/k_0\right)\right)}\mathrm{for}n\le N_a(\mathrm{EMED})$$

(25)

where δ_n depends on the integration of the transversal modal functions over the cross section. ¹ Once the excitation wave amplitudes $A_n^{+}$ are known, the mode matching system of equations is solved in step 3. Then, step 4 implies the computation of the acoustic incident/transmitted power at the end sections of the inlet/outlet ducts ¹

$$\begin{array}{c}\multicolumn{1}{c}{\mathrm{Po}{t}_{\mathrm{inc}}={\int }_S\frac{1}{2}\mathrm{Re}\left(P_{\mathrm{inc}}{U}_{\mathrm{inc}}^{*}\right)\mathrm{dS}=\frac{\pi }{{\rho }_0\omega }\sum _{n=0}^{N_a}{k}_{A,n}{\left|A_n^{+}\right|}^2}\\ {\int }_0^{R_1}{\mathrm{J}}_0{\left({\alpha }_{n}r/R_1\right)}^{2}r\mathrm{dr}\hfill \end{array}$$

(26)

$$\begin{array}{c}\multicolumn{1}{c}{\mathrm{Po}{t}_{\mathrm{trans}}={\int }_S\frac{1}{2}\mathrm{Re}\left(P_{\mathrm{trans}}{U}_{\mathrm{trans}}^{*}\right)\mathrm{dS}=\frac{\pi }{{\rho }_0\omega }\sum _{n=0}^{N_e}{k}_{E,n}{\left|E_n^{+}\right|}^2}\\ {\int }_0^{R_3}{\mathrm{J}}_0{\left({\alpha }_{n}r/R_3\right)}^{2}r\mathrm{dr}\hfill \end{array}$$

(27)

where Re denotes the real part of a complex number and * is used for the complex conjugate. Finally, in step 5 the transmission loss is evaluated as ²

$$\mathrm{TL}\left(\mathrm{dB}\right)=10\log \left(\frac{\mathrm{Po}{\mathrm{t}}_{\mathrm{inc}}}{\mathrm{Po}{\mathrm{t}}_{\mathrm{trans}}}\right)$$

(28)

Acoustic modelling of the ATD monolith

Figure 2 shows examples of monoliths and schemes of the associated flow pattern. Catalytic converter capillaries have upstream and downstream ends open, while particulate filter capillaries are characterized by channel ends alternatively plugged to force the flow through the porous walls.^17–21

Figure 2.

Monolith and flow pattern: (a) catalytic converter (CC) and (b) particulate filter (PF).

Both types of ATD can be modelled by means of a transfer matrix T.^18,26 In particular, for the catalytic converter, the transfer matrix T can be expressed as^18,22

$$\begin{array}{c}\multicolumn{1}{c}{\mathbf{T}=\left(\begin{array}{cc}{T}_{11}^m\hfill & T_{12}^m\hfill \\ \multicolumn{1}{c}{T_{21}^m}& T_{22}^m\hfill \end{array}\right)=\left(\begin{array}{cc}\cos \left(k_{\mathrm{cc}}{L}_C\right)\hfill & \frac{\mathrm{j}{\rho }_{\mathrm{cc}}{c}_{\mathrm{cc}}\sin \left(k_{\mathrm{cc}}{L}_C\right)}{{\varphi }_{\mathrm{cc}}}\hfill \\ \multicolumn{1}{c}{\frac{\mathrm{j}{\varphi }_{\mathrm{cc}}\sin \left(k_{\mathrm{cc}}{L}_C\right)}{{\rho }_{\mathrm{cc}}{c}_{\mathrm{cc}}}}& \cos \left(k_{\mathrm{cc}}{L}_C\right)\hfill \end{array}\right)}\end{array}$$

(29)

where ϕ_cc is the monolith porosity, L_C the length of the monolith capillary ducts (see Figure 1) and $k_{\mathrm{cc}}=\omega /c_{\mathrm{cc}}$ the wavenumber. The effective density ρ_cc and speed of sound c_cc can be expressed respectively as

$${\rho }_{\mathrm{cc}}={\rho }_0(1+\frac{\mathrm{Re}{s}_{\mathrm{cc}}{\varphi }_{\mathrm{cc}}}{\mathrm{j}\omega {\rho }_0}{G}_c\left(s_w\right))$$

(30)

and

$$c_{\mathrm{cc}}=\frac{c_0}{\sqrt{(1+\frac{\mathrm{Re}{s}_{\mathrm{cc}}{\varphi }_{\mathrm{cc}}}{\mathrm{j}\omega {\rho }_0}{G}_c\left(s_w\right))(\gamma -(\gamma -1)F)}}$$

(31)

Res_cc being the monolith steady flow resistivity, γ is the ratio of specific heats and s_w is the shear wave number calculated as

$$s_w=\alpha \sqrt{\frac{8\omega {\rho }_0}{\mathrm{Re}{s}_{\mathrm{cc}}{\varphi }_{\mathrm{cc}}}}$$

(32)

where α depends on the geometry of the capillary cross section. Also, F in equation (31) is expressed as

$$F=\frac{1}{1+\frac{\mathrm{Re}{s}_{\mathrm{cc}}{\varphi }_{\mathrm{cc}}}{\mathrm{j}\mathrm{Pr}\omega {\rho }_0}{G}_c\left(\sqrt{\mathrm{Pr}}{s}_w\right)}$$

(33)

where $\mathrm{Pr}=\mu C_p/\kappa$ is the Prandtl number, μ the viscosity, C_p the specific heat at constant pressure and κ the thermal conductivity. ¹⁸ To conclude, in equations (30), (31) and (33) G_c(s_w) is given by

$$\begin{array}{c}\multicolumn{1}{c}{G_c\left(s_w\right)=\frac{-\frac{s_w}{4}\sqrt{-\mathrm{j}}\frac{{\mathrm{J}}_1\left(s_w\sqrt{-\mathrm{j}}\right)}{{\mathrm{J}}_0\left(s_w\sqrt{-\mathrm{j}}\right)}}{1-\frac{2}{s_w\sqrt{-\mathrm{j}}}\frac{{\mathrm{J}}_1\left(s_w\sqrt{-\mathrm{j}}\right)}{{\mathrm{J}}_0\left(s_w\sqrt{-\mathrm{j}}\right)}}}\end{array}$$

(34)

where J₀ and J₁ are Bessel functions of the first kind and zeroth and first order. Further information of the acoustic model can be found in Selamet et al. ¹⁸

Regarding the particulate filter monolith, only some general aspects are presented in this work for the sake of brevity (the reader is referred to Allam and Åbom ¹⁹ for further details). As indicated previously, the acoustic modelling assumes plane wave propagation in the capillaries, which can be computed through transfer matrices involved in equations (7)–(9), (17) and (18), relating the acoustic fields at both sides of the monolithic region. For the case of a particulate filter, if acoustic pressures and volume flows are considered,^19,21 this yields

$$\left\{\begin{array}{c}{P}_B\left(r,z=L_B\right)\\ S_B{U}_B\left(r,z=L_B\right)\end{array}\right\}=T\left\{\begin{array}{c}{P}_D\left(r,z\prime =0\right)\\ S_D{U}_D\left(r,z\prime =0\right)\end{array}\right\}r\in S_B\equiv S_D$$

(35)

where the corresponding transfer matrix T can be expressed as

$$\mathbf{T}=\left(\begin{array}{cc}{T}_{11}^m\hfill & T_{12}^m\hfill \\ \multicolumn{1}{c}{T_{21}^m}& T_{22}^m\hfill \end{array}\right)={\mathbf{T}}_{\mathbf{IN}}{\mathbf{T}}_{\mathbf{I}}{\mathbf{T}}_{\mathbf{II}}{\mathbf{T}}_{\mathbf{III}}{\mathbf{T}}_{\mathbf{OUT}}$$

(36)

T_IN and T_OUT being the four-pole transfer matrices corresponding to the area changes between region B and the monolith inlet, and the monolith outlet and region D, respectively (see Figure 1), T_I and T_III are the transfer matrices associated with the short channels with impermeable walls located at the inlet and outlet of the filter, with area S and length l_plug corresponding to that of the plugs of the adjacent capillaries (see grey plugs in Figure 2(b)), and finally T_II is the transfer matrix of the particulate filter central section. ¹⁹ It is worth noting that, since the length l_plug is usually very short (less than 0.01 m), it is possible to simply add its acoustic effect as an end correction (mass plug) to the IN and OUT sections. In the presence of mean flow with small Mach number, T_IN+I and T_III+OUT can expressed respectively as ¹⁹

$$\begin{array}{c}\multicolumn{1}{c}{{\mathbf{T}}_{\mathbf{IN}\mathbf{+}\mathbf{I}}=\left(\begin{array}{cc}1\hfill & \frac{{\rho }_0{c}_0{M}_{\mathrm{IN}}}{S_B}(\frac{1}{m_{\mathrm{IN}}^2}-1)+\mathrm{j}\frac{{\rho }_0\omega l_{\mathrm{plug}}}{NS}\hfill \\ \multicolumn{1}{c}{0}& 1\hfill \end{array}\right)}\end{array}$$

(37)

$$\begin{array}{c}\multicolumn{1}{c}{{\mathbf{T}}_{\mathbf{III}\mathbf{+}\mathbf{OUT}}=\left(\begin{array}{cc}1\hfill & 2\frac{{\rho }_0{c}_0{M}_{\mathrm{OUT}}}{S_D}(1-\frac{1}{m_{\mathrm{OUT}}})+\mathrm{j}\frac{{\rho }_0\omega l_{\mathrm{plug}}}{NS}\hfill \\ \multicolumn{1}{c}{0}& 1\hfill \end{array}\right)}\end{array}$$

(38)

where m_IN and m_OUT are the open area ratios at the inlet and outlet, M_IN and M_OUT are the mean flow Mach numbers in the central chambers B and D, respectively and N is the corresponding number of open channels (approximately half the total number of channels N_t of the particulate filter monolith).

Obtaining the transfer matrix T_II of the particulate filter central section is the most laborious part of the model presented in Allam and Åbom. ¹⁹ The convective acoustic wave equations for two neighbouring narrow channels are solved, including viscous and thermal losses. Calculation of the propagation constants is carried out through an eigenvalue problem that is solved by iteration, using the no-flow roots as starting values (note that Mach number in typical filters is very small, M < 0.1). Finally, from the roots and the associated eigenvectors, the sound field in the filter section can be computed, leading to the acoustic transfer matrix. In the present work, the model of Allam and Åbom ¹⁹ has been implemented, and the corresponding results will be shown in subsections ‘Geometry 1: N_a = 1’ and ‘Geometry 2 (doubled radii): N_a = 4’ for some particular configurations.

Regarding the effect on the attenuation of the presence of mean flow and the corresponding gradients, Allam and Åbom ¹⁹ indicate that the influence is very small. From a practical point of view, the gradients can be replaced by the use of the average value, giving an error of the order of 0.3 dB for the monolith attenuation. Also, the main impact of the mean flow is related to the area changes in the inlet and outlet, modelled through matrices T_IN and T_OUT defined in equations (37) and (38).

Results and discussion

Validation

Figure 3 shows some pictures of a modular catalytic converter prototype with circular cross section. The main radial dimensions are R₁ = R₃ = 0.0258 m for the inlet/outlet ducts and R₂ = 0.1275 m for the chambers and monolith. Regarding the lengths, the values L_B = L_D = 0.0475 m and L_C = 0.075 m are considered in the computations and experimental measurements (see Figure 1 for details).

Figure 3.

Catalytic converter prototype: (a) external view and (b) open prototype and internal structure.

Cold flow hypotheses with T = 25°C are retained through the computations and laboratory conditions, the properties for the air being c₀ = 346.1 m/s and ρ₀ = 1.184 kg/m³. Regarding the monolith acoustic model, the following values are considered: resistivity Res_cc = 1500 Pa s/m², porosity ϕ_cc = 0.88, geometrical factor α = 1.14 (triangular metallic capillaries), dynamic viscosity μ = 1.837 × 10⁻⁵ Pa s, thermal conductivity κ = 0.0255 W/(m K) and specific heat C_p = 1006.4J/(kg K). The detailed model for the four-pole matrix computation can be found in the bibliography.^18,22

Figure 4 shows a comparison of the results computed by the proposed mode matching method and also through a FE formulation,^18,22 as well as experimental measurements carried out in the facilities of the research centre with a set-up based on the transfer function method. ² As it can be observed, both modelling techniques exhibit an excellent agreement, with undistinguishable TL curves in all the frequency range, while benchmarking against experimental attenuation is also very good, although the usual discrepancies appear due, among other things, to inaccuracies of the experimental measurements for long wavelengths at low frequencies (ineffective anechoic termination, insufficient distance between microphones, etc.), geometric deviations of the tested prototype and approximations associated with the acoustic model of the monolith.

Figure 4.

TL of an ATD incorporating a catalytic converter monolith: ___, mode matching method; +++, FEM; o o o, experimental measurement.

It is worth noting that, in the previous results, plane waves are required at the microphone locations to guarantee a correct experimental measurement, the first cut-on frequency being 8180.8 Hz in the inlet/outlet ducts. ² Multimodal propagation conditions exist, however, in the central chambers B and D, with one higher order mode (N_m = 1 with a cut-on frequency of 1688.6 Hz) in the frequency range considered in Figure 4. Experimental measurement of multidimensional incident modal fields at the inlet duct is a challenge beyond the scope of the current work. As indicated in the literature, to investigate the effect of mode mixture in the incident wave, it is convenient to start by using theoretical models, as these allow one to remove experimental uncertainty in order to focus on the physical processes involved. ¹⁴ Therefore, the conclusions obtained to date await more accurate representations of the sound source, as well as supporting experimental evidence, ¹⁵ the latter not yet widely available in the literature, at least to the best of the authors’ knowledge.

Acoustic characteristics of the monolith: catalytic converter and particulate filter

Regarding the catalytic converter, two axisymmetric configurations are considered for computation purposes. In the first case, geometry 1 has inlet and outlet ducts defined by the radii R₁ = R₃ = 0.125 m, while the chambers and monolith are characterized by the radius R₂ = 0.25 m and the lengths L_B = L_D = 0.1475 m and L_C = 0.175 m (see Figure 1 for details). In the second case, geometry 2 is defined by the dimensions R₁ = R₃ = 0.25 m, R₂ = 0.5 m, L_B = L_D = 0.475 m and L_C = 0.75 m. It is worth noting that in this latter case radii have been doubled and thus more higher order modes will propagate. Similar to subsection ‘Validation’, cold flow hypotheses with T = 25°C are retained through the computations, the properties for the air being c₀ = 346.1 m/s and ρ₀ = 1.184 kg/m³. Regarding the monolith acoustic model, the following values are assumed: resistivity Res_cc = 1500 Pa s/m², porosity ϕ_cc = 0.88, geometrical factor α = 1.14, dynamic viscosity μ = 1.837 × 10⁻⁵ Pa s, thermal conductivity κ = 0.0255W/(m K) and specific heat C_p = 1006.4 J/(kg K). The detailed model for the four-pole matrix computation can be found in the bibliography.^18,22

As far as the particulate filter is concerned, the dimensions of geometries 1 and 2 are kept. For the model presented in Allam and Åbom, ¹⁹ the following values and parameters are taken into account: plug length l_plug = 0.005 m, total number of channels (per square metre) N_t = 3.1 × 10⁵, channel width d_h = 1.44 × 10⁻³ m, wall thickness h_t = 3.55 × 10⁻⁴ m, wall permeability σ_w =2.5 × 10⁻¹³ m², thickness of soot loading h_soot =h_t/10 and soot permeability σ_soot = 1.5 × 10⁻¹⁴ m². Also, the mean flow Mach number values M_IN = 0, M_IN = 0.01 and M_IN = 0.02 are considered in region B for illustration purposes (the corresponding values in the inlet duct A are M = 0, M = 0.04 and M = 0.08, respectively) ¹⁹ . Finally, the values for temperature, dynamic viscosity, thermal conductivity and specific heat are the same as for the previous case of the catalyst.

Geometry 1: N_a = 1

For the conditions associated with this case, the incident modal field has a maximum of one higher order mode (N_a = 1) in the inlet duct, its cut-on frequency being 1688.6 Hz, ² while propagation of up to four higher order modes (N_m = 4) is expected in the central chambers B and D. As it can be seen in Figure 5, the mode mixture rules modify the attenuation of the catalytic converter above this cut-on frequency value; for this particular problem, when incident energy propagates in higher order modes this seems to have a positive impact on the acoustic performance. As observed by Kirby and Lawrie ¹¹ for large dissipative silencers, the three different multimode excitations considered for the incident acoustic field yield comparable results in general, although EMA predictions are lower than both EMED and EMP predictions. For this configuration, the maximum attenuation difference found between the multimodal excitation models and the PW prediction is 3.2 dB (about 30% of discrepancy) at approximately 1690 Hz. The presence of acoustic energy dissipation within the capillary ducts of the monolith leads to a relatively high attenuation, without the undesirable pass bands with low sound attenuation found in purely reactive configurations. In addition, the acoustic behaviour of the device may exhibit resonant characteristics and therefore the ATD transmission loss can benefit from some high attenuation peaks at particular frequencies.

Figure 5.

TL of an ATD incorporating a catalytic converter monolith for different modal excitations: ___, PW; – – –, EMA; – · – ·, EMP; ·········, EMED.

Figure 6 shows the impact of the mode mixture rules on the attenuation of a particulate filter with the same dimensions as the previous catalytic converter (geometry 1) and the properties defined in subsection ‘Acoustic characteristics of the monolith: catalytic converter and particulate filter’. Once again, the differences of acoustic behaviour above the cut-on frequency and the importance of a correct modelling of the propagation conditions are highlighted. It is worth noting that a change in the trend is found for this device, since now EMA predictions are higher than EMED and EMP results over a wide frequency range.

Figure 6.

TL of an ATD incorporating a particulate filter for different modal excitations: ___, PW; – – –, EMA; – · – ·, EMP; ·········, EMED.

Figure 7 shows the TL curves for the previous configuration of particulate filter, considering EMED mixture and three different values of mean flow Mach number in region B given by M_IN = 0, M_IN = 0.01 and M_IN = 0.02, ¹⁹ the corresponding values in the inlet duct A being M = 0, M = 0.04 and M = 0.08, respectively. As can be observed, the presence of mean flow leads to a slight TL increase due to the acoustic energy dissipation associated with the monolith expansion/contraction as well as its porous walls (note that the flow is forced to pass through the permeable walls of the particulate filter). This trend is consistent with the results shown in Allam and Åbom. ¹⁹

Figure 7.

Influence of Mach number on TL of an ATD incorporating a particulate filter for EMED mixture: ___, M_IN = 0; – – –, M_IN = 0.01; – · – ·, M_IN = 0.02.

Geometry 2 (doubled radii): N_a = 4

In this case, the incident acoustic field has a maximum of four higher order modes (N_a = 4), the cut-on frequencies being 844.3 Hz, 1545.8 Hz, 2241.6 Hz and 2935.7 Hz respectively. ² Now, propagation of up to eight higher order modes (N_m = 8) takes place in the central chambers B and D. As shown in Figure 8 for a catalytic converter, the main trends observed in Figure 5 are also found here, the same comments being applicable in general. Anyway, it is worth noting that the larger dimensions of the monolith increase the acoustic energy dissipation, thus leading to higher attenuation levels compared with Figure 5. Now, the relative discrepancies among the TL given by the different models seem to have decreased slightly in the high frequency range. In particular, the maximum attenuation difference found between the multimodal excitation models and the PW prediction is 6.7 dB (about 20% of discrepancy) at approximately 1830 Hz.

Figure 8.

TL of an ATD incorporating a catalytic converter monolith for different modal excitations and doubled radii: ___, PW; – – –, EMA; – · – ·, EMP; ·········, EMED.

Finally, Figure 9 shows the attenuation curves for a particulate filter. Again, it can be observed the evident influence of the mode mixture rules, but now there is some alternation in the frequency intervals where EMED and EMP have a higher TL than EMA, and vice versa. As expected, the larger dimensions of this filter compared with the one associated with Figure 6 lead to higher TL values. Also, as for the case of the catalytic converter (Figure 8), the relative differences among the acoustic attenuation at high frequencies given by the different approaches seem to have decreased slightly compared with geometry 1.

Figure 9.

TL of an ATD incorporating a particulate filter for different modal excitations and doubled radii: ___, PW; – – –, EMA; – · – ·, EMP; ·········, EMED.

Conclusions

An accurate and efficient computational method has been presented for the acoustic modelling of large ATD with a monolith, considering both a catalytic converter as well as a particulate filter. The model is based on the analytical mode matching method and assumes 3D wave propagation in the inlet/outlet regions (ducts and chambers) and 1D wave propagation along the small capillaries of the central monolith reactor. In large ATD, plane wave (PW) propagation in the inlet/outlet ducts no longer holds for the frequency range of interest, and energy incident on the device propagating in higher order modes may have a strong impact on the acoustic attenuation performance. The specification of the incident modal field requires a knowledge of the source characteristics, but this is currently an open field of research and there is a lack of accurate representations as well as experimental validation. In addition to the PW model, several mode mixtures have been implemented: equal modal amplitude (EMA), equal modal power (EMP) and equal modal energy density (EMED). As expected, a strong influence of the mode mixture technique on the attenuation has been found. The three different multimode conditions adopted for the incident sound field can yield predictions that are consistent with the literature for some particular configurations of catalytic converter, in the sense that at higher frequencies the acoustic performance of the exhaust device generally improves when energy from the source is carried by higher order modes (especially if dissipative effects are present in the analysis, such as those taking place within the capillary ducts of the monolith). However, the opposite trend can also be found for other configurations of particulate filter and thus, a plausible conclusion is that an accurate transmission loss estimation in large ATD requires the consideration of multimodal excitation, but the results are strongly dependent on the particular configuration under analysis and the mode mixture rule chosen. Although literature seems to suggest that EMED is probably the most suitable model, this depends on the noise source characteristics and further work is required.