An optimal local active noise control method based on stochastic ﬁnite element models (cid:73)

A new method is presented to obtain a local active noise control that is optimal in stochastic environment. The method uses numerical acoustical modeling that is performed in the frequency domain by using a sequence of ﬁnite element discretizations of the Helmholtz equation. The stochasticity of domain geometry and primary noise source is considered. Reference signals from an array of microphones are mapped to secondary loudspeakers, by an oﬀ-line optimized linear mapping. The frequency dependent linear mapping is optimized to minimize the expected value of error in a quiet zone, which is approximated by the numerical model and can be interpreted as a stochastic virtual microphone. A least squares formulation leads to a quadratic optimization problem. The presented active noise control method gives robust and eﬃcient noise attenuation, which is demonstrated by a numerical study in a passenger car cabin. The numerical results demonstrate that a signiﬁcant, stable local noise attenuation of 20–32 dB can be obtained at lower frequencies ( < 500 Hz) by two microphones, and 8–36 dB attenuation at frequencies up to 1000 Hz, when 8 microphones are used.


Introduction
Active noise control (ANC) [1][2][3][4] is a well-established noise control technique nowadays. ANC is based on artificially produced secondary, opposite-phase anti-sound, and it has proven to be effective in low frequency noise attenuation. Passive methods such as absorbing and insulating acoustic materials are commonly used to attenuate high frequency noise, but in many applications they tend to be less effective and too heavy and bulky, if used in low frequency noise attenuation. Thus, ANC methods are of special interest in applications such as passenger cars and aircrafts, where low frequency noise is present and excessive weight has to be avoided.
The basic idea of ANC goes back to Paul Lueg's patent in 1936 [5]. The topic did not receive much attention until the rapid development of digital technology in 1980's, which enabled the research and development of ANC controllers and necessary circuitry. Since then, ANC techniques have been of a considerable research interest. There have been plenty of ANC studies in aircraft [6][7][8][9] and car cabin applications [10][11][12][13][14][15][16][17][18][19]. Aircraft propeller and car engine are typical low frequency noise sources that are considered relatively easy to control by ANC techniques, as engine rotation can be used as a reference signal to the system that controls secondary source actuators. Nowadays, there are commercial ANC systems available from several manufacturers for such applications.
The low frequency sources in cabin cavity are mainly due to structural vibration from engine and tires [20]. Especially structure-borne noise is low frequency: the mechanical vibratory noise from tires is mainly below 1 kHz and main frequencies of engine noise are below 500 Hz. Furthermore, low frequency noise is often amplified by resonances of the cabin structure. Road noise, that is generated by the interaction between wheel and road surface, is another significant noise source in cars, but it has been considered more difficult for ANC methods, mainly due to lack of a proper reference signal. ANC methods for road noise cancelation have been studied in, e.g. [17][18][19].
Numerical tools are nowadays popular in acoustical investigations and they have been used to study ANC methods. However, it is not easy to assess the effectiveness of ANC in a complicated three-dimensional domain.
Accurate modeling of acoustics can be a formidable task [21], particularly due to complicated noise sources and boundary impedances. Finite/boundary element method (FEM/BEM) can be used to simulate acoustical fields in ANC systems [22][23][24][25], and these methods have also been combined with a optimization methods to optimize for example sensor/actuator locations and the acoustic source strengths [26][27][28][29][30][31].
In typical ANC applications, there are unpredictable changes that occur in the acoustical space. In a passenger car, for example, driver, passengers, and parts of the machinery move. Stochastic domains can be conveniently used to model such random geometrical changes. There are plenty of publications on partial differential equations (PDEs) with stochastic coefficients, but PDEs in stochastic domains have had much less attention. In [32], a mapping from random to fixed domain is used to transform the problem to be one with stochastic coefficients. Fictitious domain approach is used in [33][34][35] and an extended finite element method is employed in [36] to treat stochastic domains.
The expected value of sound pressure can be conveniently computed in a stochastic domain by integrating numerically the product of sound pressure and probability density function over sampled set of possible domains. Thus, a functional of stochastic PDE solution can be computed directly, without any approximation of the solution to the actual stochastic PDE. This is a non-intrusive approach which allows using a solution method for non-stochastic problems without any modification. In [37], this approach was used to develop a method to study the performance of a local noise control that is optimal in a stochastic domain and in [31], the method was further used to find optimal secondary source locations for such local ANC system. However, the method could be used only for performance assessment purposes. A practical implementation of the method in a ANC system was not possible, particularly due to the unrealistic assumption of constant time-harmonic noise source.
In this paper, we propose a novel ANC method for enclosed cavities. The method is based on [37], but the optimization of secondary source signals is now reformulated such that an arbitrary number of reference microphones is used to adapt optimal ANC to prevailing acoustic field. This means that the system remains optimal even if changes in phase and amplitude of the noise occur. Moreover, a stochastic model of the noise source is considered, which allows location and properties of the noise source to vary to some extent. This should also make the model less sensitive to inaccuracy in noise source modeling. Most ANC systems employ error microphones to adapt the control system to the prevailing acoustical field. The presented method does not use actual error microphones, instead we could say that it employs purely virtual error microphones based on the acoustical model. Due to the stochastic geometry in the acoustical model, the locations of these virtual error microphones are allowed to be stochastic. This enables an efficient, local active noise control which is robust in highly random environment.
The contribution of this paper is the novel method to obtain an optimal local ANC in stochastic environment. The acoustical model used in the numerical example of Section 5 is not realistic, although the boundary parameter values and geometries should be indicative, and thus the results should be as well. The acoustical model can be interchanged without changing the method that is described here. Accurate modeling of the acoustics is out of the scope of this paper. This paper is organized as follows. In Section 2, a mathematical model of sound propagation, the Helmholtz partial differential equation, and a numerical method to solve it are briefly presented. In Section 3, the active local noise control in a stochastic domain with a stochastic noise source is formulated as a quadratic optimization problem. In Section 4, the details are given on how to implement the ANC system in practice. In Section 5, the performance of optimal ANC in a three-dimensional car cabin is studied numerically. In Section 6, concluding remarks are given.

