Ray-tracing semiclassical low frequency acoustic modeling with local and extended reaction boundaries

TL;DR

This paper validates the ray-tracing semiclassical (RTS) method for low frequency acoustic modeling with complex boundary conditions, demonstrating its effectiveness in realistic room acoustics scenarios involving frequency-dependent materials.

Contribution

The study extends the RTS method to handle local and extended reaction boundaries, including porous layers, for more accurate low frequency acoustic modeling in realistic environments.

Findings

RTS accurately reproduces Green's function for complex boundaries.

RTS results closely match finite element method in pressure and decay curves.

RTS effectively models low frequency sound fields with diverse boundary conditions.

Abstract

The recently introduced acoustic ray-tracing semiclassical (RTS) method is validated for a set of practically relevant boundary conditions. RTS is a frequency domain geometrical method which directly reproduces the acoustic Green's function. As previously demonstrated for a rectangular room and weakly absorbing boundaries with a real and frequency-independent impedance, RTS is capable of modeling also the lowest modes of such a room, which makes it a useful method for low frequency sound field modeling in enclosures. In practice, rooms are furnished with diverse types of materials and acoustic elements, resulting in a frequency-dependent, phase-modifying absorption/reflection. In a realistic setting, we test the RTS method with two additional boundary conditions: a local reaction boundary simulating a resonating membrane absorber and an extended reaction boundary representing a porous layer backed by a rigid boundary described within the Delany-Bazley-Miki model, as well as a combination thereof. The RTS-modeled spatially dependent pressure response and octave band decay curves with the corresponding reverberation times are compared to those obtained by the finite element method.

Figures

Tables

Table 1. Table 1: The octave-band-averaged diffuse absorption coefficients for the three BCs.

	${\overline{\alpha }}_{\text{diff}}$
BC	63 Hz	125 Hz
real impedance	$0.12$	$0.12$
resonance	$0.27$	$0.15$
extended reaction	$0.13$	$0.24$

Equations

$$p(\mathbf{r},t)=A\exp [j(\mathbf{k}\cdot \mathbf{r}-2\pi ft)],$$

$$B_i=r{B}_{i-1},$$

$$r(\theta )=\frac{z\cos \theta -z_0}{z\cos \theta +z_0},$$

$${\alpha }_{\text{diff}}=2{\int }_0^{\pi /2}\alpha (\theta )\cos \theta \sin \theta \mathrm{d}\theta ,$$

$$z_m(f)=r_m-j\left(2\pi f{m}_S-{\rho }_0{c}_0\cot k{d}_m\right)$$

$$z_m(f)=r_m-j\frac{2\pi m_S}{f}\left(f^2-f_m^2\right),$$

$$f_m=\frac{c_0}{2\pi }\sqrt{\frac{{\rho }_0}{m_S{d}_m}}$$

$$z_p=z_0{g}_z\left(\frac{f}{\sigma }\right),k_p=k{g}_k\left(\frac{f}{\sigma }\right),$$

$$g_z(x)$$

$$g_k(x)$$

$$z_e(f,\theta )$$

$$-\left({\nabla }^2+k^2\right)G(\mathbf{r},{\mathbf{r}}_0,k)=\delta (\mathbf{r}-{\mathbf{r}}_0),$$

$$-\left({\nabla }^2+k^2\right)G_p(\mathbf{r},{\mathbf{r}}_0,k)=({\dot{\varphi }}_m{)}_0\delta (\mathbf{r}-{\mathbf{r}}_0).$$

$$({\dot{\varphi }}_m{)}_0=\sqrt{8\pi {\rho }_0{c}_{0}P}.$$

$$G_p(\mathbf{r},{\mathbf{r}}_0,k)=\sqrt{8\pi {\rho }_0{c}_{0}P}G(\mathbf{r},{\mathbf{r}}_0,k).$$

$$G(\mathbf{r},{\mathbf{r}}_0,k)=\frac{1}{\pi R^{2}N}\sum _l^L{B}_{l,i}{d}_l\exp \left[jk{d}_l\right],$$

Full text

I Introduction

Room acoustics simulations have become a valuable design tool for acousticians to predict room acoustic parameters, model room modal shapes, optimize the positioning of public address systems etc. Notwithstanding their increasingly versatile applicability, these tools should be still used carefully. In fact, it is not in all circumstances that relevant results are obtained and caution as well as thorough understanding of the underlying principles of each modeling technique is required from the user. According to their fundamental physical assumptions, the simulation methods used in acoustics are primarily divided into wave-based methods and geometrical methods.

