Thin layers of porous media are widely adopted in sound absorption and noise control applications due to their compact arrangement. Particularly, porous media with low flow resistivity exhibit complex, non-local reaction behavior. Therefore, sound field prediction above these media is computationally challenging. This is due to singularities and branch points in the expression for the reflection coefficient. This paper introduces a framework based on the direct discrete complex image method to analyze the sound field above a rigid-backed, non-locally reacting porous sample. In contrast to traditional complex image applications, the proposed framework avoids the extraction of the quasi-static term and the poles from the reflection coefficient. Instead, the reflection coefficient—including its singularities—is directly approximated in terms of a series of complex exponentials, whose coefficients are determined with the matrix pencil method. This study analyzes the sound field above two test samples made of melamine and rockwool. The sound field computation is efficient and accurate in the near- and in the far-field. Moreover, predicted specific impedances agree well with experimental in situ impedance measurements. The proposed framework serves for more applications including object detection in multilayered porous grounds or sound propagation prediction in layered atmosphere.

Porous materials are an effective tool in noise and reverberation control. Sound propagation inside and above such materials generally depends on their complex pore structure. Simplified locally reactive material behavior can be considered for samples with high flow resistivity, i.e., σ ≥ 30 000 Ns/m4.1 However, many materials in engineering applications, e.g., foams or fibrous materials, are characterized by low flow resistivity and therefore tend to possess non-local reaction behavior. Consequently, the sample's surface acoustic impedance and the sound absorption coefficient depend on the angle of incidence of an impinging sound wave.1,2

There are well-known methods to calculate the sound field caused by a harmonic monopole above a locally reactive surface such as the model of Nobile and Hayek3 and the pole subtraction method proposed by Chien and Soroka.4 Di and Gilbert5 presented an exact formulation for local and extended reaction semi-infinite grounds based on the Laplace transform of the reflection coefficient. Li et al.6 stated that a non-locally reactive sample of finite thickness with rigid backing can be approximated by a semi-infinite ground with extended reaction behavior, if the sample is thick enough with respect to its intrinsic wave number. However, the need for lightweight and compact design often demands for thin layers. Consequently, non-locally reactive material behavior has to be considered in these cases.

Li et al.6 and Dragonetti and Romano7 tried to overcome the dependence of the surface impedance of a rigid-backed layer with the incidence angle by modifying the well-known models for locally reactive material behavior.3–5 The literature shows that this is well accomplished by the modified models, but all of them lack accuracy at low frequencies and when the source is placed near to the sample under investigation, since the spherical character of the radiated wave dominates under these circumstances.7,8 Yet, such a source-receiver setup is of particular interest for in situ surface impedance estimation.1,7,9

Dragonetti and Romano7 showed that the formulation of Allard et al.10 is the most accurate model for sound field computations above a non-locally reactive sample backed by a rigid wall, since it is based on the velocity potential fields above and inside the porous sample. In this case, the diffractive integral term can be evaluated by using direct numerical integration. However, the integrand has to undergo multiple mathematical transformations prior to the numerical integration,1,10 but even after mathematical manipulation, singularities remain in the integrand.1 These singularities are notably poles of the reflection coefficient, yielding surface waves in the physical space, which primarily contribute to the far-field.11–13 Accuracy and stability of numerical integration in the far-field can be increased by providing information on the integrand's singularities to the integration routine.1 Yet, finding these singularities numerically and extracting them is computationally challenging and expensive. This holds in particular for multilayered media, as the singularities cannot be determined analytically.14 

Alternatively, the fast field program15 can be used to calculate the sound field in the far-field. Li16 introduced the integrated exact solution (IES), which combines the direct numerical integration solution for the near-field with the fast field program solution for the far-field. By this, a numerically exact but computationally expensive solution can be obtained.17,18 This makes it impractical for optimization tasks, such as the inference of the surface impedance.1 

For this reason, recent research concentrates on efficient approximation methods for Sommerfeld integrals. Based on the work of Brekhovskikh and Godin,19 Liu and Li20,21 presented an improved integration scheme along the steepest descent path to obtain the sound field above a porous sample of finite thickness backed by a rigid wall or a surface impedance. This method is effective for long distances between source and receiver as well as at high frequencies. Further, Li and Tao22 proposed an accurate asymptotic solution of the diffraction integral using the modified saddle point method. By this, no singularity extraction from the diffraction integral kernel is required. Yet, this only holds for model setups considering normal incidence. Recently, Wang et al.23 and Li and Wang24 extended the application of the saddle point method in conjunction with the pole subtraction method to estimate the sound field radiated from a point source moving at constant speed above extended-reacting surfaces.

In the field of electromagnetics, the discrete complex image method (DCIM) has attracted interest as an efficient alternative to approximate Sommerfeld integrals11,17,25 by computing the integral kernel in terms of a series of damped exponentials.26 In the field of acoustics, however, few researchers have addressed the DCIM for sound field analysis with applications limited to either locally reactive behavior27–30,32 or extended reaction behavior of semi-infinite grounds.27,34,35 The initial work of Li and White27 focused on the principle of complex images to evaluate the sound field above an impedance plane and a semi-infinite ground with extended reaction behavior. Fawcett34 analyzed the acoustic scattering from a seabed, which was either modeled as a half-space with extended reaction behavior or an acousto-elastic half-space. The DCIM served to approximate the Green's function for these half-space problems and a non-linear least squares optimization was used to find the coefficients of the complex images. Fawcett35 extended his work to calculate the Green's function modeling a bounded wave guide and making it applicable to layered half-spaces. Ochmann and Piscoya29,30 demonstrated the incorporation of the DCIM for modeling an impedance boundary condition in boundary element simulations. They analyzed the sound field above a locally reactive plane and showed the advantages of the DCIM over conventional steepest descent approaches. Zheng et al.31 incorporated the principle of complex images into a fast multipole boundary element method to efficiently solve large-scale locally reactive half-space problems. Brandão and Lenzi32 reported on the advantages of the DCIM in in situ surface impedance inference for a locally reactive material. A general overview of the DCIM and its use for acoustic half-space problems can be found in the recent book of Ochmann and Piscoya.33 