Acoustic model
The time harmonic sound propagation is modeled by the Helmholtz equation where ρ is the density of the material, and c is the speed of sound. The complex-valued sound pressure p defines the amplitude and the phase of the sound pressure. The sound pressure at time t is obtained by Re e −iωt p , where ω is the angular frequency of sound and i = √ −1. A partially absorbing wall material is described by the impedance boundary conditions where η is the absorption coefficient depending on the properties of the surface material and f is sound source term on part S ⊂ ∂Ω. The value η = 1 approximates a perfectly absorbing material and the value η = 0 approximates a sound-hard material (the Neumann boundary condition). The absorption model is a simplification which is not very accurate with large values of η. Note also that the absorption coefficient can be frequency dependent with the considered frequency domain formulation. An approximate solution for the partial differential equation Eq. (1) can be obtained using a finite element method [38][39][40].

The noise control problem
An acoustical model is considered in an enclosed stochastic domain Ω(r), where r is a stochastic variable conforming to known probability density function F r (r). The total pressure field p (x, r, s, γ) is a sum of the sound pressures due to stochastic noise source and n a anti-noise sources where p 0 (x, r, s) is the sound pressure originating from a stochastic primary noise source, T is a vector of sound pressures due to secondary sources and γ = (γ 1 , . . . , γ na ) T is a complex-valued vector that contains coefficients γ i that alter the amplitudes and phases of each secondary source. The noise source p 0 is a function of stochastic variable s that conforms to a known probability density function F s (s).
To attenuate noise in a control volume, we will find an optimal linear form between reference signal measurements and secondary source amplitudes. A complex-valued vectorm, which contains perturbed noise measurements, is defined asm where m = (m 1 , . . . , m nm ) T is measured sound pressure at n m measurement points and e = (e 1 , . . . , e nm ) T is a normally-distributed variable containing corresponding measurement errors. We assume the real and imaginary parts of e i are normally distributed with mean values µ = µ = 0 and standard deviations σ 2 and σ 2 . Let us also define n a × n m complex-valued matrix C such that Substituting this from Eq. (3) yields The linear mapping C is optimized such that noise is minimized in a control volume denoted by V C (r) ⊂ Ω(r). To shorten the presentation, let us denote the integral of g (e) times the normal probability density function f n x; µ, σ 2 with the mean value µ and the standard deviation σ 2 as Let us now denote the ith column vector of the matrix C as Our noise measure is defined as where overline symbol· denotes element-wise complex conjugate, superscript · H denotes the Hermitian conjugate. Let us note thatˆ∞ for all normal distribution probability density functions f n x; µ, σ 2 . We obtain now where δ ij is the Kronecker delta function. The expected value of the noise measure is given by An approximation of the integral (13) is chosen as the objective function for optimization. It is given by the numerical quadrature where the pairs r j , w r j give the quadrature points and weights for r and the pairs (s k , w s k ) give the quadrature points and weights for s. Our intermediate optimization problem reads In order to give the objective function in a convenient matrix form, the following notations are introduced: where (12)). Furthermore, we define the following matrix and vectors where c i is a vector corresponding to the ith row of the matrix C. By expanding terms and by using the notations in Eqs. (16) and (17), the objective function in Eq. (14) can be expressed in a form In the case when there are no constraints, the linear mapping C is given by the optimality condition ∇ c J = 0. This leads to the system of linear equationsĀ which has the solution A real-valued formulation for Eq. (18) is given in Appendix A. More specific details for the noise control formulation are given in [41].  Figure 1: A block diagram illustrating the presented ANC system. It is assumed here, that the sensor microphones are located such that feedback is not significant compared to the primary noise.