Wave-based methods solve the complete wave equation for a given geometry and boundary conditions (BCs) and are therefore in principle able to model all properties of waves — interference, diffraction, and scattering. In practice these effects are crucial at low frequencies, where boundary element methods (BEM) kawai:2007 , finite element methods (FEM) easwaran:1995 , and finite difference methods (FDM) kowalczyk:2008 are most commonly used. Unfortunately, at higher frequencies wave-based methods become over demanding in terms of processing time and memory usage (desmet:2002, ; harari:1992, ) although efforts have been devoted to extend their computationally effective frequency range (e.g., mehra:2012 ).

On the other side, geometrical methods assume that sound propagates from the source in straight lines, implying that diffraction/scattering effects are not included by definition. Generally, the geometrical modeling techniques are limited to large rooms/sufficiently high frequencies, where the sound field is not significantly shaped by individual room modes. A special, rather restricted but practically widely used subgroup of these methods are energy methods, which consider solely the propagation of sound energy implying that also interference effects are disregarded. These methods are used in commercial room acoustics simulation engines (e.g., ODEON Naylor , CATT-Acoustic Dalenback ). The energy methods are complemented by the so-called phased geometrical methods in which phase information is included.

Various implementations of geometrical methods exist and their further differentiation is commonly performed on the basis of their stochastic nature (ray-tracing method krokstad:1968 ), the presence of virtual sources (image source method allen:1979 ) or other criteria. Each implementation has its advantages and disadvantages for specific types of problems; for recent reviews cf., e.g., refs. vorlander:2013 ; savioja:2015 ; pompei:2015 .

Generally, the inclusion of the phase information improves the quality of the results obtained at low frequencies. Suh and Nelson showed suh:1999 that the geometrically modeled impulse response agrees better with the experiment if the phase is included. Various aspects of phased image source implementations have been independently studied and compared with other modeling methods vorlander:1989 ; jeong:2008 ; dance:1995 ; aretz:2014 ; lam:2005 . We (Prislan et al.) prislan:2016 have been concentrating on reproducing individually resolved resonances well below the Schröder frequency schroder:1954 by implementing the ray-tracing semiclassical (RTS) method which reproduces the Green’s function, $G$, rigorously within the semiclassical (eikonal) approximation. Our study showed that the results obtained by RTS can correctly model the modal structure of a rectangular room including its lowest modes, which is decisive for using the method at low frequencies.

To achieve relevant results with room acoustics simulations it is also important to fully understand the implemented BC. In geometrical modeling the BC is described in terms of the boundary impedance which determines the reflection coefficient of the boundary as if it were an extensive surface. It is common to assume that the boundaries are locally reacting, in which case the surface impedance is independent of the angle of incidence of the sound wave rossing:2007 . For locally reacting boundaries, the importance of the BC in a phased geometrical method was demonstrated by Jeong jeong:2012 .

Experiments on material samples klein:1980 ; tamura:1995 showed that the assumption of local reaction is not accurate for a layer of porous material covering a rigid base, in which case the boundary impedance evidenced angle dependence. In such case the extended reaction assumption is preciser and improves the theoretical and experimental agreement of sound pressure level predictions franzoni:2003 . Moreover, the significance of introducing extended reaction boundaries into computer models has been demonstrated – using phased beam tracing, the influence on the steady state sound pressure level in frequency bands was shown by Hodgson and Wareing Hodgson:2008 and on several acoustic parameters (strength, reverberation time and rapid speech transmission index) by Yousefzadeh and Hodgson yousefzadeh:2012 . Recently, Yasuda et al. yasuda:2016 showed by wave-based modeling that the extended reaction alters the response in particular at low frequencies, which is also the frequency range of our interest.

In the previous study prislan:2016 we have shown that in a rectangular room with a rather clean BC the lowest room modes can be correctly modeled by the RTS modeling technique. The main focus of the study was firstly to adopt the semiclassical theory and establish a connection of this elaborated concept with acoustic modeling, and secondly to test the behavior of the RTS technique by comparing the modeled $G$ with its analytical solution. Weakly damped boundaries were assumed in the simulation, corresponding to a real frequency- and angle-independent surface impedance (local reaction boundary) identical for all bounding surfaces.

