Journal of Sound and Vibration 332 (2013) 4560–4571
Contents lists available at SciVerse ScienceDirect
Journal of Sound and Vibration journal homepage: www.elsevier.com/locate/jsvi
A simple absorbing layer implementation for transmission line matrix modeling Gwenaël Guillaume, Judicaël Picaut n LUNAM Université, IFSTTAR, AME, LAE, CS 4, F-44344 Bouguenais, France
a r t i c l e in f o
abstract
Article history: Received 14 November 2011 Received in revised form 22 March 2013 Accepted 1 April 2013 Handling Editor: R.E. Musafir Available online 1 May 2013
An absorbing layer formulation for transmission line matrix modeling is proposed. The approach consists in attenuating the incident pulse propagating toward the absorbing layer only, using an attenuation factor which gradually decreases as the sound wave propagates along the absorbing medium. The formulation of the damping function followed by the attenuation factor along the absorbing layer is depicted and discussed. The efficiency of the present formulation is validated by comparison with another absorbing layer model and virtual boundary conditions proposed in the literature. Numerical simulations are also given in order to evaluate the effects of both the attenuation factor and the depth on the absorbing layer efficiency. As expected, results are consistent with absorbing layer implementation in other numerical methods; firstly, the attenuation at the entrance of the absorbing layer must be gentle, and secondly the efficiency increases with the layer depth. Lastly, it is shown that the unwanted reflection seems to vanish over the time when increasing the layer depth, meaning that reflections continuously occur within the absorbing layer and not on the geometrical limits of the absorbing layer. Although the approach is dedicated to outdoor sound propagation modeling (only an example on urban acoustics application is given), the proposed formulation of absorbing layers can be applied in other domains of acoustics. However, its application in shielded areas should be avoided because unwanted reflections due to an insufficient attenuation can be significant in comparison with the ambient noise in such quiet environments. & 2013 Elsevier Ltd. All rights reserved.
1. Introduction The growing increase of noise annoyance in cities and their impact on resident health make major interest in the knowledge of sound propagation mechanisms in urban areas. Many numerical approaches can be employed for acoustical modeling. Most of them are limited considering moving and time-dependent sound sources and time-varying parameters like micrometeorological conditions. Thus, time-domain methods have attracted a lot of attention during the last decade [1,2], due to the increase of computational resources and, overall, to their ability to model more realistic propagation conditions. Among these, the Transmission Line Matrix (TLM) method consists in an inherent discrete representation of wave propagation. Initially developed by Johns and Beurle in electromagnetism [3], the TLM method was extended to acoustics by Saleh and Blanchfield in the beginning of the 1990s [4]. Recently, a few applications for outdoor sound propagation was proposed [5,6], in particular for sound propagation over a porous ground [7–9]. Tsuchiya showed that the frequency-dependency of the atmospheric absorption can be taken into account including digital filters into the TLM numerical scheme [10]. Dutilleux
n
Corresponding author. Tel.: +33 2 40 84 57 89; fax: +33 2 40 84 59 92. E-mail address:
[email protected] (J. Picaut).
0022-460X/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.jsv.2013.04.003
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
4561
introduced wind-induced sound speed gradients through the definition of the effective sound speed, which allows to account for micrometeorological effects [11]. Nevertheless, difficulties are encountered when modeling free-space propagation due to the lack of efficient absorbing boundary conditions. Thus, the computational domain is usually enlarged for outdoor applications in order to avoid unwanted reflections. A powerful direct TLM implementation of Bérenger's perfectly matched layer (PML) [12] was proposed by Dubard and Pompéi within the electromagnetism context [13,14]. PMLs were applied to most of the numerical methods in acoustics [15], but no similar approach is available for TLM modeling of sound propagation yet. Instead, absorbing conditions for TLM models in acoustics are described by means of absorbing boundaries [16,17,8] and absorbing layers [18,8], but their efficiency as well as their implementation is still not satisfactory. As an alternative approach, a simple implementation of the absorbing layers concept, inspired by De Cogan et al. [18], is proposed in this paper, based on an intuitive formulation of absorbing layers. No rigorous PML implementation is formulated because the particle velocity would then be required, which can be avoided for TLM modeling in acoustics. Firstly, the implementation of absorbing conditions that are proposed in the literature is described and an alternative formulation is given (Section 2). The efficiency of the proposed formulation is then evaluated in Section 3 through a quantitative comparison with other artificial absorbing conditions, for several parameters of the absorbing layers and for several sound waves incident angles on the computational domain limit. The efficiency of the proposed formulation for long duration simulations and for realistic environment modeling is also discussed. 2. Absorbing boundary conditions for TLM modeling 2.1. TLM method principle Based on Huygens' principle, the TLM method consists of physically modeling undulatory phenomena through both a spatial and a temporal inherent discretization. Each volume element of the discrete propagation medium is represented by node exchanging pressure pulses with its neighbors through transmission lines. Thus, the discrete propagation medium can be seen as a transmission-line network linking nodes to each other. Inhomogeneities and dissipation in the propagation medium are contained in each node, and connection laws between neighboring nodes determine the sound propagation in the medium. In practice, these connection laws express, in a two-dimensional problem, the relation between scattered pulses v t Sði 7 1;j 7 1Þ leaving nodes ði7 1; j 7 1Þ neighboring to the node (i,j) by the transmission line v, at time t, and incident pulses u tþΔt I ði;jÞ
coming through the transmission line u to the node (i,j) at time t þ Δts (Fig. 1). Each node is surrounded by four neighboring nodes, meaning that u and v are equal to 1, 2, 3 or 4, in accordance with the pulse direction, +x, −x, +y and −y respectively, when considering incident pulses, and −x, +x, −y and +y for scattered ones [6]. Transmission lines are defined by a characteristic impedance Z ¼ ρ0 c, where ρ0 is the mass density of air and c ¼ Δl=Δt is the velocity in the network defined by the spatial step Δl and on the time increment Δt between two successive nodes (see Fig. 1). In order to ensure that sound propagation occurs with an apparent velocity equal to the classical adiabatic speed pffiffiffiof sound in air c0, it is necessary to correct the celerity c introduced in the TLM approach, by considering for instance c ¼ 2c0 pffiffiffi in a homogenous and non-dissipative medium, meaning that the characteristic impedance is changed to Z ¼ 2ρ0 c0 . 2.2. TLM artificial absorbing boundary conditions 2.2.1. Absorbing boundaries Absorbing boundaries consist of damping terminations at the limit of the computational domain in order to avoid unwanted reflections in the propagation domain of interest. Such boundary conditions should require less computational resources than the absorbing matched layers (see Section 2.2.2), since they are only applied at the limit of the propagation
Fig. 1. Incident tþΔt I uði;jÞ and scattered t Svði 7 1;j 7 1Þ pulses at a node (i,j) in the transmission-line network. Representation of the absorbing layer at the right side limit of the computational domain.
4562
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
medium and not within an artificial layer extending the propagation domain. Boundary operators. El-Masri et al. proposed an approach based on the method developed by Saguet in electromagnetism [19,16]. It consists in determining the acoustic pressure that should be reflected by the boundary, and, then, in removing the corresponding pressure from the effective sound field [20]. In this method, the acoustic pressure on the computational domain limit is obtained by means of Taylor's series expansion. As presented later in Section 3.3.2, this approach induces instabilities. Saguet assigned them to the presence of micro-pulses within the propagation medium after the wavefront passing, which does not fulfill the pressure field continuity condition. These micro-pulses correspond to spurious evanescent and propagating modes due to the periodical spatial sampling imposed by the mesh [21]. Thus, Saguet recommends to apply a pressure reflection coefficient equal to zero on the boundary after the wavefront reflection to avoid this annoying effect [22]. Consequently, this absorbing boundary condition cannot be used with continuous noise sources. Portí and Morente considered Higdon's operators [23] that allow the optimization of the boundary for several incidence angles [17]. This technique gives satisfying results provided that boundaries are located far from diffusing structures. However, this approach is efficient only for some incidence angles and around the design frequency. Moreover, problems of stability also appear using these kinds of artificial boundary conditions and digital filtering is then required [24]. In addition, Chien et al. showed that Higdon's condition and the one based on Taylor's expansion can be mathematically expressed by the same equation [25]. Subsequently, only the first boundary operators approach will be considered. Real impedance boundary condition. As suggested by Hofmann and Heutschi [8], an absorbing boundary can be implemented by identifying the boundary with a purely real impedance condition. The authors indicate that the choice of an impedance condition such as Z ¼ Z 0 (Z 0 ¼ ρ0 c0 being the characteristic impedance of air) leads to unsatisfying results pffiffiffi for oblique incidences. They advice in this case us to use Z ¼ 2Z 0 (i.e. the nominal impedance of transmission lines). However, this last condition fails for sound waves at normal incidence. As with Higdon's formulation, the adaptation of the impedance condition requires to take the angle of incidence into account, which makes this approach inapplicable for reliable simulations.
2.2.2. Absorbing matched layers The perfectly matched layers method, proposed by Bérenger, consists in surrounding the propagation domain with an anisotropic medium which gradually attenuates the sound waves [12]. Nowadays, this technique is widely used in numerical skills in order to avoid the retro-propagation of unwanted reflections toward the medium of interest. Many applications of this method were carried out for time-domain finite-difference (FDTD) modeling [26–28]. No rigorous implementation of Bérenger's perfectly matched layers was proposed for TLM method applications in acoustics. Indeed, even if the TLM model can be related to the finite-difference approach [9], the modeling of a PML in a TLM model starting from the PML equations is not straightforward, since the TLM method relies on connection laws between nodes instead of discrete equations. Alternatively, two methods can be found in the literature. Dissipative medium. Hofmann and Heutschi suggested to consider a dissipative medium with a gradual dissipative term in order to create a straightforward absorbing layer [8]. According to their analysis, the layer depth should be generously dimensioned, i.e. about 500 nodes, which represents about 50 wavelengths λ for simulations defined with an usual spatial step of Δl ¼ λ=10. In addition to the substantial computational burden of this technique, qualitative studies within the framework of the present study have shown that the practical implementation of this absorbing layer is not easy because of the gradual modification of the nodal reflection and transmission coefficients due to the progressive variation of the dissipative term. This approach is no longer considered in the following. Matched connection laws. De Cogan et al. proposed to modify the classical nodal connection laws of the TLM model within the absorbing layer, by introducing an attenuation factor increasing gradually as the sound wave propagates along this absorbing medium. This method is called graded perfectly matched load [18]. Thus, the two-dimensional connection laws within the absorbing layer are expressed as (Fig. 1) 1 tþΔt I ði;jÞ
¼ F ði;jÞ t S2ði−1;jÞ ;
(1a)
2 tþΔt I ði;jÞ
¼ F ði;jÞ t S1ðiþ1;jÞ ;
(1b)
3 tþΔt I ði;jÞ
¼ F ði;jÞ t S4ði;j−1Þ ;
(1c)
4 tþΔt I ði;jÞ
¼ F ði;jÞ t S3ði;jþ1Þ ;
(1d)
where F ði;jÞ ≡Fðdði;jÞ Þ is an attenuation factor that depends on the distance dði;jÞ between the node (i,j) in the absorbing layer and the computational domain limit (expressed in terms of number of nodes). De Cogan et al. proposed to consider the following attenuation factor: " 2 # −dði;jÞ 1−exp ; (2) B
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
4563
between 0 and 1, where B is a decay constant that determines the rate of transition. The choice of an attenuation factor that gradually increases along the absorbing layer ensures to minimize unwanted reflections within the layer. However, no application is given in Ref. [18]. Proposed formulation. It is important to remind that the basis of Bérenger's PML implementation consists in including an artificial absorbing medium just before the limit of the computational domain, acting on the normal components of the acoustic observables (i.e. the acoustic pressure and the particle velocity). On the other hand, in the method suggested by De Cogan's et al., all components within the absorbing layer are attenuated. The proposed formulation combines both approaches, by attenuating the incident pulses within the layer propagating toward the computational domain limit only (i.e. “one-way” attenuation). This corresponds to Eq. (1a) for a layer located at the west side of the domain, Eq. (1b) for the east side limit, Eq. (1c) for the north of the computational domain and Eq. (1d) in the south limit case. This method slightly differs from Bérenger's approach since, in this last one, the sound field is attenuated in both directions of each cartesian axis at once (e.g. +x and −x simultaneously), while, in the proposed approach, only the field components toward a given direction (+x or −x) are attenuated. Qualitatively (it will be presented later), the modification of one of the connection laws only, within the absorbing layer effectively leads to a better minimization of unwanted reflections, instead of using four or two connection laws simultaneously. This effect can be explained by the fact that the pressure field continuity condition is almost verified in the case of a “one-way” attenuation, while an attenuation in all directions simultaneously breaks more suddenly this condition. The effect can also be seen as an unexpected impedance jump, according to the TLM concept. In the proposed approach, the form of the attenuation factor F ði;jÞ is equivalent to the one suggested by De Cogan et al. (see Eq. (2)), which offers the advantage to swiftly attenuate the incident pulses at the beginning of the layer (with F ði;jÞ tending to 1 when dði;jÞ is equal to the layer depth eAL ), and to gently vanish the pulses at the vicinity of the computational domain limit (F ði;jÞ tending to 0 when dði;jÞ tends to 0). However, it must be reminded that a null reflection coefficient on the limit can generate instabilities. Finally, the attenuation factor F is chosen in such a way that it varies from Fðdði1 ;j1 Þ ¼ 0Þ ¼ ϵ;
(3)
for nodes of discrete coordinates ði1 ; j1 Þ located at the computational domain limit (Fig. 1), where the term ϵ tends to zero, to FðdðiN ;jN Þ ¼ eAL Þ ¼ 1;
(4)
for the nodes of discrete coordinates ðiN ; jN Þ on the interface between the propagation medium and the absorbing layer, N being the number of nodes inside the absorbing layer. In the last equation, eAL is the layer thickness, which can be expressed in terms of number of nodes as eAL ¼
λN λAL ; Δl
(5)
where the parameter NλAL represents the equivalent depth of the absorbing layer in wavelength unit λ (e.g. NλAL ¼ 5 corresponds to a layer thickness equivalent to 5 wavelengths). In order to satisfy conditions (3) and (4), using a similar function as the one suggested by De Cogan et al., the attenuation factor can be chosen as " Fðdði;jÞ Þ ¼ ð1 þ ϵÞ−exp
2
−dði;jÞ B
# ;
(6)
with 0 o ϵ≤1 and where the decay constant is given by B¼−
e2AL : ln ϵ
(7)
This expression of the attenuation factor is equal to 1 for dði;jÞ ¼ eAL , which ensures the pressure field continuity between both media, and tends to ϵ when dði;jÞ tends to zero. Fig. 2 presents the attenuation factor for several values of ϵ, i.e. different values of the decay constant B. Results are given as a function of the distance d from the computational domain limit, expressed in terms of node number. It can be observed that the slope becomes more gentle at the beginning of the absorbing layer when reducing the value of ϵ (reducing the decay constant B reciprocally), which should limit unwanted reflections thanks to a better respect of the pressure field continuity condition during the sound propagation inside the absorbing layer. Although other damping functions could be formulated in order to fulfill the conditions (3) and (4), the proposed formulation (6) already shows an appropriate behavior both on its bounds and during the slope. However, the study of other damping functions could also be investigated. This formula involves the lower frequency of the simulation through the layer thickness eAL given at Eq. (5), which represents a critical parameter as it strongly affects the absorbing layer efficiency.
4564
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
Fig. 2. Attenuation factor F ði;jÞ as a function of the distance dði;jÞ from the limit of the computational domain, for several values of ϵ, in a layer depth equivalent to 50 nodes.
3. Numerical results 3.1. Principle and geometries 3.1.1. Reflection error The efficiency of the absorbing layer can be characterized by comparing the sound propagation in a computational domain surrounded by artificial absorbing boundary conditions (i.e. in the “virtual free-field” case), with the equivalent propagation in the free-field case. In practice, the acoustic pressure pðtÞ that is obtained when the propagation medium is delimited using absorbing conditions, is compared with the equivalent free-field acoustic pressure pff ðtÞ obtained with a larger computational domain (i.e. in order to eliminate unwanted reflections in the time scale of interest). Then, a reflection error relative to free-field is defined in dB, by [29] errorðx; yÞ ¼ 10 log10
∑Tt ¼ 0 jpff ðx; y; tÞ−pðx; y; tÞj2 ∑Tt ¼ 0 jpff ðx; y; tÞj2
;
(8)
where T is the total simulation duration. A mean error over the whole propagation medium can also be evaluated by averaging the reflection error on x- and y-directions. Ideally, no unwanted reflection occurs for a perfect absorbing layer, implying that pðtÞ ¼ pff ðtÞ, and a reflection error tending to −∞. In practice, since unwanted reflections are present in the absorbing layer, and boundaries do not produce sound energy, the reflection error is strictly negative. 3.1.2. Geometry Two geometries are considered in this paper in order to evaluate the absorbing layers efficiency. For the first geometry given in Fig. 3, absorbing layers of depth eAL are set downstream from the upper and from the right side propagation medium limits. The point source S is located at the bottom left side corner of the propagation medium of dimensions Dx Dy at the junction of two perfectly reflective boundaries. This geometry allows the simulation of a “freefield” domain of dimensions ð2 Dx Þ ð2 Dy Þ that would be completely surrounded by absorbing layers, the sound source being in the center of the computational domain. The two excitation frequencies are considered for the simulations presented further. For simulations performed at 100 Hz (Sections 3.2, 3.3.1 and 3.4), the dimensions of the propagation medium are Dx ¼ Dy ¼ 10 m. In order to estimate the reflection error, a corresponding reference free-field propagation medium is also defined, with dimensions Dx ¼ Dy ¼ 85 m such as no reflection can return within the domain of interest until the simulations end, for a simulation duration of t ¼0.25 s. For simulations carried out at 1 kHz (Sections 3.4 and 3.5), the dimensions of the propagation medium and the reference free-field domain are Dx ¼ Dy ¼ 5 m and Dx ¼ Dy ¼ 75 m respectively, for a simulation duration of t¼0.21 s. The second geometry, given in Fig. 4, is useful for evaluating the efficiency of absorbing conditions as a function of the incidence angle of sound waves (Section 3.3.2). An absorbing layer is set downstream from the upper limit regarding the propagation medium of size Dx Dy , with Dx ¼ 700 m and Dy ¼ 365 m. The source is placed on the left boundary that
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
4565
Fig. 3. Computational domain used to evaluate the absorbing layers efficiency. The two-dimensional computational domain is made up of a propagation medium with dimensions Dx Dy extended by absorbing layers of depth eAL at the upper and right sides of the domain of interest. The sound source is located at the lower left corner at the junction of two perfectly reflective boundaries, in order to simulate a two-dimensional domain of size ð2 Dx Þ ð2 Dy Þ entirely surrounded by absorbing layers.
Fig. 4. Computational domain designed in order to estimate the efficiency of artificial absorbing boundary condition as a function of the incidence angle on the absorbing boundary. The source is placed on the left boundary, which is perfectly reflecting, in order to enlarge the propagation domain according to the x-axis, from Dx Dy to ð2 Dx Þ Dy . Receivers are distributed on each node from the source to the right side, along the x-axis, and are associated with incident angles θ on the upper interface between both media.
artificially increases the dimension of the computational domain along the x-axis to ð2 Dx Þ Dy . The reference “free-field” domain is enlarged along the y-axis to Dy ¼ 710 m. All simulations are performed on a duration t ¼2 s with a sound source located at a distance Δ′ ¼ 355 m from the lower limit of the computational domain (i.e. at a distance Δ ¼ 10 m from the upper limit). Due to the substantial computational requirement for the simulations with this geometry, all numerical simulations are carried out with a Gaussian pulse source at 100 Hz, as detailed in [9, Eq. (37)], in order to ensure relevant results over a large range of incident angles θ.
3.2. Effect of the layer depth and of the parameter ϵ A first numerical study is carried out in order to evaluate the effect of the parameter ϵ and of the layer depth N λAL on the absorbing condition efficiency. The first geometry (Fig. 3) is used with a Gaussian sound source at 100 Hz. Fig. 5 gives the reflection error as a function of the absorbing layer depth as NλAL ranges from 1 to 10, for several values of the parameter ϵ. As expected, a substantial improvement of the layer efficiency is observed while reducing the parameter ϵ and while enlarging the artificial absorbing layer, which is due to a better respect of the pressure field continuity condition as the sound waves progress inside the absorbing layer. However, it appears that the efficiency improvement of the absorbing layer is less and less significant as ϵ decreases. Finally, the choice of ϵ ¼ 10−5 and NλAL ¼ 5 seems to be a good compromise between
4566
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
Fig. 5. Effect of the parameter ϵ on the mean reflection as a function of the absorbing layer thickness N λAL .
absorbing layer efficiency and computational burden. Except when explicitly mentioned, all numerical simulations presented in the following are carried out with this set of parameters. 3.3. Comparison between absorbing conditions This section presents a comparison of the absorbing conditions presented in Section 2.2, with the proposed formulation. 3.3.1. Comparison with De Cogan's formulation The gradually matched layer implementation proposed by De Cogan et al. (Section 2.2.2) is compared with the proposed implementation of a “one-way” attenuation in the absorbing layer using the first geometry (Fig. 3), both with the same attenuation factor F ði;jÞ . It must be noted that preliminary numerical simulations have shown that the parameter ϵ has no effect when considering De Cogan's absorbing layer implementation. Fig. 6 presents the acoustic pressure pðtÞ obtained at a receiver M located at a distance of 1 m both from the right and the lower sides of the propagation domain. Annotations 1 and 2 on Fig. 6 stand for the direct pulse and the first unwanted pulse reflections detected using De Cogan's implementation (2) and the proposed (2′) formulations respectively. Fig. 6 shows the significant reduction of the unwanted reflected pulse using the “one-way” attenuation. The delay between the unwanted reflected pulses 2 and 2′ comes from the fact that these reflections occur inside the absorbing layer and not on the two geometrical limits of the absorbing layer (i.e on the limit of the computational domain and on the interface between the two propagation media). It is also important to observe that the present absorbing layer implementation induces a spreading of the spurious pulse that is related to a reduction of the unwanted reflection magnitude. Both absorbing condition formulations are then compared by computing the mean reflection error for several absorbing layer thicknesses N λAL ranged from 1 to 10 (Fig. 7). Here again, the present formulation leads to a reduction of the unwanted reflection magnitude in comparison with De Cogan's approach whatever the layers thickness is. The mean reflection error obtained with De Cogan's formulation remains constant (around −15 dB) independent of the layer depth, due to the generation of unwanted reflections just on the interface between the propagation medium and absorbing layers, and not within the absorbing layers as with the proposed formulation. It must be reminded that the reflection error represents the ratio of the reflected field to the incident field at the same location in the computational domain. Thus Fig. 7 means for instance that, for De Cogan's approach, the reflected field is lower than the incident field by about −15 dB so that the cumulated sound field at any observation point in the computational domain will be mainly supplied by the incident field. The improvement of the absorption layer efficiency using the proposed formulation is not significant in this situation. However, when considering a propagation medium that generates multiple diffuse reflections on the limit of the computational domain and on absorbing boundaries, the cumulative energy of the unwanted reflections becomes non-negligible in comparison with the incident energy. This is particularly true for the numerical modeling of sound propagation in urban areas, when considering for instance a street for which the open-top is “closed” by an absorbing layer (see Section 4). 3.3.2. Comparison with absorbing boundaries The two absorbing boundaries based on Taylor's series expansion of the sound pressure field on the limit of the computational domain (Section 2.2.1) and on purely real impedance conditions (Section 2.2.1) are now compared with the
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
4567
Fig. 6. Temporal pressure signals pðtÞ recorded at a node located at a distance of 1 m both from the lower and from the right boundaries of the propagation domain (Fig. 3), using De Cogan's implementation of the absorbing layer and the proposed formulation.
Fig. 7. Mean reflection error as a function of the layer thickness N λAL for De Cogan's implementation and the proposed formulation.
proposed formulation. Simulations are performed using a Gaussian sound pulse at 100 Hz in the second geometry (Fig. 4), in order to compare all artificial absorbing conditions for several incidence angles θ. Fig. 8 shows the unwanted temporal signals obtained at a normal incidence (i.e. θ ¼ 01) for the different approaches: with the absorbing boundaries defined by Taylor's series expansion pffiffiffi of the acoustic pressure on the limit of the computational domain, with a purely real impedance condition (for Z ¼ 2Z 0 and Z ¼ Z 0 ), with De Cogan's implementation and with the proposed formulation. Taylor's expansion based condition gives the worst results and produces large numerical instabilities as underlined in Section 2.2.1. However, it must be noted that the technique that consists in applying a null reflection coefficient on the limit after the passing of the wavefront has not been performed, since this approach is not relevant when considering permanent sound field or long duration simulations. Due to these instabilities, the reflection error is very important; this absorbing condition is no longer considered in the following. Fig. 9 shows the mean reflection error as a function of the incident angle θ, for a purely real impedance condition and for absorbing layers implemented with De Cogan's approach and with the proposed formulation. The purely real impedance condition gives quite interesting results. However, as explained by Hofmann and Heutschi, the impedance Z must be pffiffiffi adapted according to the angle θ (a better absorption is obtained for θ ¼ 01 with Z ¼ Z 0 than with Z ¼ 2Z 0 , while opposite results are observed for θ ¼ 451), what consequently makes difficult the use of this method, especially when considering multiple sound sources. As alreadypobserved in Fig. 8, De Cogan's approach gives similar results than the purely real ffiffiffi impedance condition defined by Z ¼ 2Z 0 . Lastly, although the proposed formulation leads to similar results than the purely
4568
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
Fig. 8. Temporal pressure signals pðtÞ recorded at normal incidence (θ ¼ 01) for absorbing boundaries defined pffiffiffi by a Taylor's series expansion of the acoustic pressure on the computational domain limit and by a purely real impedance condition (Z ¼ Z 0 and Z ¼ 2Z 0 ), compared with the proposed formulation.
Fig. 9. Mean reflection error as a function of the incidence angle θ, for absorbing boundaries defined by pffiffiffia Taylor's series expansion of the acoustic pressure on the computational domain limit and by a purely real impedance condition (Z ¼ Z 0 and Z ¼ 2Z 0 ), and for absorbing layers using De Cogan's implementation and the proposed formulation.
pffiffiffi real impedance condition with Z ¼ 2Z 0 for θ o 151, it leads to better results when the whole incidence angles range is considered. Moreover, as it can be observed in Fig. 9, the mean reflection error remains quite constant at about −40 dB when θ is below 651. However, for larger angles, the mean reflection error increases for all the tested absorbing conditions. It must be pointed out that, unlike the PML approach, the present absorbing layer implementation is not relevant for attenuating the sound field components parallel to the computational domain limit, which explains the efficiency declines at grazing incidence. Then, the proposed absorbing boundary condition is more satisfying than all other virtual boundary formulations. However, it must be pointed out that the attenuation is much lower than the one obtained with current PML implementation within other numerical models like FDTD.
3.4. Stability Simulations have been carried out using the first geometry (Fig. 3) in order to verify the stability of the proposed absorbing layer implementation for long duration simulations. Numerical simulations have been performed using a sine
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
4569
Fig. 10. Reflection error spectrum for two absorbing layer depths, both with N λAL ¼ 5, based on the upper and lower values of the frequency range (f max ¼ 1 kHz and f min ¼ 100 Hz respectively).
wave source at both frequencies 100 Hz and 1 kHz, for a duration of 30 s and 5 s, respectively. It has been observed that the temporal signals remain very stable in the whole propagation medium for all the duration of the simulations. 3.5. Variation of efficiency with frequency The influence of the frequency on the absorbing layer efficiency was also studied in order to determine the ability of the present formulation for broadband simulations. The reflection error spectrum, relative to free-field, is computed using the following relation: errorðωÞ ¼ 10 log10
PSD½pff ðx; y; tÞ−pðx; y; tÞ ; PSD½pff ðx; y; tÞ
(9)
where PSD is the power spectral density function. The average error is also computed on the frequency range of interest. Numerical simulations presented further have been carried out using the first geometry, with a Gaussian sound pulse excitation of central frequency 1 kHz. Fig. 10 shows the reflection error spectrum computed on the frequency range of interest [100,1000] Hz, at the receiver N located at a distance of 1 m from the lower limit of the propagation domain, and at the last node before the absorbing layer along the x-axis. Simulations have been performed for two layer depths (i.e. for two reference frequencies), both with NλAL ¼ 5. The first one is sized by adjusting the thickness of the layer from the upper frequency f max ¼ 1 kHz (i.e. the minimal wavelength) of the frequency range, while the second case is based on the lower frequency f min ¼ 100 Hz (i.e. the maximal wavelength). As expected, the improvement of the efficiency of the absorbing layer is significant when adjusting the absorbing region thickness with the maximal wavelength, leading to very good results for the whole frequency range. The average reflection error spectrum over the frequency range [100,1000] Hz is also computed for several layer depths NλAL ranging from 1 to 10 (Fig. 11). It can be observed that the efficiency of the absorbing layer is almost constant, whatever the layer depth is, here again when absorbing region thickness is based on the maximal wavelength. As expected, the layer efficiency is equivalent for N λAL ¼ 10 with f ¼1 kHz than for NλAL ¼ 1 with f¼ 100 Hz, since the layer depths are equal in both cases. 4. Urban acoustics application A realistic case of sound propagation within a street section is considered to illustrate the efficiency of the proposed formulation. The same geometry as depicted in Fig. 4 is used. The street has width and height of Dx ¼ 6 m and Dy ¼ 9 m, respectively, that correspond to a street with an aspect ratio of 0.7. A Gaussian sound source (100 Hz) is located at Δ′ ¼ 1 m above the ground (i.e. x¼ 0). Receivers are located at each node along a parallel axis to the ground, passing through the source. As previously suggested in Section 3.3.1, this kind of geometry generates multiple reflections on boundaries; the absorbing boundaries are then particularly involved. Equivalent sound levels at receivers are calculated on a duration T ¼1 s according to the relation " # 1 ∑Tt ¼ 0 jpðx; y; tÞj2 Leq;T ¼ 10 log10 ; (10) T p20
4570
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
Fig. 11. Average reflection error on the whole frequency range [100,1000] Hz as a function of the absorbing layer thickness N λAL . The layer depth is based either on the upper frequency (f max ¼ 1 kHz) or the lower frequency (f min ¼ 100 Hz) of the frequency range.
Fig. 12. Equivalent sound pressure levels at every node along the street transversal axis located at Δ′ ¼ 1 m high above ground. The façades and the ground are perfectly reflecting. The open-top is modeled by a real impedance condition (Z ¼ Z 0 ) and by absorbing layers following De Cogan's approach as well as the proposed formulation.
with p0 the reference sound pressure. Numerical simulations have been performed using perfectly reflecting boundary conditions on ground and façades, and several absorbing conditions for the open-top of the street. The “free-field” open-top configuration is obtained by enlarging the propagation domain along the y-dimension (Dy ¼ 200 m). Fig. 12 shows the efficiency of the proposed formulation in comparison with other approaches. Results confirm that wavefronts reach the virtual limit (i.e. the absorbing boundary or the layer) with multiple angles of incidence, since the equivalent sound pressure level along the street width is around 3 dB higher with the purely real boundary condition (i.e. Z ¼ Z 0 ) than with the proposed absorbing layer implementation using NλAL ¼ 5 (see Section 3.3.2 and Fig. 9). In addition, although the proposed formulation of absorbing layers gives satisfying results with a thickness of N λAL ¼ 5, with an error lower than 1 dB, the enlarged absorbing layer (i.e. N λAL ¼ 10) leads to sound pressure levels which are much more similar to the reference case. 5. Conclusion An absorbing layer formulation for TLM modeling of sound propagation is presented in this paper. The proposed formulation is based on De Cogan's suggestion of absorbing layer implementation, for which the matched connection laws
G. Guillaume, J. Picaut / Journal of Sound and Vibration 332 (2013) 4560–4571
4571
are restricted to the incident pulse propagating toward the absorbing layer only. The proposed approach improves results significantly, compared with other absorbing layer and absorbing boundary models. Numerical simulations have been carried out in order to estimate the effect of the attenuation factor, the layer depth and the incidence angle on the absorbing layer efficiency. It is shown that a decay constant computed with ϵ ¼ 10−5 with a sufficiently wide absorbing layer (e.g. a depth equivalent to 5 or more wavelengths) guarantees a significant reduction of unwanted reflections on the limit of the computational domain. Moreover, the present formulation ensures a good efficiency for oblique angles except for grazing ones, because the method does not stand for attenuating the sound field components at grazing incidence. As expected, results also show that, for a given frequency range of interest, the depth of the absorbing layer must be designed considering the lower frequency limit in order to ensure a significant attenuation of unwanted reflections over the whole corresponding frequency range. However, the present formulation does not verify a priori Bérenger's PML equations. The efficiency of the proposed formulation is then much lower than the one expected with “real” PML as used in FDTD methods. Moreover, the proposed implementation results in a significantly expensive computational algorithm, in terms of CPU-time as well as in terms of memory storage. However, the formulation leads to acceptable efficiency for broadband frequencies and oblique incidences, which is of particular interest for long range outdoor sound propagation, like in urban areas, even if its application for shielded areas for example has to be avoided. In addition, although the present approach is mainly dedicated to outdoor sound propagation applications, the proposed formulation of absorbing layers can be applied in other domains of acoustics. References [1] D. Heimann, Three-dimensional linearised Euler model simulations of sound propagation in idealised urban situations with wind effects, Applied Acoustics 68 (2) (2007) 217–237. [2] T. Van Renterghem, D. Bottledooren, Numerical evaluation of sound propagating over green roofs, Journal of Sound and Vibration 317 (3–5) (2008) 781–799. [3] P.B. Johns, R.L. Beurle, Numerical solution of two dimensional scattering problems using a transmission line matrix, Proceedings of the IEE 118 (9) (1971) 1203–1208. [4] A.H.M. Saleh, P. Blanchfield, Analysis of acoustic radiation patterns of array transducers using the TLM method, International Journal of Numerical Modelling: Electronic Networks, Devices and Fields 3 (1990) 39–56. [5] G. Guillaume, J. Picaut, G. Dutilleux, G. Gauvreau, Use of the transmission line matrix method for the sound propagation modelling in open-space, Proceedings of the 13th International Symposium on Long Range Sound Propagation, Proceedings on CD-ROM, Lyon, France, 2008, pp. 417–426. [6] G. Guillaume, J. Picaut, G. Dutilleux, Use of the transmission line matrix method for the sound propagation modelling in urban area, Proceedings of the ASA-EAA-CFA Joint congress, Acoustics'08, Paris, France, 2008, pp. 1901–1906. [7] U. Kristiansen, K. Jezzine, TLM model for sound propagation above ground, Proceedings of the 29th International Congress and Exhibition on Noise Control Engineering (Inter-Noise 2000), Nice, France, 2000. [8] J. Hofmann, K. Heutschi, Simulation of outdoor sound propagation with a transmission line matrix method, Applied Acoustics 68 (2) (2007) 158–172. [9] G. Guillaume, J. Picaut, G. Dutilleux, B. Gauvreau, Time-domain impedance formulation for transmission line matrix modelling of outdoor sound propagation, Journal of Sound and Vibration 330 (26) (2011) 6467–6481. [10] T. Tsuchiya, Numerical simulation of sound wave propagation with sound absorption in time domain, Proceedings of the 13th International Congress on Sound and Vibration, Vienna, Austria, 2006. [11] G. Dutilleux, Applicability of TLM to wind turbine noise prediction, Proceedings of the Second International Congress on Wind Turbine Noise, Lyon, France, 2007. [12] J.-P. Bérenger, A perfectly matched layer for the absorption of electromagnetic waves, Journal of Computational Physics 114 (2) (1994) 110–117. [13] J.L. Dubard, D. Pompéi, Simulation of Bérenger's perfectly matched layer with a modified TLM node, IEE Proceedings—Microwaves, Antennas and Propagation (Part H) 144 (3), 1997, pp. 205–207. [14] J.L. Dubard, D. Pompéi, Optimization of the PML efficiency in 3-D TLM method, IEEE Transactions on Microwave Theory and Techniques 48 (7) (2000) 1081–1088. [15] J. Diaz, P. Joly, A time domain analysis of PML models in acoustics, Computer Methods in Applied Mechanics and Engineering 195 (2006) 3820–3853. [16] S. El-Masri, X. Pelorson, P. Saguet, P. Badin, Vocal tract acoustics using the transmission line matrix (TLM) method, Proceedings of the International Conference on Spoken and Language Processing, vol. 2, Philadelphia, USA, 1996, pp. 953–956. [17] J.A. Portí, J.A. Morente, TLM method and acoustics, International Journal of Numerical Modelling: Electronic Networks, Devices and Fields 14 (2001) 171–183. [18] D. De Cogan, W.J. O'Connor, S. Pulko, Transmission Line Matrix in Computational Mechanics, CRC Press, London, England, 2005 ISBN 978-0-415-32717-6. [19] P. Saguet, TLM method for the three dimensional analysis of microwave and MM-wave structures, Proceedings of the International Workshop of the German IEEE M7T/AP Chapter on CAD Oriented Numerical Techniques for the Analysis of Microwave and MM-Wave Transmission-Line Discontinuities and Junctions, Stuttgart, Germany, 1991. [20] D.T. Murphy, J. Mullen, Digital waveguide mesh modelling of room acoustics: improved anechoïc boundaries, Proceedings of the Fifth International Conference on Digital Audio Effects, Hamburg, Germany, 2002. [21] J. Nielsen, Spurious modes of the TLM-condensed node formulation, IEEE Microwave and Guided Wave Letters 1 (8) (1991) 201–203. [22] P. Saguet, Analyse des milieux guidés—la méthode MTLM (Guided medium analysis: MTLM method), PhD Thesis, Institut National Polytechnique de Grenoble, France, 1985. [23] R.L. Higdon, Absorbing boundary conditions for difference approximations of the multi-dimensional wave equation, Mathematics of Computation 47 (176) (1986) 437–459. [24] C. Eswarappa, W.J.R. Hoefer, One-way equation absorbing boundary conditions for 3-D TLM analysis of planar and quasi-planar structures, IEEE Transactions on Microwave Theory and Techniques 42 (9) (1994) 1669–1677. [25] Z. Chen, M. Ney, W.J.R. Hoefer, Absorbing and connecting boundary conditions for the TLM method, IEEE Transactions On Microwave Theory and Techniques 41 (11) (1993) 2016–2024. [26] F.Q. Hu, On absorbing boundary conditions for linearized Euler equations by a perfectly matched layer, Journal of Computational Physics 129 (1) (1996) 201–219. [27] Q.-H. Liu, J. Tao, The perfectly matched layer for acoustic waves in absorptive media, Journal of the Acoustical Society of America 102 (4) (1997) 2072–2082. [28] G.-X. Fan, Q.H. Liu, A well-posed PML absorbing boundary condition for lossy media, IEEE Antennas and Propagation Society International Symposium 3 (2001) 2–5. [29] Y.-L. Sheu, P.-C. Li, Simulations of photoacoustic wave propagation using a finite-difference time-domain method with Berenger's perfectly matched layers, Journal of the Acoustical Society of America 124 (6) (2008) 3471–3480.