Implementation of active noise control
To implement the presented ANC system in a practical environment, it is necessary first to perform off-line computations to produce the matrix C for each frequency. In Fig. 1, the ANC system is illustrated in a block diagram. The following steps are performed for each considered frequency.
1. Compute the sound fields in the control volume V c due to the primary noise source p 0 and the secondary noise sources, vector p, using the finite element model. If the control volume V c is reduced to a set of points (i.e. virtual error microphones), the reciprocity principle can be used, see Appendix B. 2. Compute the sound field due to the noise source p 0 on the reference microphone locations, i.e. vector m, using the finite element model. The reciprocity principle can be used. 3. Assemble the matrix A and vector b in Eq. (17) and solve the system of linear equations, Eq. (19). The resulting matrix C is stored to a database.
The matrix database is stored in the ANC controller. The following procedure is performed in the controller circuit during ANC operation in frequency domain.
1. Complex-valued pressure amplitudes (i.e. the phase and amplitude of the signal) in reference microphones, the vector m (see Eq. (6)), are obtained for each considered frequency. 2. For each considered frequency, the complex-valued pressure amplitude of each secondary source, the vector γ in Eq. (3), is obtained by multiplying m with the pre-computed matrix C, as in Eq. (5) with e = 0.

Noise control in a car cabin
As a numerical example, we will study the presented local active noise control in BMW 330i car interior (see Fig.  2a). The height of studied car interior is 1.3 m, width 1.8 m and length is 3.1 m. It is assumed here, that the driver is alone in otherwise empty car. Thus, the cabin interior excluding the driver is the computation domain Ω (r) for acoustical sound pressure field.
The objective of the noise control is to minimize noise in driver's ears. Thus, we define control volume V C (see Eq. (8)) as a set of two points where x e l (r) and x er (r) are the coordinates of the left and right ears, respectively; these points can be interpreted as virtual error microphones. The noise measure equation (8) Various driver's properties alter sound propagation inside the car cabin. Especially posture and head position affect the sound pressure that is experienced in ears and these parameters lead to the stochasticity of the computation domain. As in [37], the domain stochasticity variable r = (r 1 , r 2, r 3 ) T consists of three parameters: r 1 is driver's sideways bending angle, r 2 is forward bending angle, and r 3 is head rotation angle to left/right (see Fig. 2b); all these variables are degrees by dimension. The probability density function F r for the variable r is given by a piecewise trilinear function defined by the nodal values at 5 3 = 125 points, and elsewhere by trilinear interpolation, thus F r (r) is composed of 2 3 = 8 trilinear functions. The integral in Eq. (13) is approximated by a three-dimensional generalization of the trapezoidal quadrature rule.
The primary noise source is modeled as a vibrating quadrangular surface behind the leg room, which is a simplification resembling the real primary noise source (see Fig. 2a). The amplitude f (see Eq. (2)) of the quadrangular noise source is given by a bilinear function specified by the corner values f (x i , y i ) = s i , i = 1 . . . 4, where s i is the amplitude coefficient for the ith corner of the quadrangular noise source surface located at (x i , y i ) defining the noise source stochasticity variable s = (s 1 , s 2 , s 3 , s 4 ) T .
The probability density function F s for the variable s is given by a piecewise quadrilinear function defined by the nodal values at 3 4 = 81 points, and elsewhere by quadrilinear interpolation, thus F s (s) is composed of 2 4 = 16 quadrilinear functions. The integral in Eq. (13) is approximated by a four-dimensional generalization of the trapezoidal quadrature rule.

Numerical study 5.2.1. Model definition
For the stochastic variable r defining the domain Ω, we will use sampling r ∈ R = R 1 × R 2 × R 3 , where On the boundaries of hyper-quadrangles defined by sets R and S, probability density functions are set to zero, i.e. F r (r) = 0 and F s (s) = 0, from which it follows that we do not need to create samples for these sampling points. In addition, the probability function values at the center point of the hyper-quadrangles are set so that the integral of the probability density function over the probability space is one. The locations of reference signal measurement points and secondary sources for ANC are presented in Fig. 2a and they are labeled as follows: L1...L8 for left side measurement points, R1...R8 for right side measurement points and A1...A6 for six circle-shaped secondary sources. The absorbency coefficients (see Eq. (2)) were set per material as listed in Table 1.
To solve the Helmholtz equation in Eq. (1) with the finite element method, a collection of meshes consisting of tetrahedra and triangles were generated with Ansys ICEM CFD. Standard linear Lagrangian basis functions were used in finite elements. The solutions were computed with a damped Helmholtz preconditioner technique described in [39,40]. Each mesh corresponds to a different driver posture and they were generated so that there are at least 10 nodes per wavelength at the highest studied frequency (1000 Hz). The total number of meshes needed for domain sampling is 5 3 = 125 which is the number of parameter combinations (r 1 , r 2 , r 3 ).
The study was performed in the frequency range 10-1000 Hz with 10 Hz steps. This means that 100 frequencies are sampled. By employing the reciprocity principle (see Appendix B) a sound source is placed in each of 16 measurement points and 2 ears. The acoustic model is solved for all 125 sampled driver's postures for each 18 sources. Thus, a discrete Helmholtz equation is solved 125 × 100 × 18 = 225000 times.