The aim of the present study is to validate the RTS method in a practically more relevant environment, i.e., by using more realistic BCs. Besides i) the BC used in the previous study, two additional BCs are implemented: ii) a frequency-dependent local reaction BC of a resonating membrane absorber and iii) an extended reaction BC of a porous boundary layer covering a rigid base described within the Delany-Bazley-Miki model miki:1990 . These particular BCs were chosen because they are well known in theoretical acoustics. On equal basis, alternative BC models can be implemented into RTS as long as the reflection coefficient can be formulated.

The results obtained by RTS are compared with FEM which is considered reliable at low frequencies. In a rectangular room with the same proportions as in the previous study prislan:2016 , the comparison is made for the three BCs applied homogeneously and on top of that, for their more realistic inhomogeneous distribution. The comparison includes the magnitude of the obtained pressure response and the octave band decay curves with the corresponding reverberation times. The pressure response was computed up to 300 Hz, since at higher frequencies the FEM method would not produce correct results in acceptable time limits for the models under study.

The article is organized as follows. In Sec. II the BCs are introduced and in Sec. III the most important aspects of the simulation methods are presented. In Sec. IV, the studied cases are defined and the results are presented, compared and discussed. Main conclusions are drawn in Sec. V.

II Boundary conditions

A Specular reflection coefficient

As we are concerned with the behavior of the method at lower frequencies, only specular reflections are considered. Diffuse (stochastic) reflections are not included although they are known to be significant in geometrical modeling techniques vorlander:1995 ; savioja:2015 . A link between the diffuse sound field and its geometrical interpretation was described in ref. prislan:2014 , together with experimental tests that indicated that the complexity of the room boundaries (which produces diffuse reflections) is relevant only at higher frequencies.

In geometrical modeling techniques, the BC is implemented via the specular reflection coefficient $r$ kuttruff:2007 , which is the ratio $B/A$ of complex sound pressure amplitudes of the reflected wave, $p^{\prime}({\bf r},t)=B\exp[j(\mathbf{k}^{\prime}\cdot\mathbf{r}-2\pi f\,t)]$, and the incident wave,

[TABLE]

where $|{\bf k}|=|{\bf k}^{\prime}|\equiv k=2\pi f/c_{0}$ is the wavenumber, $f$ is the frequency, $c_{0}$ the speed of sound in air and $j$ the imaginary unit. Generally, $r$ is a complex quantity inducing both a magnitude as well as phase change of the reflected sound wave. As each sound ray is multiply reflected, its amplitude $B_{i}$ after the $i$-th reflection is calculated recursively as

[TABLE]

with $B_{0}=1$.

The reflection coefficient is completely determined by an effective impedance $z$ of the boundary kuttruff:2007 and the angle of incidence $\theta$ of the ray with respect to the boundary normal,

[TABLE]

where $z_{0}=\rho_{0}c_{0}$ is the air impedance and $\rho_{0}$ the density of air ($\rho_{0}=1.204$ kg/m3 and $c_{0}=343$ m/s is considered throughout this study).

Under the local reaction assumption the BC includes only the normal velocity component fahy:2001 and the effective impedance is just the specific impedance of the bounding material, which is angle-independent. In contrast, in the case of an extended reaction boundary, its effective impedance depends on the frequency and the angle of incidence, as discussed in Sec. D.

The reflection coefficient is directly related to the (energy) absorption coefficient $\alpha(\theta)=1-|r(\theta)|^{2}$. The half-space solid angle average of $\alpha(\theta)$ is the diffuse sound field absorption coefficient $\alpha_{\rm diff}$,

[TABLE]

which is known as the Paris formula mechel:2008 . For common materials only $\alpha_{\textrm{diff}}$ is standardly specified (e.g., measured according to the standard ISO354 ISO354 ), which opens the discussion jeong:2008 of how to determine the reflection coefficient $r$ required in phased geometrical modeling.

The computational efficiency of the RTS method is not compromised if $z$ depends on frequency and angle of incidence. In fact, in this study we put this under a direct test and introduce just these dependencies, arising naturally in the implementation of the three different BC described in detail in the following.

B Real impedance boundary