All these complex image approaches rely on the traditional DCIM, for which the quasi-static term, i.e., the limiting value of the reflection coefficient for the intrinsic wave number in air tending to zero, limk00R̃ (contributing to the near-field) and the poles of the reflection coefficient (contributing to the far-field) have to be identified and extracted from the Sommerfeld kernel.14 This method is impractical since these poles, which lead to surface waves in the physical space, are difficult to determine analytically or numerically. Further, errors occur in the near-field, due to the singularity at the origin of the Hankel function.11,13,17 Another drawback of the traditional DCIM is the use of polynomial methods, in particular Prony's method,36,37 to determine the coefficients for the complex images. Prony's method shows unstable behavior with an increasing number of complex images.26 Yet, a multilayered ground demands for an increased number of complex images due to its challenging reflection behavior,35 which also holds for a rigid-backed, non-locally reacting sample of finite thickness.

This work introduces a robust and computationally efficient framework to analyze the sound field above a non-locally reacting layer with rigid backing. It is based on the direct DCIM14 and uses the matrix pencil method26,38 (MPM) to calculate the coefficients for the complex images. To the best of the authors' knowledge, this is the first time the DCIM has successfully been applied to predict the sound field above a rigid-backed, non-locally reacting sample of finite thickness.

This paper is structured as follows: Sec. II presents the general equations for the sound field above a non-locally reactive sample of infinite radial extent backed by a rigid wall. Further, the fundamentals of the direct DCIM are demonstrated. Particular attention is given to the MPM, which is the backbone of the direct DCIM. Results of numerical studies and experimental validation of the proposed method are displayed in Sec. III. Section IV concludes this article and outlines further applications of the proposed framework.

The outlined theory here considers the sound field above a rigid-backed, non-locally reactive sample of finite thickness and infinite radial extent (Sec. II A) and the procedure to estimate the sound field with use of the direct DCIM (Sec. II B). Moreover, the process to calculate the coefficients of the complex images by using the MPM is demonstrated (Sec. II C).

A radiating monopole above a partly reflective plane of infinite radial extent is a widely used model to estimate the sound absorption properties of porous layers in situ.1,7,27,39–41 Figure 1 shows the model setup for a porous layer of finite thickness d above a rigid backing. The source is located at rs=(xs,ys,zs) in the upper, homogeneous half-space made of air and radiates with angular frequency ω. Assuming a harmonic time dependence of eiωt with i=1, the sound pressure at a prescribed receiver position r=(xr,yr,zr) can be expressed by27 

p(r,rs)=SG(r,rs),
(1)

given a source strength S. A unit source strength of S = 1 kg/s2 is assumed for all further analyses. For the sake of convenience, the problem is reformulated in cylindrical coordinates19 and the radial distance between source and receiver r=(xrxs)2+(yrys)2 is introduced.

FIG. 1.

Model setup of a harmonic monopole source above a rigid-backed, non-locally reactive porous sample of finite thickness d.

FIG. 1.

Model setup of a harmonic monopole source above a rigid-backed, non-locally reactive porous sample of finite thickness d.

Close modal

The specific half-space Green's function G can be stated as27 

G(r,rs)=eik0D1D1+i0R̃(kr)eiKz,0(zr+zs)krJ0(krr)Kz,0dkr,
(2)

with the first term representing the direct contribution from the source to the receiver as a free-field Green's function with the direct distance D1=r2+(zrzs)2 between source and receiver. The intrinsic wave number of air is k0=ω/c0, with c0 being the speed of sound in air. The zeroth order Bessel function is represented by J0. The second term in Eq. (2) accounts for the reflected wave. With the presented method, there is no need to separate this term into a reflective (quasi-static) term and a diffractive integral term. Instead, the integrand remains in its spectral domain representation,42 which is expressed in terms of the vertical wave number in air

Kz,0=+k02kr2,
(3)

and the radial wave number kr. The vertical and the radial wave number can be interpreted as a decomposition of the wave number vector in vertical and radial coordinate directions. Formally, this manifests in the Sommerfeld identity43 

eik0|rrs||rrs|=i0e±iKz,0(zrzs)Kz,0krJ0(krr)dkr,
(4)

which represents the decomposition of a spatial spherical wave into a plane wave propagating in vertical direction and a cylidrical wave propagating in radial direction in the wave number domain. The negative (positive) sign in the integrand is selected if the source is located above (below) the receiver.27 

The contribution and mathematical complexity of the reflected wave depends on the plane wave reflection coefficient R̃(kr). For a single layer of finite thickness d backed by a rigid wall, it can be expressed as6 

R̃(kr)=(iKz,0ρ0)(Kz,1ρ1)tan(Kz,1d)(iKz,0ρ0)+(Kz,1ρ1)tan(Kz,1d).
(5)

In this equation, ρ0 and ρ1 represent the density of air and the complex-valued density of the porous medium, respectively. The vertical wave number in the porous medium is

Kz,1=+k12kr2,
(6)

which depends on the complex-valued intrinsic wave number k1 of the porous medium. Introducing the reflection coefficient [Eq. (5)] in the integral term in Eq. (2), there is no closed-form solution for the Green's function.1 The integrand is mathematically challenging, due to its poles and zeros,1 the numerous discontinuities in the tangent function and the highly oscillating behavior of the Bessel function.17 This raises the need to manipulate the integrand to make it numerically tractable with conventional methods, e.g., numerical integration.

The direct DCIM as introduced by Yuan et al.14 promises to overcome the need to extract the quasi-static term and the poles from the reflection coefficient R̃(kr) [Eq. (5)], which is directly approximated by a series of damped exponentials

R̃(kr)n=1NanebnKz,0,
(7)

with expansion coefficients an and bn. Inserting the approximated reflection coefficient [Eq. (7)] into the integrand in Eq. (2) yields

G(r,rs)=eik0D1D1+in=1Nan0eiKz,0(zr+zs+ibn)krJ0(krr)Kz,0dkr.
(8)

The integral in the reflected wave term can be integrated exactly using the Sommerfeld identity [Eq. (4)], which leads to a closed-form expression for the specific half-space Green's function

G(r,rs)=eik0D1D1+n=1Naneik0DnDn.
(9)

Each summand represents a spherical wave with complex-valued amplitude an, which is radiated from an image source located at the complex-valued distance Dn=r2+(zr+zs+ibn)2 from the receiver. As the coordinates of the complex images solely differ in their values for bn, all complex images are located at the same real-valued position, but can be interpreted as a line source along the additional imaginary axis. This makes this formulation inherently convenient to include the cylindrical wave nature of surface waves as a contribution from a line source given the total number of image sources N in Eq. (9) is sufficient.14,44