Expected noise attenuation and measurement error sensitivity
There are several error sources in the presented ANC method that are due to inaccuracy of microphones and loudspeakers and the model itself. To obtain good noise control, the ANC method should not be too sensitive to such errors.
The sensitivity of the performance to measurement error in the model is adjusted by the standard deviation parameters σ 2 i (see Eq. (16)). These parameters are set proportional to the absolute value of the measured value |m i | as σ 2 i = ε |m i |. This is based on the assumption that the measurement error is proportional to the amplitude of measurement value. It was observed by numerical experiments that the error can be underestimated, thus the coefficient ε was set small, ε = 10 −4 . Due to the quadratic formulation, even small standard deviation values σ 2 i decrease the error sensitivity significantly. If the values σ 2 i were not underestimated, the results would indeed have very low error sensitivity, but the expected noise attenuation would be poor.
The influence of the measurement error e i is studied by considering the measurement error present in the real and imaginary parts of each measurement m i . The error in this study is due to coefficients k ∼ N 0, σ 2 and k ∼ N 0, σ 2 that are used to multiply the real and imaginary parts of the error, i.e. the perturbed measurement ism where e i = k m + i k m . The noise attenuation is sampled with 100 different measurement errors. In Fig. 3 23) and (24). Fig. 3 shows a typical noise attenuation obtained by ANC when the data from 2 and 4 measurement points are used. The standard deviation and min/max curves reveal that the noise control is robust with respect to stochastic features: the worst sampled values produce still good, around 20 dB reduced noise attenuation at low frequencies.
It can be seen from perturbed noise attenuation curves that 4% measurement error can cause 10-30 dB deterioration of the noise control at low frequencies.
In Fig. 4, the performance of ANC with different number of measurement points is studied. In Fig. 4a it is seen that increasing the number of measurement points enhances the noise control. Fig. 4b shows the same results as Fig. 4a, but with 4% measurement error. With errors in measurements, the noise control is slightly less effective, which is expected. The enhancement gained by increased number of measurement points is also reduced. In Fig. 5, each line of Fig. 4a are plotted as a separate figure together with lines that correspond to 1%, 2% and 4% measurement errors. In addition to these lines, the corresponding noise attenuation with the largest sampled measurement error is plotted. The figure clearly shows that the error sensitivity is very reasonable, and even with the large (4%) measurement error, the noise attenuation is very good. It can be also concluded that a large number of measurements is not required to get very good attenuation, but with increasing the number, it is possible to further enhance the noise control performance.
The influence of the measurement error control variables σ i can be seen by comparing Fig. 6 to Figs. 5a and 5b, when the optimal matrix is computed without measurement error term (σ i = 0 in Eq. (16)). When the number of measurement points used is at least 4, the error sensitivity causes poor results, even with small 1% measurement error.

Conclusions
A local active noise control (ANC) method that is optimal in stochastic environment was presented. The ANC method uses modeling data that is obtained by acoustical simulations using the finite element method in stochastic computational domain. The noise control obtained by the ANC method is not sensitive to measurement errors or changes in noise source or geometry.
The efficiency of the ANC method was demonstrated numerically in a car interior with a driver in varying postures. Sound attenuation of 50-60 dB was obtained at lowest frequency (10 Hz). The error sensitivity of the ANC method was studied and it was found that with 1-4% error, 30-46 dB expected attenuation was obtained at low frequencies, and expected attenuation was more than 20 dB up to 500 Hz. It can be concluded, that the numerical results are very good. The results could be further enhanced by choosing optimal locations for sensor and secondary source locations.
The acoustical finite element model (especially boundary impedances, noise source model and geometries) used in the numerical example is demonstrative and is not accurate nor realistic. Thus, the results obtained in the numerical example are not reproducible in real experiments. However, the model should be suggestive and the results can be used to estimate the efficiency of the ANC method. To implement the ANC method in real environment, more accurate acoustical model is needed. The acoustical model used in this paper can be interchanged without changing the method that is described here.