The BC previously tested with RTS prislan:2016 , i.e., weakly damped boundaries with a real impedance $z_{n}=60\,z_{0}$, is an example of local reaction and is used also in this study for comparison. The corresponding reflection and absorption coefficients follow from Eqs. (3) and (4) and are frequency-independent as also shown in Fig. 1. This BC is not particularly realistic, although rather low and frequency-independent absorption coefficients are in practice measured for massive structures with smooth finishing (e.g., painted concrete walls).

C Resonance absorber

This BC models a resonating sound absorber such as of a membrane type, characterized by its complex impedance cox:2009

[TABLE]

consisting of specific resistance $r_{m}$, specific inertance $-j2\pi fm_{S}$ of the membrane mass per area $m_{S}$ and specific stiffness $j\rho_{0}c_{0}\cot{kd_{m}}$ of the air-filled cavity with depth $d_{m}$. For low frequencies ($kd_{m}\ll 1$), $\cot kd_{m}\approx 1/kd_{m}$ and the impedance Eq. (5) is simplified to

[TABLE]

where

[TABLE]

is the resonance frequency of the absorber.

The parameters used in the simulation are $m_{S}=10$ kg/m${}^{2},\,d_{m}=0.14$ m and $r_{m}=15\,z_{0}$ resulting in the resonance frequency of $f_{m}=50.6$ Hz. In Fig. 1, the corresponding deep in the reflection coefficient and peak in the absorption coefficient can be seen at the resonance frequency.

D Extended reaction boundary

A layered structure of the boundary leads to the extended reaction BC. In our case we use a $d_{p}=5$ cm thick layer of porous material covering a perfectly rigid wall.

An empirical model of rigid-frame porous materials was introduced by Delany and Bazley in 1970s delany:1970 . The model was later modified and generalized by Miki miki:1990 and is today known as the Delany-Bazley-Miki model. In our implementation, we follow its recent description by Yasuda et al. yasuda:2015 .

In this model, the specific impedance of the bulk porous material $z_{p}$ and the wavenumber $k_{p}$ (dispersion relation) are given by

[TABLE]

where $\sigma$ [N s/m4] is the linear specific airflow resistivity and the functions $g_{z}(x)$ and $g_{k}(x)$ are

[TABLE]

with the constants $a_{z}=0.07,b_{z}=-0.107,q_{z}=0.632,a_{k}=0.109,b_{k}=-0.160,q_{k}=0.618$. Note that the signs of $b_{z}$ and $b_{k}$ are opposite to those of ref. yasuda:2015 , in accord with our definition of the plane-wave time factor, Eq. (1).

The effective impedance of such an extended reaction boundary, formed by a rigid wall covered with a layer of the described porous material with thickness $d_{p}$, can be shown to be yasuda:2015

[TABLE]

and depends on the frequency $f$ as well as on the angle of incidence $\theta$. To complete the picture, $\sqrt{g_{k}^{2}-\sin^{2}{\theta}}/g_{k}=\cos\theta^{\prime}$, where $\theta^{\prime}$ is the (complex) angle of refraction in the porous layer.

The value $\sigma=500$ Ns/m4 was used in the simulation and the corresponding reflection coefficient Eq. (3) and diffuse field absorption coefficient Eq. (4) are plotted in Fig. 1. It can be seen that in the considered frequency range the reflection coefficient monotonically decreases while the absorption coefficient monotonically increases with frequency.

III Methods

A Pressure response function

The Green’s function $G({\bf r},{\bf r}_{0},k)$ (the response to a dimensionless point source of strength 1 located at ${\bf r}_{0}$) of the Helmholtz equation prislan:2016 satisfies

[TABLE]

while the physical sound pressure response function $G_{p}({\bf r},{\bf r}_{0},k)$, i.e., the pressure response to a pointlike sound source of strength $(\dot{\phi}_{m})_{0}$, satisfies

[TABLE]

The amplitude $(\dot{\phi}_{m})_{0}$ of the rate of change of the radiative mass flux, i.e., the radiative mass acceleration generated by the sound source, is directly connected with the acoustic power $P$ emitted by the source,

[TABLE]

This relates the Green’s function $G$, which is furnished by RTS, with the sound pressure response function $G_{p}$:

[TABLE]

B RTS method