Equation (9) serves for an efficient evaluation of the Green's function. Remark that the only approximation error occurs in the approximation of the reflection coefficient [Eq. (7)], while all other operations are analytically exact. Since the reflection coefficient is directly approximated without any terms extracted, Eqs. (7)–(9) differ from their complements in traditional DCIM formulations.27–30,32,34,35

In order to compute the coefficients an and bn, it is necessary to sample the reflection coefficient in a real-valued variable.14 In traditional DCIM applications, a sampling path in the complex Kz,0-plane27 

Kz,0|C=k0[it+(1tT0)],0tT0
(10)

is used. The reflection coefficient has to be evenly sampled in the real-valued variable t along this path to allow for the use of Prony's method, while the end of this path T0 should be chosen larger than the branch points to avoid singularities. This is not necessary for the direct DCIM, as the MPM (Sec. II C) can estimate the behavior beyond T0 given sufficient information on the reflection coefficient's singularities is included in the sampling path.14 For this reason, a modified sampling path

Kz,0|C12=γk0[it+1],0tT0γ,0<γ1,
(11)

with the scaling parameter γ has been proposed by Sarkar et al.14 for the direct DCIM.

The reflection coefficient should be sampled evenly along the modified path such that all possible singularities (in the multilayered case) are equally accounted for during the sampling and the reconstruction process using the MPM. Figure 2 qualitatively compares two variants of the modified sampling path for the direct DCIM, C2 for γ = 1 and C12 for 0<γ<1 in Eq. (11), with the path used in traditional DCIM applications C [Eq. (10)]. The modified path runs parallel to the imaginary axis, which simplifies the relationship between the variable t and the vertical wave number in air Kz,0 compared to traditional DCIM applications.

FIG. 2.

Comparison of sampling paths of the traditional discrete complex image method (gray solid), the direct discrete complex image method (black dashed) and the extended range direct discrete complex image method (black solid) in the complex Kz,0-plane (a) and in the complex kr-plane (b).

FIG. 2.

Comparison of sampling paths of the traditional discrete complex image method (gray solid), the direct discrete complex image method (black dashed) and the extended range direct discrete complex image method (black solid) in the complex Kz,0-plane (a) and in the complex kr-plane (b).

Close modal

Figure 3 shows the reflection coefficient of the melamine foam sample (see Table I) at 1000 Hz sampled along the modified path C2, i.e., γ=1.0 [Fig. 3(a)] and along a shifted modified path C12 with γ=0.05 [Fig. 3(b)]. Smaller values of the scaling parameter γ move the modified path closer to the imaginary axis, which leads to more variation in the sampled reflection coefficient. This correlates to more information on the reflection coefficient's singularities in the samples used in the MPM, which is pivotal to account for the contribution of surface waves in particular in the far-field. By this, accurate sound field estimations can be achieved for higher distances between source and receiver. In other words, a smaller value of γ extends the accuracy of the direct DCIM in the far-field.14 However, shifting the modified path closer to the imaginary axis comes at the cost that the sampling path does no longer start at the origin of the complex kr-plane [Fig. 2(b)]. Consequently, a compensation sampling path along the real axis—C11—has to be introduced if γ<1. Therefore, the procedure to calculate the coefficients for the complex images has to be slightly modified.14 The modified calculation procedure includes the following steps:

FIG. 3.

Real part (solid line) and imaginary part (dashed line) of the reflection coefficient [Eq. (7)] of the melamine foam sample (Table I) for normal incidence at 1000 Hz sampled along the modified sampling path [Eq. (11)] for γ = 1 (a) and γ=0.05 (b).

FIG. 3.

Real part (solid line) and imaginary part (dashed line) of the reflection coefficient [Eq. (7)] of the melamine foam sample (Table I) for normal incidence at 1000 Hz sampled along the modified sampling path [Eq. (11)] for γ = 1 (a) and γ=0.05 (b).

Close modal
TABLE I.

Material parameters for test samples under investigation.

Flow resistivityPorosityTortuosityViscous char. lengthThermal char. length
σ (Ns/m4)ϕ (—)α (—)Λ (m 106)Λ (m 106)
Melamine foam 9207 0.99 1.0 300 600 
Rockwool 11 817 0.97 1.0 113 122 
Flow resistivityPorosityTortuosityViscous char. lengthThermal char. length
σ (Ns/m4)ϕ (—)α (—)Λ (m 106)Λ (m 106)
Melamine foam 9207 0.99 1.0 300 600 
Rockwool 11 817 0.97 1.0 113 122 
  1. Select the scaling parameter γ according to the desired maximum radial distance between source and receiver.

  2. Generate M evenly distanced samples of the reflection coefficient along the modified path C12 and approximate the reflection coefficient with a first series expansion
    R̃(kr)|C12n=1NanebnKz,0|C12=n=1NAneBnt|C12
    (12)

    in terms of the coefficients in the transformed space An and Bn.

  3. Set the primary truncation threshold tol and apply the MPM (Sec. II C) to the samples of Eq. (12) to get the coefficients an, bn and the number of primary complex image sources N; note that N < M.

  4. Generate M samples of the residual along the compensation path C11 (Fig. 2) and approximate the residual in terms of a second series expansion
    R̃(kr)|C11n=1NanebnKz,0|C11n=1NanebnKz,0|C11.
    (13)
  5. Set the compensation truncation threshold tol and apply the MPM to the samples of Eq. (13) to calculate the coefficients an,bn and the number of compensation complex images N; note that N<M.

  6. Hence, for 0<γ<1, the expression for the Green's function [Eq. (9)] changes to
    G(r,rs)=eik0D1D1+n=1Naneik0DnDn+n=1Naneik0DnDn,
    (14)

    with Dn=r2+(z+zS+ibn)2.

In the following, Eq. (14) is denoted as the extended range solution, in which the total number of complex images is the sum of the number of primary image sources N (step 3) and compensation image sources N (step 5).

For applications in the near-field, it suffices to use the modified path C2 (Fig. 2), i.e., γ = 1 in Eq. (11). In this case, no compensation path is required and the two-path procedure reduces to an approximation of the reflection coefficient along a single path. By this means, only steps 1–3 are relevant and the Green's function can then be efficiently evaluated using Eq. (9). We recommend the direct DCIM procedure for near-field analyses in the range of k0r<10, as proposed by Li and White.27 