The RTS method resides on the semiclassical approximation of quantum mechanics and in our previous publication prislan:2016 we have demonstrated that this concept can be transferred also to acoustics. The resulting acoustic RTS method is a frequency domain method, the main advantage of which is that it yields directly the Helmholtz equation Green’s function $G$ of the room.

From the technical perspective it is a ray-tracing method, where sound rays are propagated from the source ${\bf r}_{0}$ in random directions, upon hitting a boundary they are reflected according to reflection rules, and are in the end summed up in a point (region) of interest. In our case, the rays are emitted in random directions with isotropic probability distribution. At the boundary they are specularly reflected and the pressure amplitude factor changes as defined by Eq. (2), with the reflection coefficient $r$ corresponding to the BC.

The segments of a multiply reflected ray constitute a trajectory. If a trajectory crosses the observation region centered at $\mathbf{r}$, it qualifies as a path, i.e., a valid connection between the source and the observation region. The Green’s function $G$ is numerically constructed as prislan:2016

[TABLE]

where

• •

$N$ is the total number of sound rays emitted from the point souce at ${\bf r}_{0}$,

• •

$R$ is the radius of the spherical observation region centered at $\mathbf{r}$,

• •

$L$ is the total number of paths $l$ with lengths $d_{l}$, that are identified by the simulation,

• •

$B_{l,i}$ is the complex amplitude of the sound ray of path $l$ after $i$-th reflection.

The factor $d_{l}$ represents the density of paths for the three-dimensional case of flat or piecewise flat boundaries, where (de)focusing effects are absent prislan:2016 .

C FEM method

The same room and BCs were also modeled using the PDE FEM solver COMSOL. There we used a monopole source with acoustic power of 1 W to obtain the complex pressure response $G_{p}$ to a source with unit acoustic power. Taking into account the connection Eq. (15), $G_{p}$ and the response $G$ obtained by RTS can be directly compared.

IV Results and discussion

The dimensions of the room were $a\times b\times c=4.215\,\textrm{m}\times 3.647\,\textrm{m}\times 3\,\textrm{m}$ as in our previous study (prislan:2016, ). To exclude that any conclusions would be based on the particularity of a chosen position, the simulation was performed for two pairs of source/observation points, Fig. 2, designated as corner and inner positions.

Rather than fixing the number of reflections of each ray, we introduce a low-amplitude cutoff: a ray sequence is terminated when $|B_{i}|<10^{-6}$ independent of the frequency. For the chosen geometry and absorption strengths this gives typically between 500 and 900 reflections per sound ray.

The size of the observation region was $R=0.2$ m and the number of emitted sound rays was $N=6\,758\,400$. Without parallelization, this simulation would be running for $\sim$ a month on a single core processor of a desktop computer (no code optimization). In our case the code was parallelized and was running on a cluster for 4 days. We however emphasize that practically valuable results are obtained already with drastically less rays and reflections, reducing the computational resources to a fraction.

Four BCs were tested: real impedance (Sec. B), resonance (Sec. C), extended reaction (Sec. D) and mixed. In the first three cases the same BC is applied to all boundaries (homogeneous BC), whereas in the mixed case the three BC are distributed as shown in Fig. 2.

Adjusting the rigid wall position of the extended reaction boundary such as to compensate for the 5 cm layer of porous material, the inner room dimensions were the same for all studied cases. In this way the rays were reflected from the same planes in all cases and only one RTS simulation was needed, in which the four studied BCs were computed simultaneously. To the contrary, for each BC a separate FEM simulations had to be performed in COMSOL.

The pressure response $G_{p}$ was calculated with both methods at the same frequencies in 0.05 Hz steps up to 300 Hz. The FEM meshing was performed automatically with the maximum allowed mesh spacing set to 0.34 m. This limit granted at least 5 mesh points per wavelength up to 200 Hz. An additional verificatory simulation with a much finer mesh did not deliver any observable changes.

A Pressure response

In Fig. 3, the magnitudes of the pressure responses obtained with RTS and FEM for the four BCs and both positions are presented in dB scale up to 300 Hz. We observe the following.

• •

For all studied cases the normalization of $G_{p}$ is correct. In most cases also the peaks and deeps of $|G_{p}|$ coincide well in frequency, level and width.

• •

A general exception is present for the lowest room modes (the peaks below 70 Hz) which are systematically slightly underestimated by RTS. This phenomenon was observed and discussed already in the previous study prislan:2016 and seems to be an inherent feature of the semiclassical approximation.