The MPM is the backbone of the direct DCIM, since it enables the direct approximation of the reflection coefficient [Eq. (7)]. Compared to Prony's method, the MPM is a single-step process to calculate the expansion coefficients an and bn (step 3) and in a second function call an and bn (step 5) in case the extended range solution is used.26 In the following, the procedure for the MPM is outlined for process step 3, but it holds in similar manner for step 5 in which all primary variables have to be replaced by their compensation counterparts, i.e., MM and toltol.

The MPM takes M evenly distanced samples of the reflection coefficient generated along the modified path R̃(kr)|C12=s[j] with j=1,2,,M (or M evenly distanced samples of the residual along the compensation path C11) as an input. In an initial step, these samples are arranged in a rectangular Hankel matrix of dimension (ML)×(L+1)

Y=(s[1]s[2]s[L]s[L+1]s[2]s[3]s[L+1]s[L+2]s[ML]s[ML+1]s[M1]s[M]),

with L denoted as the pencil parameter. We set the pencil parameter to L=(M/2)1, according to the recommendation of Sarkar et al.26,45 Next, a truncated singular value decomposition of the Hankel matrix Y is performed to determine the dominant singular values fulfilling the condition σj>max(σj)×tol. Only the columns of the diagonal matrix containing these dominant singular values in descending order and their corresponding right singular vectors are kept in the filtered matrices Σf and Vf, while all other singular values are discarded. This truncation is an essential step in the direct DCIM, since the number of dominant singular values defines the number of primary complex images N (or compensation images N), i.e., N=rank(Σf). Then, the matrix pencil Y2λY1 is formed with the submatrices Y1=UΣfV1fH and Y2=UΣfV2fH. The matrices V1fH and V2fH result by eliminating the last row and the first row of Vf, respectively.26 The matrix U contains all left singular vectors of the full singular value decomposition. Finally, the dominant signal poles zn are obtained by solving the general eigenvalue problem

zn=eigval(Y1Y2),
(15)

where Y1 is the Moore-Penrose pseudoinverse of the first submatrix, which is used to guarantee a stable and efficient eigenvalue analysis, as both submatrices are ill-conditioned.38 The damping coefficients in the deformed space

Bn=ln(zn)Δt
(16)

follow directly from the signal poles and the sampling interval Δt. The unknown amplitudes in the deformed space An follow from the least squares problem

(z10z20zN0z11z21zN1z1M1z2M1zNM1)(A1A2AN)=(s[1]s[2]s[M]).
(17)

Once the coefficients in the transformed space An and Bn are found, retransformation yields the original coefficients for the complex images

an=AneiBnandbn=iBnk0.
(18)

Note that the set of coefficients an, bn (and an,bn as well) is independent of the source and of the receiver position, as stated by Li et al.27 By this, any source-receiver combination can efficiently be evaluated with the same set of coefficients for the complex images. Moreover, the direct DCIM enables efficient evaluation of the particle velocity with the same set of coefficients by using the linearized Euler equation. Thus, the specific impedance can be evaluated at any receiver position above the porous sample without any additional effort.

Algorithm 1 summarizes the proposed framework to calculate the specific half-space Green's function based on the direct DCIM. It should be mentioned that the presented method is only valid, when source and receiver are located above the porous ground. However, the presented framework can be easily extended to analyze sound propagation within the porous layer. For this purpose, the derivative of the Sommerfeld identity is required.14 

Algorithm 1:

Pseudocode of the direct discrete complex image method.