• •

A noisy behavior of RTS can be observed in frequency zones where the magnitude of $G_{p}$ is lower. This might be an indication of a noise floor that is not of statistical nature, but rather of numerical origin, i.e., it is not readily lowered by using more rays. A possible source of such numeric noise could be the finite floating point accuracy of determining the reflection angles and coordinates of boundary impacts. Although relatively small for double precision, the error increases with every reflection and propagation and is amplified by two or three orders of magnitude over the whole paths. Since the paths are long but at the same time confined to the room, the rather small relative error of a long path can translate into a substantial error relative to the size of the room.

• •

Due to the higher absorption of the extended reaction boundaries above 120 Hz (see Fig. 1), $|G_{p}|$ is lower in this frequency range in comparison with the other BCs. Therefore, again, more noise in $G_{p}$ computed by RTS can be observed, which is especially evident for the inner position. Moreover, for this BC the curves appear more erratic even for $|G_{p}|$ as high as $\sim 100$ dB, which is not observed for the other BCs. As the extended reaction BC is more involved, i.e., it requires more algebraic operations and function evaluations to determine the reflection coefficient of each reflection, this supports the finite floating point accuracy as a candidate responsible for this error.

• •

Moreover, a noisy behavior of RTS for the extended reaction BC is observed also in the frequency range below 40 Hz, most probably owing to the lack of absorption there ($\alpha_{\textrm{diff}}<0.05$, Fig. 1). As this noise is present only for very low absorption, it is not perceived as a practical limitation of the RTS method.

Further investigations where we modeled more reflections have revealed that this noise is reduced only if also more rays are emitted at the same time. To explain this observation, we should bear in mind that rays diverge during propagation and hence their population density decreases. For the noise to average out, it is required that portions of the phase space (prislan:2016, ) that significantly contribute to $G$ are adequately represented. Therefore, in case of low absorption when the amplitude of rays remains higher even after many reflections and thus still significantly contributes to $G$, more rays must be emitted to keep their population density adequate also in these distant parts of the paths.

• •

For the mixed BC, in the corner position the agreement between both methods is equally good as for the other BCs and no peculiarities are observed. In contrast, in the inner position deviations are severer. In particular the deeps between 160 Hz and 250 Hz do not coincide. We find this phenomenon interesting, as it is not systematic, i.e., it occurs only in one position and in a limited frequency range. Performing a FEM simulation on a finer mesh on one hand, and testing the influence of RTS parameters (number of reflections, number of emitted sound rays, radius of the observation region) on the other did not bring any evident improvement. At this stage we conclude that this peculiarity should be systematically addressed in further more specific studies – by testing more positions and more variety in the distribution of different BCs.

In case of the extended reaction BC the peaks occur at lower frequencies in comparison to the other two BCs. This frequency shift results from the fact that i) for the extended reaction BC the linear size of the room without the inside porous layer is $\sim 3$ % larger and ii) the sound speed in the porous layer is lower ($\sim 1$ % at 10 Hz up to $\sim 8$ % at 300 Hz). Both effects add up in lowering the resonance frequencies. The frequency shift is the smallest for the lowest modes and increases for higher modes. It is most evident in the central position for the peaks in the 90-170 Hz range. In the mixed BC case, the shift is less visible.

B Decay curves and reverberation time

The impulse response was calculated as an inverse Fourier transform of the octave-band-filtered (class 1 filters in agreement with EN 61260-1:2014 EN61260 ) pressure response. The two octave bands considered are indicated in Fig. 3. Next, the 0.8 s window of the squared impulse response was reverse integrated schroeder:1965 and the decay curve was obtained, Fig. 4. By a linear fit to the -5 dB to -25 dB part of the decay curve, $T_{20}$ was calculated as defined by the standard ISO3382 ISO3382 , serving as an estimate for the reverberation time.

For comparison the reverberation times were calculated also by the Sabine equation using absorption coefficients as octave band averages of $\alpha_{\textrm{diff}}$ given by the Paris formula Eq. (4). The so obtained average absorption coefficients $\bar{\alpha}_{\rm diff}$ for each BC are listed in Table 1.

The decay curves and reverberation times for both simulation points are presented in Fig. 4. The following can be observed.

• •