SET scaling parameter γ according to desired maximum radial distance 
GENERATEM samples of reflection coefficient along modified  path C12 [Eq. (11)
CALL MPM to compute coefficients an, bn and number of primary  complex images N 
IF0<γ<1THENextended range solution 
  GENERATEM samples of residual [Eq. (13)] along compensation  path C11 
  CALL MPM to compute coefficients an,bn and number of  compensation complex images N 
  COMPUTE Green's function [Eq. (14)
ELSEnear-field solution 
  COMPUTE Green's function [Eq. (9)
END IF 
SET scaling parameter γ according to desired maximum radial distance 
GENERATEM samples of reflection coefficient along modified  path C12 [Eq. (11)
CALL MPM to compute coefficients an, bn and number of primary  complex images N 
IF0<γ<1THENextended range solution 
  GENERATEM samples of residual [Eq. (13)] along compensation  path C11 
  CALL MPM to compute coefficients an,bn and number of  compensation complex images N 
  COMPUTE Green's function [Eq. (14)
ELSEnear-field solution 
  COMPUTE Green's function [Eq. (9)
END IF 

The proposed direct DCIM framework is applied to three test scenarios. First, the accuracy and efficiency of the direct DCIM are numerically validated for normal incidence. Second, numerical prediction of the specific impedance is validated against in situ impedance measurements in a hemi-anechoic chamber. Third, the extended range solution of the proposed method is analyzed numerically. Computations are conducted on a machine with eight physical cores operating at 3.30 GHz and with 96 Gb RAM.

Each scenario treats the sound field above a rigid-backed, non-locally reactive sample of finite thickness with source and receiver located at a height of zs=0.33 m and zr=0.015 m above the surface of the sample. The density and the sound speed of air are ρ0=1.21 kg/m3 and c0=343 m/s. A melamine foam sample and a rockwool sample are considered. Rigid frame porous material behavior is assumed and the intrinsic wave number of the sample k1 and its complex density ρ1 are calculated with the Johnson-Champoux-Allard model.46 Since this paper focuses on the capability of the direct DCIM to predict the sound field above the samples, this simpler material model is chosen, even though more elaborate material models have proven a better fit with respect to impedance tube measurements.1 Table I lists the parameters characterizing the considered samples, whose measurement is detailed in the work of Brandão et al.1 

The first two analyses deal with normal incidence of the radiated sound wave. In these cases, there is no radial offset between source and receiver, hence r = 0. Therefore, the near-field solution process for the direct DCIM with γ = 1 is used.

1. Numerical validation

A melamine foam sample of thickness d = 40 mm is modeled with material properties as listed in Table I. The solution with direct numerical integration (DNI) as presented by Brandão et al.1 is used as reference. The reflection coefficient is sampled along the modified path C2 with a sampling rate of Δt=0.1 in the interval 0<t<T0=7.5. The truncation tolerance to obtain the dominant singular values is set to tol=1010. The sound field is analyzed in a frequency range from 100 to 10 000 Hz in steps of 1 Hz.

Figure 4 shows the dimensionless sound pressure, i.e., the total sound pressure p normalized against the free-field sound pressure pff=S(eik0D1/D1) above the melamine foam sample obtained with the direct DCIM (dotted) and the reference solution (dashed) across the considered frequency range. Its real part is displayed in Fig. 4(a) and the imaginary part in Fig. 4(b). Reflection at the sample's surface and non-local acoustic material behavior lead to interference effects above the sample. The direct DCIM solution yields nearly an identical sound pressure compared to the reference solution across the whole frequency range. The MPM reveals a maximum number of N = 11 primary image sources to obtain the presented results. No compensation sources are calculated, since we use the near-field solution (see Algorithm 1) in case of normal incidence (r = 0). The proposed direct DCIM framework shows high efficiency in terms of computation time compared to the reference solution (DNI). Given one source-receiver setup, the sound pressure across the whole frequency range as displayed in Fig. 4 is obtained in 25.2 s using the direct DCIM approach compared to 750 s using direct numerical integration.

FIG. 4.

Comparison of real (a) and imaginary part (b) of the dimensionless sound pressure above the melamine foam sample of 40 mm thickness calculated with direct numerical integration (dashed) and with the direct discrete complex image method (dotted).

FIG. 4.

Comparison of real (a) and imaginary part (b) of the dimensionless sound pressure above the melamine foam sample of 40 mm thickness calculated with direct numerical integration (dashed) and with the direct discrete complex image method (dotted).

Close modal

Since, the direct DCIM solution and the reference solution are visually indistinguishable (Fig. 4), the accuracy of the proposed framework is further assessed with the normalized error in the total sound pressure

ε(fq)=|{pDNI(fq)pDCIM(fq)}|{pDNI}RMSandε(fq)=|{pDNI(fq)pDCIM(fq)}|{pDNI}RMS
(19)

individually for real and imaginary part at a discrete frequency fq. The absolute error in the sound pressure is normalized with the root mean square across the entire frequency range of either the real part {pDNI}RMS=1Nfq=1Nf({pDNI(fq)})2 or the imaginary part {pDNI}RMS=1Nfq=1Nf({pDNI(fq)})2 of the total sound pressure, with Nf representing the total number of discrete frequencies considered.

The normalized error is displayed in Fig. 5(a) for the real part and in Fig. 5(b) for the imaginary part. Both errors show a slight increase for frequencies above 2500 Hz, but remain below a value of 105 for the entire frequency range. Actually, the DCIM solution oscillates around the reference solution, which is a consequence of the approximation in terms of a series of complex exponentials. Yet, the accuracy of the proposed method is validated and its efficiency makes it a viable tool in optimization problems, e.g., in surface impedance estimation or material characterization.

FIG. 5.

Normalized errors for real (a) and imaginary part (b) in the total sound pressure above the melamine foam sample of 40 mm thickness.

FIG. 5.

Normalized errors for real (a) and imaginary part (b) in the total sound pressure above the melamine foam sample of 40 mm thickness.

Close modal

2. Experimental validation

Experimental validation is conducted for a 40 mm thick melamine foam sample with a sample size of 0.6 m × 0.6 m and a 50 mm thick rockwool sample with a sample size of 1.2 m × 0.6 m with properties listed in Table I. Numerical estimations of the specific impedance Zm above the samples are compared to in situ impedance measurements from previous studies.1,9 A free-field calibrated pressure-particle velocity sensor was used for measurements in a hemi-anechoic chamber with a volume of 100 m3 and a cut-off frequency of 150 Hz for the absorbing wedges. 30 sequential measurements with high reproducibility have been conducted for each sample under normal incidence given a source height of zs=0.33 m and a receiver height of zs=0.015 m. The reader is referred to the works of Brandão et al.1,9 for detailed information on the measurement setup. Regarding the direct DCIM framework, all parameters remain the same as in the previous section, i.e., γ=1.0,Δt=0.1,T0=7.5, and tol=1010.

Figure 6 shows the real (left column) and imaginary parts (right column) of the normalized specific impedance determined above the melamine foam sample [Figs. 6(a) and 6(b)] and the rockwool sample [Figs. 6(c) and 6(d)]. The measured impedances are displayed as the solid black line, the calculated impedances using the direct DCIM framework are displayed in gray, and the reference solution (DNI) is plotted as a black dashed line. The specific impedance shows similar trends for both materials. The numerical models predict the specific impedance well for the melamine foam sample. However, they underestimate the real part for the rockwool sample at low frequencies and show lesser fluctuation than the experimental data at higher frequencies. These deviations can be attributed difficulties in the experimental determination of the macroscopic parameters (Table I) and to the assumption of rigid-frame material behavior assumed with the Johnson-Champoux-Allard model.46 Better predictions, in particular in the low frequency regime might be obtained by using more elaborate material models, such as the Johnson-Champoux-Allard-Lafarge model,47 which has been identified as the most appropriate model for the materials analyzed,1 or the Biot-Allard model.48 These models include low-frequency limits of the dynamic material parameters or consider elastic behavior of the porous structure. Moreover, edge effects affect the measured data, in particular in the low frequency regime, since the measured samples are of finite size.49,50 In contrast, the numerical setup considers samples of infinite radial extent. Hence, deviations in the low frequency regime also stem from model errors. Again, for the specific impedance, the direct DCIM yields indistinguishable results compared to the numerical integration, which requires integration on both pressure and particle velocity and, therefore, can become computationally expensive. The presented framework shows good performance in estimating the specific impedance and it is inherently designed to serve as the core of a robust impedance inference algorithm.

FIG. 6.

Real and imaginary parts of the normalized specific impedance calculated with the direct discrete complex image method (solid gray) and with direct numerical integration (dashed black) compared with experimental data from in situ measurements (Ref. 1) (solid black) for a melamine foam test sample (upper row) and a rockwool test sample (lower row).

FIG. 6.

Real and imaginary parts of the normalized specific impedance calculated with the direct discrete complex image method (solid gray) and with direct numerical integration (dashed black) compared with experimental data from in situ measurements (Ref. 1) (solid black) for a melamine foam test sample (upper row) and a rockwool test sample (lower row).

Close modal

The proposed direct DCIM framework intrinsically also allows for fast and accurate sound field predictions in the far-field, as the coefficients for the complex images are independent of the source and the receiver position. The user only has to adjust the scaling parameter γ<1 customized to the associated application and the method follows the steps outlined in Algorithm 1 to obtain an accurate approximation of the Green's function.

The sound field above the rigid-backed melamine foam sample of 40 mm thickness is analyzed for varying radial distance between source and receiver. As in the previous analyses, the source and the receiver are located at a height of zs=0.33 m and zr=0.015 m above the surface of the layer. Studies are conducted at 100, 1000, 5000, and 10 000 Hz. Here, the IES according to Li16 serves as a reference to determine the accuracy of the direct DCIM approximation. The end of the modified path for the direct DCIM framework is set to T0=8.5 to include more information on the singularities. The sampling rate and the truncation threshold on the modified path C12 are Δt=0.01 and tol=1010. Note that the results for 100 Hz are obtained with a sampling rate of Δt=0.1 on C12, which allows one to analyze the effect of the number of samples on the calculation time. The residual [Eq. (13)] is sampled with M=50 samples along the compensation path C11 and a truncation threshold of tol=103 is prescribed for the MPM on the samples of the residual. Only a few compensation sources are required to approximate the behavior of the residual, since the residual does not contain a lot of variation as the reconstruction with the primary image sources does account for most of the reflection coefficient's variation.14 For this reason, the number of samples of the residual M and the truncation threshold tol can be set lower than for the primary step (step 3).

Figure 7 shows the estimated excess attenuation EA=20log10|p/pff| for a variation of the radial distance between source and receiver. The results for 100 Hz are displayed in Fig. 7(a), for 1000 Hz in Fig. 7(b), for 5000 Hz in Fig. 7(c), and for 10 000 Hz in Fig. 7(d). Sound field analysis is conducted for scaling parameters γ=1.0 (blue), γ=0.4 (magenta), and γ=0.05 (orange). Additionally, the IES reference solution is displayed (black dashed). The excess attenuation shows nearly no variation for small radial distances, higher variation with dips at 5000 and 10 000 Hz in an intermediate range and decreases rapidly for larger distances between source and receiver. The direct DCIM framework yields accurate results in the near-field for all considered scaling parameters. The single path solution (γ=1.0) deteriorates the earliest from the reference solution at each frequency. The solution with γ=0.4 yields accurate results up to higher radial distances, but only the solution with scaling parameter γ=0.05 gives accurate results across the entire range. Table II lists the number of primary and compensation complex images N and N, the maximum radial distance max(r)5%, when the direct DCIM solutions deteriorate 5% relative to the reference solution and the total calculation time τ for all direct DCIM solutions displayed in Fig. 7.

FIG. 7.

(Color online) Excess attenuation at 100 Hz (a), 1000 Hz (b), 5000 Hz (c), and 10 000 Hz (d) for a variation of the radial distance between source and receiver for the integrated exact solution (black dashed) and different variations of the direct discrete complex image method γ=1.0 (blue), γ=0.4 (magenta), and γ=0.05 (orange).

FIG. 7.

(Color online) Excess attenuation at 100 Hz (a), 1000 Hz (b), 5000 Hz (c), and 10 000 Hz (d) for a variation of the radial distance between source and receiver for the integrated exact solution (black dashed) and different variations of the direct discrete complex image method γ=1.0 (blue), γ=0.4 (magenta), and γ=0.05 (orange).

Close modal
TABLE II.

Results of the extended range analysis.

f in HzγNNmax(r)5%τ in s
100 1.0 10 n.a. 11.6λ 0.254 
 0.4 12 85.2λ 3.95 
  0.05 20 >408.8λ 5.66 
1000 1.0 11 n.a. 7.6λ 2.87 
 0.4 12 62.4λ 6.09 
  0.05 15 >2283.6λ 685.0 
5000 1.0 12 n.a. 55.6λ 2.92 
 0.4 13 111.6λ 6.24 
  0.05 15 >7301.9λ 674.0 
10000 1.0 13 n.a. 54.0λ 2.9 
 0.4 14 190.8λ 6.15 
  0.05 17 >8057.9λ 710.0 
f in HzγNNmax(r)5%τ in s
100 1.0 10 n.a. 11.6λ 0.254 
 0.4 12 85.2λ 3.95 
  0.05 20 >408.8λ 5.66 
1000 1.0 11 n.a. 7.6λ 2.87 
 0.4 12 62.4λ 6.09 
  0.05 15 >2283.6λ 685.0 
5000 1.0 12 n.a. 55.6λ 2.92 
 0.4 13 111.6λ 6.24 
  0.05 15 >7301.9λ 674.0 
10000 1.0 13 n.a. 54.0λ 2.9 
 0.4 14 190.8λ 6.15 
  0.05 17 >8057.9λ 710.0 

The proposed framework yields accurate results for larger distances between source and receiver if a smaller scaling parameter is chosen. Reducing the scaling parameter means that more primary image sources N are required to account for the increased variation in the sampled reflection coefficient. Further, smaller values of the scaling parameter lead to larger sampling ranges, 0tT0/γ, according to Eq. (11). Since the sampling rate Δt is held constant, this results in an increased number of samples M for the reflection coefficient along the modified path C12. Therefore, the dominant singular values have to be determined from a larger Hankel matrix and a larger least squares problem has to be solved for the amplitudes, which explains increased computation times. Coarser sampling can significantly reduce the calculation time (see results for 100 Hz in Table II), but this comes at the cost of a decreased maximum radial distance. The number of compensation image sources N is rather unaffected by the choice of the scaling parameter and only few compensation images (N5 in most cases) suffice for accurate predictions. This has similarly been observed by Yuan et al.14 Note that the frequency hardly affects the calculation time. By adjusting the scaling parameter, the proposed framework can be used for sound propagation predictions into the far-field, e.g., in underwater acoustics above a layered seabed or sound propagation in layered atmosphere. Further, we stress that the computation of an, bn and an,bn is independent of the position of sources and receivers and that the evaluation of Eqs. (9) or (14) is straightforward. The operator only has to choose the scaling parameter according to a desired maximum radial distance. The values for max(r)5% in Table II serve as an orientation for such tasks.

This article introduces the direct discrete complex image method to the field of acoustics. Compared to traditional complex image applications, neither the quasi-static term nor the poles have to be extracted from the reflection coefficient prior to the application of the complex image method. This avoids the cumbersome identification and extraction of singularities of the reflection coefficient. Instead, the reflection coefficient is directly sampled along a modified sampling path and the coefficients of the complex images are calculated with the matrix pencil method, which includes an automated selection of the number of image sources. With these attributes, the direct discrete complex image method allows for a simplified application compared to traditional discrete complex image method applications. Further, the proposed method allows for the analysis of computationally challenging scenarios, such as the sound field above a rigid-backed, non-locally reacting sample of infinite radial extent, which is validated against experimental in situ impedance measurements for a melamine foam and a rockwool sample. Accurate sound field predictions are obtained in the near- and in the far-field. The scaling parameter in the modified sampling path, which the user has to prescribe a priori, decides up to what distance the approximation yields accurate results. Smaller values for the scaling parameter enable the analysis to larger distances between source and receiver. The proposed method extends the use of the complex image method in the field of acoustics and makes it a viable alternative for computationally efficient predictions of the sound field above layered media, since the reflection coefficient is directly encountered without the need to extract any singularities. Future work includes the incorporation of the direct discrete complex image method in a robust impedance inference framework and half-space sound field analyses using the boundary element method.

1.
E.
Brandão
,
P.
Mareze
,
A.
Lenzi
, and
A. R.
da Silva
, “
Impedance measurement of non-locally reactive samples and the influence of the assumption of local reaction
,”
J. Acoust. Soc. Am.
133
(
5
),
2722
2731
(
2013
).
2.
R.
Dragonetti
,
C.
Ianniello
, and
R. A.
Romano
, “
Measurement of the resistivity of porous materials with an alternating air-flow method
,”
J. Acoust. Soc. Am.
129
(
2
),
753
764
(
2011
).
3.
M. A.
Nobile
and
S. I.
Hayek
, “
Acoustic propagation over an impedance plane
,”
J. Acoust. Soc. Am.
78
(
4
),
1325
1336
(
1985
).
4.
C.
Chien
and
W.
Soroka
, “
Sound propagation along an impedance plane
,”
J. Sound Vib.
43
(
1
),
9
20
(
1975
).
5.
X.
Di
and
K. E.
Gilbert
, “
An exact Laplace transform formulation for a point source above a ground surface
,”
J. Acoust. Soc. Am.
93
(
2
),
714
720
(
1993
).
6.
K. M.
Li
,
T.
Waters-Fuller
, and
K.
Attenborough
, “
Sound propagation from a point source over extended-reaction ground
,”
J. Acoust. Soc. Am.
104
(
2
),
679
685
(
1998
).
7.
R.
Dragonetti
and
R. A.
Romano
, “
Considerations on the sound absorption of non locally reacting porous layers
,”
Appl. Acoust.
87
,
46
56
(
2015
).
8.
R.
Dragonetti
and
R. A.
Romano
, “
Errors when assuming locally reacting boundary condition in the estimation of the surface acoustic impedance
,”
Appl. Acoust.
115
,
121
130
(
2017
).
9.
E.
Brandão
,
E.
Tijs
,
A.
Lenzi
, and
H.-E.
Bree
, “
A comparison of three methods to calculate the surface impedance and absorption coefficient from measurements under free field or in situ conditions
,”
Acta Acust. Acust.
97
,
1025
1033
(
2011
).
10.
J.-F.
Allard
,
W.
Lauriks
, and
C.
Verhaegen
, “
The acoustic sound field above a porous layer and the estimation of the acoustic surface impedance from free-field measurements
,”
J. Acoust. Soc. Am.
91
(
5
),
3057
3060
(
1992
).
11.
S.
Teo
,
S.
Chew
, and
M.
Leong
, “
Error analysis of the discrete complex image method and pole extraction
,”
IEEE Trans. Microwave Theory Tech.
51
(
2
),
406
413
(
2003
).
12.
J. F.
Allard
,
M.
Henry
, and
V.
Gareton
, “
Pseudo-surface waves above thin porous layers (L)
,”
J. Acoust. Soc. Am.
116
(
3
),
1345
1347
(
2004
).
13.
E.
Simsek
,
Q. H.
Liu
, and
B.
Wei
, “
Singularity subtraction for evaluation of Green's functions for multilayer media
,”
IEEE Trans. Microwave Theory Tech.
54
(
1
),
216
225
(
2006
).
14.
M.
Yuan
,
T. K.
Sarkar
, and
M.
Salazar-Palma
, “
A direct discrete complex image method from the closed-form Green's functions in multilayered media
,”
IEEE Trans. Microwave Theory Tech.
54
(
3
),
1025
1032
(
2006
).
15.
T.
Richards
and
K.
Attenborough
, “
Accurate FFT-based Hankel transforms for predictions of outdoor sound propagation
,”
J. Sound Vib.
109
(
1
),
157
167
(
1986
).
16.
K. M.
Li
, “
Predicted attenuation of sound in a rigid-porous ground from an airborne source
,”
J. Acoust. Soc. Am.
123
(
3
),
1352
1363
(
2008
).
17.
F.
Ling
and
J.-M.
Jin
, “
Discrete complex image method for Green's functions of general multilayer media
,”
IEEE Microwave Guided Wave Lett.
10
(
10
),
400
402
(
2000
).
18.
K. M.
Li
and
S.
Liu
, “
Fast asymptotic solutions for sound fields above and below a rigid porous ground
,”
J. Acoust. Soc. Am.
130
(
3
),
1103
1114
(
2011
).
19.
L. M.
Brekhovskikh
and
O. A.
Godin
,
Acoustics of Layered Media II
(
Springer
,
Berlin
,
1999
), Vol.
10
, pp.
1
32
.
20.
S.
Liu
and
K. M.
Li
, “
Efficient computation of the sound fields above a layered porous ground
,”
J. Acoust. Soc. Am.
131
(
6
),
4389
4398
(
2012
).
21.
K. M.
Li
and
S.
Liu
, “
Propagation of sound from a monopole source above an impedance-backed porous layer
,”
J. Acoust. Soc. Am.
131
(
6
),
4376
4388
(
2012
).
22.
K.
Li
and
H.
Tao
, “
A modified saddle point method for predicting sound fields above a non-locally reacting porous medium
,”
Wave Motion
51
(
2
),
222
239
(
2014
).
23.
Y.
Wang
,
K. M.
Li
,
D.
Dragna
, and
P.
Blanc-Benon
, “
On the sound field from a source moving above non-locally reacting grounds
,”
J. Sound Vib.
464
,
114975
(
2020
).
24.
K. M.
Li
and
Y.
Wang
, “
On the three-dimensional sound fields from a moving monopole source above a non-locally reacting ground
,”
J. Acoust. Soc. Am.
147
(
4
),
2581
2596
(
2020
).
25.
Y. L.
Chow
,
J. J.
Yang
,
D. G.
Fang
, and
G. E.
Howard
, “
A closed-form spatial Green's function for the thick microstrip substrate
,”
IEEE Trans. Microwave Theory Tech.
39
(
3
),
588
592
(
1991
).
26.
T. K.
Sarkar
and
O.
Pereira
, “
Using the matrix pencil method to estimate the parameters of a sum of complex exponentials
,”
IEEE Antennas Propag. Mag.
37
(
1
),
48
55
(
1995
).
27.
Y. L.
Li
and
M. J.
White
, “
Near-field computation for sound propagation above ground—Using complex image theory
,”
J. Acoust. Soc. Am.
99
(
2
),
755
760
(
1996
).
28.
M.
Ochmann
, “
The complex equivalent source method for sound propagation over an impedance plane
,”
J. Acoust. Soc. Am.
116
(
6
),
3304
3311
(
2004
).
29.
R.
Piscoya
and
M.
Ochmann
, “
Acoustical Green's function and boundary element techniques for 3D half-space problems
,”
J. Comput. Acoust.
25
(
04
),
1730001
(
2017
).
30.
R.
Piscoya
and
M.
Ochmann
, “
Implementation of different formulas for special Green's functions in a 3D half-space bem
,”
J. Theor. Comput. Acoust.
27
(
01
),
1850051
(
2019
).
31.
C.-J.
Zheng
,
W.-Y.
Liu
,
Y.-B.
Zhang
,
C.-X.
Bi
,
H.-F.
Gao
, and
H.-B.
Chen
, “
Simulation of sound propagation over an infinite impedance plane by using a fast multipole BEM
,”
J. Theor. Comput. Acoust.
28
(
2
),
2050020
(
2020
).
32.
E.
Brandão
and
A.
Lenzi
, “
An iterative approach based on prony's method to calculate the surface impedance of acoustic materials measured in situ
,” in
SAE Technical Paper
,
SAE International
(
2010
).
33.
M.
Ochmann
and
R.
Piscoya
,
Theory and Application of Acoustic Sources Using Complex Analysis
(
Springer
,
Singapore
,
2021
), Vol.
1
.
34.
J. A.
Fawcett
, “
Complex-image approximations to the half-space acousto-elastic Green's function
,”
J. Acoust. Soc. Am.
108
(
6
),
2791
2795
(
2000
).
35.
J. A.
Fawcett
, “
A method of images for a penetrable acoustic waveguide
,”
J. Acoust. Soc. Am.
113
(
1
),
194
204
(
2003
).
36.
R.
Prony
, “
Essai exp érimental et analytique sur les lois de la dilatabilité des fluids élastiques, et sur celles de la force expansive de la vapeur de l'eau et de la vapeur de l'alkool, à différentes temperatures
,”
J. L'École Polytech.
1
,
24
76
(
1795
).
37.
T.
Sauer
, “
Prony's method in several variables
,”
Numerische Mathematik
136
,
411
438
(
2017
).
38.
Y.
Hua
and
T. K.
Sarkar
, “
Matrix pencil method for estimating parameters of exponentially damped/undamped sinusoids in noise
,”
IEEE Trans. Acoust. Speech Signal Process.
38
(
5
),
814
824
(
1990
).
39.
R.
Lanoye
,
G.
Vermeir
,
W.
Lauriks
,
R.
Kruse
, and
V.
Mellert
, “
Measuring the free field acoustic impedance and absorption coefficient of sound absorbing materials with a combined particle velocity-pressure sensor
,”
J. Acoust. Soc. Am.
119
(
5
),
2826
2831
(
2006
).
40.
M.
Müller-Trapet
,
P.
Dietrich
,
M.
Aretz
,
J.
van Gemmeren
, and
M.
Vorländer
, “
On the in situ impedance measurement with pu-probes—Simulation of the measurement setup
,”
J. Acoust. Soc. Am.
134
(
2
),
1082
1089
(
2013
).
41.
M.
Li
,
W.
van Keulen
,
E.
Tijs
,
M.
van de Ven
, and
A.
Molenaar
, “
Sound absorption measurement of road surface with in situ technology
,”
Appl. Acoust.
88
,
12
21
(
2015
).
42.
K. A.
Michalski
and
J. R.
Mosig
, “
Multilayered media Green's functions in integral equation formulations
,”
IEEE Trans. Antennas Propag.
45
(
3
),
508
519
(
1997
).
43.
A.
Sommerfeld
, “
Problems of radio
,” in
Partial Differential Equations in Physics
, edited by
A.
Sommerfeld
(
Academic Press
,
New York
,
1949
), Chap. VI, pp.
236
289
.
44.
G.
Taraldsen
, “
The complex image method
,”
Wave Motion
43
(
1
),
91
97
(
2005
).
45.
Y.
Hua
and
T. K.
Sarkar
, “
Generalized pencil-of-function method for extracting poles of an EM system from its transient response
,”
IEEE Trans. Antennas Propag.
37
(
2
),
229
234
(
1989
).
46.
Y.
Champoux
and
J.-F.
Allard
, “
Dynamic tortuosity and bulk modulus in air-saturated porous media
,”
J. Appl. Phys.
70
(
4
),
1975
1979
(
1991
).
47.
D.
Lafarge
,
P.
Lemarinier
,
J. F.
Allard
, and
V.
Tarnow
, “
Dynamic compressibility of air in porous structures at audible frequencies
,”
J. Acoust. Soc. Am.
102
(
4
),
1995
2006
(
1997
).
48.
J. F.
Allard
and
N.
Atalla
, “
Biot theory of sound propagation in porous materials having an elastic frame
,” in
Propagation of Sound in Porous Media
(
Wiley
,
New York
,
2009
), Chap. 6, pp.
111
135
.
49.
K.
Hirosawa
,
K.
Takashima
,
H.
Nakagawa
,
M.
Kon
,
A.
Yamamoto
, and
W.
Lauriks
, “
Comparison of three measurement techniques for the normal absorption coefficient of sound absorbing materials in the free field
,”
J. Acoust. Soc. Am.
126
(
6
),
3020
3027
(
2009
).
50.
E.
Brandão
,
A.
Lenzi
, and
J.
Cordioli
, “
Estimation and minimization of errors caused by sample size effect in the measurement of the normal absorption coefficient of a locally reactive surface
,”
Appl. Acoust.
73
(
6
),
543
556
(
2012
).