The obtained dynamic range is adequate to determine T20 everywhere with the exception of the RTS extended reaction BC inner position in the 125 Hz band. This case corresponds to the region of low $|G_{p}|$, where the highest amount of noise is present, Fig. 3. As a consequence, the reverberation time determined by RTS is highly overestimated and is thus irrelevant.

The largest discrepancy between the RTS and FEM $T_{20}$, $17\,\%$-$227\,\%$, is found for the case of the extended reaction BC. Again it is due to the noisy RTS response, Fig. 3, which restrains the dynamic range of the decay curves. In contrast, for the other BCs the average (maximum) discrepancy between the RTS and FEM reverberation times is $15\,\%$ ($31\,\%$). Particularly, for the 125 Hz band the average discrepancy is only $11\,\%$.

• •

The decay curves are most linear for the real impedance BC in the 125 Hz band. In these cases the difference in $T_{20}$ computed by FEM and RTS is only $3\,\%$ and $5\,\%$.

• •

Characteristically, for the 63 Hz band the $T_{20}$ computed by RTS is shorter in comparison with FEM for all studied cases. This is due to the underestimation of the $|G_{p}|$ peaks by RTS, Fig. 3. In fact, lowered peaks correspond to an increased relative width of the modes, which imitates an increased absorption and thus a shorter reverberation time.

• •

As expected, the Sabine reverberation time, which is by definition equal for both positions, is closer to $T_{20}$ determined by RTS and FEM for the corner position, where there are more peaks in the pressure response, Fig. 3, and the modal overlap is higher. At a higher overlap the diffuse sound field assumption – the basic assumption behind the Sabine equation – is better justified.

It thus follows that the RTS method is suitable for the determination of the reverberation time, provided the dynamic range of the decay curve is adequate. If the frequency band includes the lowest room modes, the RTS reverberation time can be underestimated by up to 30 % in absolute. Its relative behavior depending on the positions of the source and the listener nevertheless remains significant also in this case.

V Conclusion

We have performed a stark test of the acoustic RTS (ray-tracing semiclassical) method for a wide selection of realistic boundary conditions, confronting it face-to-face with the results of a full wave FEM method which is considered accurate in the studied lower frequency range.

First of all, we have shown that the pressure responses obtained by RTS and FEM agree in their normalization, i.e., there is no scaling uncertainty in the RTS-computed response. Moreover and most importantly, for the majority of the studied cases the peaks and deeps of the spatially dependent responses coincide remarkably well in frequency, level and shape. The agreement is systematic for all studied BCs with the exception of the mixed BC in the inner position, where in a limited frequency range the deeps in the response do not coincide.

Some further particularities of the RTS method were exposed – the underestimation of the peaks corresponding to the lowest modes and the presence of noise in case of low $|G_{p}|$ values as well as for very reflective boundaries. Generally we can conclude, that for the case of the rectangular room and tested BCs, RTS is an adequate method for modeling the spatially dependent pressure response.

Excluding the case of the pure extended reaction BC, the average discrepancy between FEM- and RTS-computed reverberation times was $15\,\%$. Owing to the weak absorption at low frequencies, the extended reaction BC produced a strong noise in this part of the pressure response, which restrained the dynamic range of the decay curves. When the extended reaction BC was combined with other BC types that supplied enough absorption at low frequencies, this noise problem disappeared.

Best agreement of the estimated $T_{20}$ was obtained for the real impedance BC in the 125 Hz band. In this case the FEM- and RTS-computed reverberation times differed only by 3 % and 5 % for both simulation positions, respectively. Overall, the study showed that RTS is an adequate method for predicting the reverberation time in rectangular rooms at low frequencies. At the same time, the importance of the adequate dynamic range of the computed decay curves was exposed.

In conclusion, our study validates the RTS method for a set of practically relevant boundary conditions, thereby further promoting it for low-frequency geometrical acoustic modeling in rectangular rooms. Some uncertainties in the case of inhomogeneous BCs were revealed, which should be systematically addressed in the future. Future studies should also validate and characterize the RTS method for general room geometries, including in particular (de-)focusing effects of curved boundaries.

Acknowledgments

The Authors acknowledge partial support of the Slovenian Research Agency (Grant J1-7435), the Slovenian Ministry of Higher Education, Science and Technology and the European Regional Development Fund.