Introducing passive acoustic filter in acoustic based condition monitoring: Motor bike piston-bore fault identification

Mechanical Systems and Signal Processing 70-71 (2016) 932–946

Contents lists available at ScienceDirect

Mechanical Systems and Signal Processing journal homepage: www.elsevier.com/locate/ymssp

Introducing passive acoustic filter in acoustic based condition monitoring: Motor bike piston-bore fault identification D.P. Jena a,n, S.N. Panigrahi b a b

Department of Industrial Design, National Institute of Technology Rourkela, Orissa, 769008, India School of Mechanical Sciences, Indian Institute of Technology Bhubaneswar, Bhubaneswar 751 013, India

a r t i c l e in f o

abstract

Article history: Received 25 October 2014 Received in revised form 30 May 2015 Accepted 27 September 2015 Available online 23 October 2015

Requirement of designing a sophisticated digital band-pass filter in acoustic based condition monitoring has been eliminated by introducing a passive acoustic filter in the present work. So far, no one has attempted to explore the possibility of implementing passive acoustic filters in acoustic based condition monitoring as a pre-conditioner. In order to enhance the acoustic based condition monitoring, a passive acoustic band-pass filter has been designed and deployed. Towards achieving an efficient band-pass acoustic filter, a generalized design methodology has been proposed to design and optimize the desired acoustic filter using multiple filter components in series. An appropriate objective function has been identified for genetic algorithm (GA) based optimization technique with multiple design constraints. In addition, the sturdiness of the proposed method has been demonstrated in designing a band-pass filter by using an n-branch Quincke tube, a high pass filter and multiple Helmholtz resonators. The performance of the designed acoustic band-pass filter has been shown by investigating the piston-bore defect of a motor-bike using engine noise signature. On the introducing a passive acoustic filter in acoustic based condition monitoring reveals the enhancement in machine learning based fault identification practice significantly. This is also a first attempt of its own kind. & 2015 Elsevier Ltd. All rights reserved.

Keywords: Passive Acoustic Filter Condition monitoring Quincke tube Helmholtz resonator Genetic Algorithm (GA) Support Vector Machine (SVM)

1. Introduction Every system is unique in itself and requires specific formulation of signal processing methodology pertinent to the system or faults under investigation [1–3]. So far, the piston-bore defect identification based on engine noise signature has not been addressed profoundly. Although various methodologies have been proposed for acoustic based condition monitoring of similar applications, suitability of such techniques are known to be highly reliant on the type of application and the quantity and quality of parameters that are monitored [4–9]. From the literature, it is also worth noticing that when expert systems are developed based on real-time data, robustness of such systems depend crucially on the nature of the signal signatures, signal parameters used to train the system and the total number of training samples [5,10–12]. The performance of supervised learning technique highly hinges on the number of parameters and numbers of training samples [2,12]. However, use of higher numbers of parameters need high computational time which is not desired. In order to achieve adequate performance of expert system using lower number of parameters an optimized filtering is required at

n

Corresponding author at: Department of Industrial Design, National Institute of Technology Rourkela, Orissa, 769008, India. E-mail address: [email protected] (D.P. Jena).

http://dx.doi.org/10.1016/j.ymssp.2015.09.039 0888-3270/& 2015 Elsevier Ltd. All rights reserved.

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

933

the time of signal acquisition. So, a passive acoustic filter has been introduced before the acoustic sensor. The corresponding enhancement in the proposed acoustic based condition monitoring has been demonstrated. In implementing real time condition monitoring systems, so far, no one has attempted to explore the possibility of implementing passive acoustic filters in acoustic based condition monitoring as a preconditioner. In an intention of introducing passive acoustic filters as preconditioner in the acoustic based condition monitoring, literature covering analysis and design of passive filters has been reviewed and summarized below. One can find detailed analogies of filtering properties of various passive acoustic filters to those of digital filters [13]. In general, parameters affecting the acoustic impedance dictate the acoustic properties of a passive filter. In recent years, numerous analytical solutions have been reported for estimating acoustic transmission loss (TL) of various types and shapes of filters [14–16]. The use of transfer matrix technique facilitates designers to design a wide range of filters with minimal acceptable error. If the filter is not a multiply connected one, on multiplying transfer matrix of each sub component the resultant transfer matrix for the whole desired filter can be achieved. For estimating acoustic TL of such a filter, it is required to have transfer matrix for each of its sub components [14,16]. After this present section of introduction, Section 2 explains the design of acoustic based expert system and its application in piston-bore fault detection. The concept of designing a passive acoustic band-pass filter and introducing in acoustic based condition monitoring has been demonstrated in Section 3. The last section summarizes some important observations of the present experimental research.

2. Development of the acoustic based expert system In the present work the most common engine defect, commonly branded as the piston-bore defect, which is also known to have direct impact on the quality of the exhaust, has been addressed. A sample of a defective piston-bore pair of a vehicle is shown in Fig. 1(a) and (b). This work has been inspired by the near super-human ability of auto mechanics in identifying this defect just by carefully listening to the sound coming from the vehicle. Therefore, in an attempt to mimic their working style, an expert system has been developed that broadly follows the stages of appropriate signal filtering, statistical learning and feature based classification. For automatic fault detection, development of a robust supervised learning system is desired. In the proposed acoustic based expert system, a supervised statistical learning mechanisms has been established. Out of many methods, artificial neural network (ANN) and support vector machine (SVM) have been widely adopted to tackle similar types of systems. With an intention to develop an expert system for automatic piston-bore fault identification using engine noise signal, an expert system based on ANN and SVM have been developed separately, and have been investigated on experimental data. The performance of six of the most preferred conventional parameters that are used in statistical learning systems to tackle similar problems, has been demonstrated. A variant of artificial neural network (ANN) and support vector machine (SVM) has been implemented to set up a robust expert system suitable for natural workshop environment. These stages of implementation are explained in detail in subsequent sections. 2.1. Multi-layered feed forward back propagation neural network Artificial neural network (ANN) is considered to have enough potential to build a fairly acceptable system model even with insufficient knowledge about the system. ANN gains knowledge through an adaptive learning process by adjusting its weights by observing a series of inputs and the corresponding outputs. This process is usually called the training of an ANN. There are various types of neural network models based on different kernel functions. Feed-forward back-propagation neural network (FBNN) structure is a widely used one in machine fault diagnosis and condition monitoring purposes [17,18,11]. In this study, a multi-layer (with 10 hidden layers) feed-forward back-propagation neural network (FBNN) is designed and used. An analogous network, with six inputs, one hidden layer and two outputs, is shown in Fig. 2(a). The general rule of selecting number of hidden layers depends on the number of inputs and outputs, i.e., it is taken to be equal to two third of the total number of inputs and outputs. However, the quality of approximation hinges highly on the

s

Fig. 1. Experiment; (a) sample faulty piston of DV-II, (b) sample faulty bore of DV-II, and (c) noise data acquisition using Bruel Kjaer 2270 .

934

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

Fig. 2. Schematic of machine learning techniques; (a) artificial neural network(ANN), and (b) support vector machine (SVM).

number of training samples. In the present case, the analysis starts with six inputs and one output. As per the mentioned thumb rule, the number of hidden layers should be 5. From experimental validation, it is also observed that change of 9–10hidden layers reduces the training error by 2%. Further increasing the number of hidden layers did not increase the efficiency significantly. Hence it is maintained at 10 for the effective analysis. In feed-forward back-propagation neural network (FBNN) structure, each propagation involves two steps, i.e., forward propagation and backward propagation. The “feed-forward” step involves in processing the developed neural network and the patterns. The “back-propagation” step involves in gaining the intelligence. The schematic of the implemented algorithm for the present situation is expressed in the following pseudo code. Pseudo-code: Multi Layer Feed Forward Back Propagation Neural Network Random initialization of the weights in the network: yi ¼ w  xi þ b Do For each output neuron q, compute feed forward Error, e ¼ di  ðb þ wxi Þ, where di is the desired output using MSE, i.e., E ¼ ð1=2NÞ

PN i

e2i

If fE r E DesiredError J Loop Iteration ¼1000 ¼ ¼True Stop Else For each output neuron q, the adapted weight and bias become P P 2P P P ðx  x~ Þðdi  d~ Þ x d  x x d i i i i i i  i i i wi ¼ P ; bi ¼ i i P 2 i

ðxi  x~ Þ

N

i

ðxi  x~ Þ2

For each non-output neuron p, the estimated back propagation error from layer Ln  1 to layer L1 is computed as: ! sl þ 1 sP l1 ðl þ 1Þ i;l i;l P ξðlÞ wlpq ξq where εðlÞ wlpq 1 εlp 1 q ¼ f LogSigmoid p ¼ ð1  εp Þεp q¼1

p¼0

ðl þ 1Þ

Update the weights in the network by: wlpq ’wlpq þ αalp δq Return (Do)

Using standard six signal statistical parameters, the inadequacy of FBNN in identifying piston-bore fault has been observed. However, such systems are highly reliant on number of training samples. 2.2. Multi-class SVM with Maximal Violation Support Vector In order to achieve further accuracy in fault identification, another statistically learning technique such as multi-class SVM with maximum violating support vector have been designed and implemented. SVM has advantages in solving complex practical problems involving non-linearity, higher dimensions, local minima and/or small training sample [12,2,19,20,1]. The purpose of classification is to provide a mapping from the feature domain to the condition domain. The implementation of SVM for fault pattern classification using sets of statistical parameters follows a predefined routine. First, a mapping of the statistical parameters (input vectors) onto one or multiple feature spaces is carried out. This mapping may be either linear or nonlinear. Afterwards it dictates the selection of a suitable kernel function. Then, within the feature spaces generated from the first step, SVM seeks an optimized partition. It constructs a hyper-surface that separates the classes significantly. SVM training looks for a globally optimized solution and avoids over-fitting. So it has the ability to deal with a large number of features. ^ Mathematically, a pattern x is given by a class yðxÞ ¼ 7 1 with a feature vector and a linear discriminate function can be expressed as [21,22] ^ yðxÞ ¼ wT ΦðxÞ þ b:

ð1Þ

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

935

The parameter w and b are obtained from training the data sets of ðyi ; xi Þ1 r i r N . By minimizing w2 =2, the SVM technique tries to find a unique separating hyper-plane (see Fig. 2(b)). Then the optimization necessity can be modeled as min pðw; bÞ ¼ w2 =2 subject to yi ðwT Φðxi Þ þ bÞ Z1

ð2Þ

This can be turned into a dual problem and can be simplified with a Lagrangian multiplier ðαÞ and can be expressed as (soft margin support vectors) max

DðαÞ ¼

n X

n   1 X y α y α Φðxi ÞT Φ xj ; 2 i;j ¼ 1 i i j j X and yα ¼0 i i i

αi 

i¼1

subject to

αi Z0

ð3Þ

The direction vector w of the optimal hyper-plane can be evaluated from the solution of α on this situation. Optimizing of only bn leads to an one-dimension solution and the corresponding decision function can be presented as 

 ^ yðxÞ ¼ wT x þb ¼



n X



yi αi Φðxi ÞT ΦðxÞ þ b

where

w ¼

i¼1

X i

yαi Φðxi Þ

ð4Þ

For a linearly non-separable case, SVM performs a nonlinear mapping of the input vector to a higher dimensional Hilbert space. The nature of the non-linear mapping is determined by the kernel function mentioned below [2,19] 8 Linear Kernel: Kðx; yÞ ¼ x:y > >    > > < Polynomial Kernel: Kðx; yÞ ¼ x:y þ 1 p ; p A 2; 3…N ð5Þ 2 2 > Radial basis Kernel: Kðx; yÞ ¼ e  jx  yj =2σ > > > : Sigmoid Kernel: Kðx; yÞ ¼ tanhν:x:y þ c The selection of kernel function depends on the quality of the input parameters, the type of classification problem under investigation and the desired accuracy of the classification results. Incorporating the Kernel function, the support vector decision function can be presented as  ^ yðxÞ ¼ wT x þb ¼

n X



yi αi Kðxi ; xÞ þb :

ð6Þ

i¼1

The solution can be achieved by a decomposition method or modified gradient projection technique. The decomposition method is opted in present situation. The decomposition method optimizes a subset of Lagrangian coefficients, i.e., αi ; i A β and excludes the remaining coefficients, i.e., αj ; j2 = β. The corresponding decomposition can be expressed as max Dðα0 Þ ¼

n X i¼1

α0i 

n   1 X y α0 y α0 K x i ; x j 2 i;j ¼ 1 i i j j

  2 = βα0i ¼ αi   A β 0 r α0 r ζ subject to i  P 0  i yi αi ¼ 0

ð7Þ

In order to achieve the optimal solution, the optimization of the variable α0 can be excluded and the simplified maximization model can be expressed as 0 1 X X 1XX ! max α0 @1 yi yj αj Ki;j A  y α y α0 K iAβ i 2i A β j A β i i j j i;j j2 =β X X subject to A β ; r α0i r ζ and yi α0i ¼  yj αj ð8Þ iAβ

j2 =β

In present problem, the implemented algorithm to achieve optimized SVM along with maximum violating support vector pair can be expressed as: Pseudo-code: SVM with Maximal Violation Support Vector Initialize the coefficients and gradients: αk -0, and, δk -1j 8 k A f1;2;3…Ng Do imaxðyi δi Þ -ijyi αi o βi and jminðyj δj Þ -ijyj αj 4 Aj If yi δi r yj δj Stop Else



βi  yi αi ; yi αi  Aj ; Ki;iyþi δKi j;j yj δ2Kj i;j 

min

-γ

δk  γ yk K i;k þ γ yk K j;k -δk j 8 K A f1…Ng

Update the coefficients and gradients: αj  yj γ -αj , and, αi þ yi γ -αi Return (Do)

936

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

Table 1 Technical specifications of the vehicle used in the experiment. Specification

Details

Engine type Ignition system Bore and stroke (mm) Clutch Drive

4 stroke single cylinder forced air-cooled Digital DC CDI 51 and 43 Pivoted clutch centrifugally operated Automatic continuously variable transmission (CVT) with torque sensor

Fig. 3. Experimentally collected noise signals and their PSD spectra; (a) noise signal from healthy vehicle (HV_Signal), (b) corresponding PSD spectrum, (c) noise signal from faulty vehicle (DV_Signala), (d) corresponding PSD spectrum, (e) noise signal from faulty vehicle (DV_Signalb), and (f) the corresponding PSD spectrum.

2.3. Experiment and data acquisition s

s

s

The present experiments were conducted with a motor bike ((TVS ), Scooty Pep  Plus ) in a company (TVS ) authorized auto-workshop. The setup including vehicle and the remaining equipment has been shown in Fig. 1(c). The specification of the vehicle is tabulated in Table 1. The experiments were carried out in the natural workshop environment with an aim to develop a system that works even in such a setting without requiring any laboratory kind of environment. The engine noise s signal has been captured using a standard hand held sound analyzer (Bruel&Kjaer 2270) at 48 kHz sampling rate. To understand the actual defect signature, no additional amplification or weighing filters were incorporated during data acquisition. As shown, the sound analyzer had been mounted on a tripod near the engine. In order to achieve a repeatable testing environment, it was ensured that the vehicle is on its main stand, engine ignition is switched on by the key and the rear wheel does not rotate during the experiment. Experiments were carried out with three different types of motor bikes of the same make and model as mentioned. With a view of comparing acoustic signal of defective and healthy vehicles, the engine noise signal from a defect free vehicle was captured for over 30 min and is referred to as “HV_Signal” hereafter. A 0.5 s sample of acoustic signal and its corresponding power spectral density (PSD) are shown in Fig. 3(a) and (b), respectively. Subsequently, two more used-vehicles (of the same make and model) were identified and predicted to have piston-bore defect. This prediction had been made by the automechanic just by hearing the engine noise signal from the vehicle. The engine noise signal from the first used-vehicle (distance travelled 34,074 km) was captured over a period of 30 min and is referred to as “DV_Signala” hereafter. A 0.5 s sample of acoustic signal and its corresponding PSD are shown in Fig. 3(c) and (d). Similarly, the engine noise signal from another used-vehicle (distance travelled 19734 Km) was captured for over 30 min and a 0.5 s sample along with its PSD are presented in Fig. 3(e) and (f), respectively. This signal is referred to as “DV_Signalb” hereafter. After the data is acquired, 350 test samples of duration 0.5 s from HV signal and 300 test samples of duration 0.5 s from DV_Signala and DV_Signalb, each, have been extracted. In total, 950 numbers of acoustic test samples of duration 0.5 s

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

937

Table 2 Mean and coefficients of variation (CV) (in %) of statistical parameters of engine noise signal samples. Parameter name (1–6)

Crest factor SNR Skewness Kurtosis Peak frequency (Hz) Peak frequency amplitude

HV_Signal

DV_Signala

DV_Signalb

Mean

CV %

Mean

CV %

Mean

CV %

4.33 0.82 0.04 2.56 85.59 0.06

1.37 1.63 290.13 9.07 22.45 22.85

5.06 0.8  0.02 2.94 74.09 0.05

4.73 1.32 810.56 8.12 64.98 27.11

5.41 0.77  0.12 3.54 42.53 0.04

0.37 1.32 133.48 7.76 72.38 14.38

Fig. 4. Testing of test samples using ANN trained by six parameters and 800 samples.

have been extracted and processed offline to establish and demonstrate the robustness and effectiveness of the proposed technique. After completion of the acoustic data acquisition, the vehicles have been dismantled by the auto-mechanic in order to reveal the actual extent and severity of the defect, as shown in Fig. 1(a).

2.4. Performance analysis of the expert system The frequency domain acoustic signals (i.e., HV_Signal, DV_Signala and DV_Signalb) indicate that, in presence of the defect, the frequency components below 300 Hz alter prominently. However, no particular frequency peak could explain the presence of fault. This may be due the facts that the faults are not identical but similar in nature. Predictably, one can notice that the vehicle acoustic signals are not stationary in nature over any particular period. In order to model the faulty system, six standard signal parameters such as crest factor (CF), SNR, skewness, kurtosis, peak frequency (Hz), peak frequency amplitude have been identified and have been estimated. The averages of such parameters have been tabulated in Table 2. As discussed earlier, with an intention to develop and establish a robust expert system, subsequently, a multi-layer (with 10 hidden layers) feed-forward back-propagation neural network (FBNN) is designed and investigated as discussed earlier. Out of the 950 acoustic signal samples, which include samples both from the defective and the non-defective vehicles, 800 random samples (300 test samples from HV_Signal and 250 test samples from DV_Signala and DV_Signalb, each) have been chosen randomly to train the ANN model. The rest 150 unused samples (50 test samples from healthy vehicle (HV_Signal), and 50 test samples from DV_Signala and DV_Signalb, each) have been used for testing and validation of the trained model. This is a supervised learning method where six statistical parameters (see Table 2, for 800 samples) have been provided as inputs to the six neurons in order to train the FBNN. Corresponding desired outputs (or targets) have been mapped as “Healthy” or “Defective” vehicle types. The stopping criteria for FBNN training is based upon the minimum error gradient magnitude. In the present case it is less than 1e-5 (refer Fig. 4). After successful training, the 150 acoustic signal samples (those were kept aside initially) have then been tested for vehicle condition identification. The healthy and the defective vehicles can be observed to be classified with an average of 94.59% and 93.16% success rate, respectively. In order to get an overview of the success percentage, the success rate for rest unused samples have been computed. From experimental analysis it is observed that only few samples are having detection probability below 80%, and may be considered as failure in the detection. Thus, this configuration of FBNN can be said to have a 91% recognition rate and 9%

938

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

Table 3 Evaluation of signal statistical parameters with respect to kernel functions of SVM using least square optimization (in %). Parameters

Linear

Quadratic

Polynomial

RBF

MLP

Parametr1 Parametr2 Parametr3 Parametr4 Parametr5 Parametr6 Parameters1–6

91.89 81.39 66.68 74.64 74.64 69.75 91.22

91.89 82.74 66.58 81.50 68.56 72.14 92.67

91.89 81.39 67.15 81.39 73.96 72.14 92.09

91.89 82.74 63.77 82.74 79.37 73.39 92.09

91.89 59.93 64.92 81.50 92.20 72.14 86.59

Table 4 Evaluation of signal statistical parameters with respect to kernel functions of SVM using sequential minimal optimization (in %). Parameters

Linear

Quadratic

Polynomial

RBF

MLP

Parametr1 Parametr2 Parametr3 Parametr4 Parametr5 Parametr6 Parameters1–6

91.89 82.74 64.45 82.07 65.85 68.30 94.02

91.89 81.39 60.86 82.74 68.56 72.14 94.02

91.89 82.74 65.70 81.50 71.26 71.36 93.45

91.89 82.74 65.70 81.39 77.34 72.14 94.4

91.89 81.39 58.32 75.62 65.85 62.42 67.722

failure rate as defined below Recognition Rate ¼

Correct no: of classification samples  100 Total no: of testing samples

ð9Þ

The achieved rate of fault detection is not satisfactory from industrial perspective. This triggers us to identify and evaluate the other supervised learning mechanism such as SVM. Hence, the same six statistical parameters, that were used to train the FBNN, have been used to train the SVM model (see Table 2). To arrive at a suitable kernel function for the existing problem, classification strength of each statistical parameter has been evaluated independently. To start with, classification has been attempted using all the statistical parameters individually. This is done using the least square (LS) optimization method for all the standard kernel functions such as (a) polynomials up to third degree, (b) Gaussian radial basis function (RBF) and (c) multilayer perceptron (MLP). Then the trained SVM has been tested with rest 150 unused acoustic signal samples. The average of healthy and defective vehicle recognition rate has been calculated. The computed recognition rate using individual parameters and using their combination are shown in Table 3. In a similar manner, the SVM model has been tested using all the statistical parameters individually and then using their combination through the sequential minimal optimization (SMO) method for all the standard kernel functions as mentioned above. The computed average recognition rate using individual parameters and their combination are tabulated in Table 4. From Tables 3 and 4 it is clear that the maximum recognition rate of 94.4% has been achieved using RBF kernel function along with sequential minimal optimization method. It is a common practice that this particular engine defect, i.e., the piston-bore fault, needs to be fixed by disassembling the engine which is both time and money consuming. Therefore, a very high identification/classification accuracy is desired before proceeding with such effort. To achieve such accuracy levels, introduction of a few more parameters can be contemplated. However, the computational time will increase on increasing number of parameters. Therefore, alternate filtering, preferably at the time of signal acquisition, is desired to enhance the performance of the expert system based on standard six statistical parameters. In order to achieve this, the idea of introducing a passive acoustic band-pass filter has been investigated in subsequent section.

3. Performance improvement using Passive Acoustic Filter In order to reduce the computational overhead in the proposed methodology, the desire is to acquire the objective acoustic signal by applying a band-pass filter. From earlier sections, it has been noticed that the useful range of frequency in the acoustic signal lies between 50 Hz and 300 Hz. Therefore the idea of introducing a passive acoustic band-pass filter for such purpose has been tried out, as shown in Fig. 5. In present work, the objective is to design a passive filter which can perform closely compared to a sophisticated digital band-pass filter.

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

939

Fig. 5. Implementation of proposed acoustic filter in acoustic based condition monitoring.

3.1. Designing of Passive Acoustic Filter The original theory on passive acoustic filters and corresponding electrical analogy has been given by Stewart [23] followed by a detailed study by Olson [13]. Over the years, many researchers have presented various analytical solutions of various types of acoustic filters [14]. In the present work, the n-element passive acoustic filter refers to filters that can be formed by combining ‘n’ numbers of acoustic filters in series. In the process of optimizing the design of such n-element acoustic filters, the sound power transmission coefficient, at, of the filter has been used, which is defined as the ratio of the transmitted acoustic power to the incident acoustic power [24]. Mathematically, the sound power transmission coefficient can be derived using acoustic transmission loss of a filter as at ¼

W tr ¼ 10TL=10 ; W in

where the transmission loss for a passive acoustic filter can be computed using the final transfer matrix as "     # Y 1 0:5 T 11 þ T 12 =Y 1 þ ðY n T 21 Þ þðY n T 22 =Y 1 Þ TL ¼ 20log   2 Yn

ð10Þ

ð11Þ

where Y1 and Yn are the characteristic impedance of the inlet and outlet duct attached to acoustic filter (for the case at hand, Y 1 ¼ Y n ¼ ρc=π r 2 ; r is the radius of the duct and ρc is the characteristic acoustic impedance; ρ is the density of air and c is the speed of the sound) and Tij are the elements of the transfer matrix defined as " # " #" # pn  1 T 11 T 12 p1 ¼ ð12Þ vn  1 T 21 T 22 v1 Here, the pi and vi are sound pressure and particle velocity, respectively. The transfer matrix of standard acoustic filters [14,13,23] used in the present work for designing the n-element filter can be found in literature [25,26]. With an assumption that the transfer matrix of each passive acoustic filter component is known, the transfer matrix of an n-element passive filter can be written as TM1  n ¼ TM1  TM 2  …  TMn

ð13Þ

Traditionally, a microphone is placed near the machine/system to acquire data for acoustic based condition monitoring [27]. However, in the present work, in order to enhance the signal quality, insertion of a passive acoustic filter in the signal path has been proposed. One end of the acoustic filter has been proposed to face the noise source (preferably normal), so that the plane acoustic waves can travel through the filter. The other end of the acoustic filter has been proposed to have an anechoic termination. Such a design has been proposed to restrict the interference of filtered acoustic signal with the environmental noise. The microphone or the sound analyzer has been mounted on the tube (before the anechoic termination) on the outlet side of the filter. The schematic of the measurement has already been shown in Fig. 5. The corresponding electrical analogy of the acoustic filter along with inlet and outlet acoustic impedances has also been shown in same figure. Mathematically, the transfer matrix for above situation can be written as

1 Z Air ð14Þ  TM1  n TM ¼ 0 1 The impedance of the open air termination ZAir can be defined as [14]

ρc 1 Z Air ¼ 2 ðkr d Þ2 þ j0:61ðkr d Þ π rd 4 unflanged

ð15Þ

940

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

Fig. 6. Acoustic properties of standard acoustic filters; (a) low-pass filter, (b) high-pass filter, (c) band-stop filter and (d) Quincke tube.

where rd is the filter inlet radius. In the present design, rd is set to 0.05 m. As discussed earlier, the simplest of the filter components have been combined in series to achieve the desired band-pass filter. It is impossible to manually decide the geometry of each component in order to arrive at the desired performance of the filter. Therefore, to achieve this, a genetic algorithm based optimization technique has been adopted. 3.2. Transfer matrix of various types of acoustic filters In present work, the objective is to introduce a band-pass passive filter to enhance acoustic based condition monitoring. Conventionally, a band-pass filter can be constructed by using a low pass filter and a high-pass filter [14,13,23]. An acoustic low-pass filter may be constructed by inserting a simple expansion chamber (SEC) in a duct. The corresponding transfer matrix can be written as 2 3 π r2 cos ðklcen Þ j ρcen c sin ðklcen Þ 5 ð16Þ TMlow  pass ¼ 4 ρc cos ðklcen Þ jπ r2 sin ðklcen Þ cen

Where rcen and lcen are the radius and length of the central chamber of the SEC. The performance of SEC as a low pass filter has been shown in Fig. 6(a). Similarly, an acoustic high-pass filter may be constructed by inserting a ‘T’ junction, or a short side branch in duct. If both the radius and the length of the side branch are smaller than the cutoff plane wave frequency of the duct, then the acoustic impedance Zhp of the side branch opening can be expressed as k ρðl þ1:5r Þ Z hp ¼ ρc þ j od 2 od 4π π rod 2

ð17Þ

where lod and rod are the length and radius of the open branch duct. The performance of a high pass filter has been shown in Fig. 6(b). The impedance Zhp acts as lumped shunt element and the corresponding transfer matrix can be written as " # 1 0 ð18Þ TMhigh  pass ¼ 1=Z hp 1 In the present work, multiple Helmholtz resonators, band-stop filters, have been proposed to augment the performance of a band-pass filter. The Helmholtz resonator behaves like a simple spring-mass system and cavity resonates on a particular frequency. Analogous to high-pass filter, the impedance Zhr also acts as lumped shunt element, when connected to a duct as side branch. The Zhr is the acoustic impedance of a Helmholtz resonator and can be expressed as [25] Z hr ¼ jZ c cotðkhÞ þ Z h 8 ρc > > < Zc ¼ S c ρc  > > : Z h ¼ S 0:0072 þjkðlneck þ 1:7  r neck Þ h

ð19Þ

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

941

where rneck and lneck are the radius and length of the resonator neck. The Sc and Sh are the cavity and neck cross sectional area. The performance of a Helmholtz resonator as a band-stop filter has been shown in Fig. 6(c). Next, an attempt has been made in introducing a Quincke tube with n-duct configuration to replace the simple expansion chamber (low pass filter) and reduce the requirement of Helmholtz resonators. Selamet and Easwaran [26] have proposed the generalized transfer matrix for n-duct configuration. The corresponding transfer matrix and transmission coefficients can be written as 2 3 T Q11 T Q12 Q 4 5 TM ¼ ð20Þ T Q21 T Q22  8 Pj ¼ n  > j ¼ 1 Aj cos ðklj Þ= sin ðklj Þ > Q > T ¼ > Pj ¼ n 11 > > > > j ¼ 1 Aj = sin ðklj Þ > > > > j ρc > Q > T 12 ¼ Pj ¼ n > > > > > j ¼ 1 Aj = sin ðklj Þ > > < P

  cos ðklj Þ j¼n Aj jX ¼n jX ¼n j ¼ 1 > sin ðklj Þ Aj Aj cos ðklj Þ j j > >    T Q21 ¼ > > > P A ρ c sin ðkl Þ ρ c sin ðklj Þ j j¼n j > j ¼ 1 j ¼ 1 > > j ¼ 1 sin ðkl Þ > > j > > >  Pj ¼ n  > > Aj cos ðklj Þ= sin ðklj Þ > Q > j ¼ 1 > T ¼ > Pj ¼ n > : 22 j ¼ 1 Aj = sin ðklj Þ

ð21Þ

The performance of a Quincke tube with n-duct configuration filter has been shown in Fig. 6(d). 3.3. Tuning of Acoustic Filter Towards this, first, a rudimentary band-pass filter has been designed by combining simple acoustic low-pass, high-pass and band-stop filter components such as a simple expansion chambers, open branch ducts, and Helmholtz resonators [23,13]. Then in order to tune the filter to our requirement, a genetic algorithm (GA) based optimization technique with an appropriate objective function and multiple design constraints has been proposed. Robustness of the proposed technique has been demonstrated by designing a band-pass filter using a 4-branch Quincke tube, a high pass filter and three Helmholtz resonators under specified constraints for the purpose of piston-bore defect identification. Steps of achieving the final design are explained in the following sections. Genetic algorithm (GA) is an evolutionary algorithm based on Darwin's “Survival of the fittest” theory, which mimics the mechanisms of natural evolution. Genetic algorithm was initially introduced by John Holland [28] and has been accepted as established practice in design optimization [29,30]. The algorithm drives a population of designs where the “n” design parameters are coded in a chromosomes. An initial randomized population, Pk, of such chromosomes is built and each individual is evaluated with the objective function. The population evolves from generation to generation (k), which gradually leads towards attaining the objective through a natural selection. The fittest individuals, from “fitness” function, are chosen for the next generation. Next generation chromosomes, P k þ 1 , are obtained from “Crossing” and “Mutation” procedures using the chromosomes of previous generation. The iterative procedure is practiced for given generations. Most often, the best chromosome from the resultant population generate an excellent solution for the optimization problem. The general process of implementing GA can be written as: Pseudo code: Genetic Algorithm (n, κ, mr) // Initialize Generation k ≔ 0; Pk ≔ a population of n randomly generated individuals; // Initialize Generation Pk: Compute fitnessðiÞ for each i A P k ; do {// Create Generation k þ 1; Step-a (Copy): Select ð1  κ Þ  n members of P k - insert in to P k þ 1 Step-b (Crossover): Select κ  n members of P k Pair them off - Produce offspring - Insert in to P k þ 1 ; Step-c (Mutate): Select mr  n members of P k þ 1 - Invert a randomly selected-bit; Step-d (Evaluate P k þ 1 ): Compute fitnessðiÞ for each iA P k þ 1 ; Step-e (Increment): k≔k þ 1; } While: Fitness of the fittest member from Pk o Desired; Return: The fittest member from Pk;

Here, n, κ, and mr are number of population, crossover fraction, and mutation rate, respectively. In optimization, most of the constraint handling methods are based on penalty functions, which penalize infeasible solutions. They solve an

942

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

unconstrained problem (on f) using the modified fitness function ! evalð x Þ ¼

( ! fð x Þ ! ! f ð x Þ þ Penaltyð x Þ

if

! x AF

otherwise

ð22Þ

! where, penalty, ( x ), is zero, if no violation occurs, and is positive, otherwise. Usually, the penalty function is based on the distance of a solution from the feasible region F. For every constraint, a family of intervals has been set which determines appropriate penalty coefficient [31]. A methodology has been established to handle the constraints as mentioned below: Implementing Reward & Penalty Algorithm Step a: Step b: Step c: Step d: Step e:

For each constraint, create several (l) levels of violations. For each level of violation and for each constraint, create a penalty coefficient, ‘Rij’, where, iA ½1…l, and jA ½1…m. Assign larger values of coefficient to corresponding higher levels of violations. Start with a random population, Pk, of individuals (feasible/infeasible). Evolve the population P k þ 1 by evaluating individuals, Pk(i), using below expression

jX ¼m ! ! 2 ! Rij f j ð x Þ evalð x Þ ¼ f ð x Þ þ

ð23Þ

j¼1

For m constraints the method requires mð2l þ 1Þ parameters in total. To establish the number of intervals for each constraint, m numbers of parameters are required and l numbers of parameters for each constraint defining the boundaries of the intervals. Another l parameters are required to represent the penalty coefficients Rij. The detail optimization algorithm has s been developed in MATLAB platform. Now, the desired band-pass filter, using a 4-branch Quincke tube, a high pass filter and three Helmholtz resonators, has been designed using the aforementioned optimization technique. A total of 31 numbers of parameters which represent the detail geometry of each filter component such as length, height, diameter including the distance between them have been optimized with multiple geometrical constants. However, the objective function is the key factor in establishing a robust GA based optimization technique. In the present case, the requirement is to have maximum possible transmission coefficient in band-pass frequency band and other frequency bands should have the minimal sound transmission coefficient. The proposed objective function can be written as OF Minimization ¼ k1  OF 1  k2 OF 2 þ k3 OF 3 ; ðk1;2;3 A ½1; 2; 3…NÞ 8 f h ¼ Lower > XLimit  1 > > > OF 1 ¼ T π ; ðMinimizationÞ > > > > fl ¼ 1 > > > > > f h ¼ Higher < X Limit OF 2 ¼ T π ; ðMaximizationÞ > > f l ¼ Lower Limit > > > > f h ¼ Duct Cutoff > X Frequency > > > OF 3 ¼ T π ; ðMinimizationÞ > > > : f l ¼ Higher Limit þ 1

ð24Þ

Here, T π ¼ at , and fl, fh are the lower and higher limit of frequency, respectively. The constants k1 ; k2 and k3 are the weights assigned to each sub-objective function. These weights have been found from numerical experimentation. In the present case the value of k1 ; k2 and k3 are set to 1,  3 and 2, respectively. Using the generated dimensions from the optimization process, a 3D model has been prepared which represent the s acoustic cavity of the filter using SOLID WORKS and has been shown in Fig. 7(a). The optimized filter coefficient has been shown in Fig. 7(c). Before fabricating the filter for experimental validation, its acoustic behavior has been evaluated using finite element method (FEM). For this purpose, the acoustic cavity, shown in Fig. 7(a), has been analyzed using a coms mercially available FEM based software, namely ANSYS . For understanding the acoustical characteristics of the filter, the acoustic transmission loss (TL) across the filter has been estimated numerically. The potential of FEM in estimating acoustic s TL using ANSYS has been discussed in subsequent chapters. The generated numerical results such as acoustic TL and the corresponding sound power transmission coefficient have been compared with analytical results in Fig. 7(d). From the numerical and the analytical results, it is worth noticing that the performance of the optimized filter is very close to required s performance. The designed filter has been manufactured using rapid prototyping machine (FORTUS400mc ) and has been s shown in Fig. 8(a). The acoustic properties of the filter have been measured using B&K transmission loss tube setup and the TL curve has been compared with the analytical and those from the numerical results (see Fig. 8(b)). A good agreement between all three technique can be observed.

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

943

Fig. 7. Optimized acoustic filter and corresponding coefficient; (a) filter cavity, (b) 3D model of acoustic filter, (c) acoustic transmission coefficient, (d) acoustic transmission loss.

Fig. 8. Proposed system; (a) designed acoustic filter mounted with microphone and anechoic terminator, (b) transmission coefficient of the designed filter.

3.4. Signal conditioning through Passive Acoustic Filter s

The designed passive acoustic filter, microphone holder, hand held sound analyzer (B&K 2270 ) and anechoic termis nation (B&K ) have been assembled as shown in Fig. 8(a). This system has been used for acquiring the engine noise at 48 kHz sampling rate. Similar to earlier exercise, no additional amplification or weighing filters have been incorporated during data acquisition. Next, the experiments have been carried out with the same vehicles. Similar to previous experiments, a total of 950 test samples of duration 0.5 s have been acquired (350 test samples from health vehicle and 300 test samples from first defective vehicle (DV_Signala0 ) and second defective vehicle (DV_Signalb0 )). The acquired samples of 0.5 s duration from the healthy 0 and defective vehicles have been shown in Fig. 9(a)–(c) and referred as HV_Signal , (DV_Signala0 ), and (DV_Signalb0 ), respectively. The corresponding PSD spectra have been shown in Fig. 10(a)–(c). Similar behavior in PSD spectrums have been observed such as prominent fluctuations in frequency components below 300 Hz for defective vehicles. Similar to earlier investigation, out of the 950 acoustic signal samples, which include samples both from the defective 0 and the non-defective vehicles, 800 random samples (300 test samples from HV_Signal and 250 test samples from 0 DV_Signala0 and DV_Signalb , each) have been chosen randomly to train the machine learning algorithm. The rest unused 0 150 samples (50 test samples from HV_Signal and 50 test samples from DV_Signala0 and DV_Signalb0 , each) have been used

944

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

0

Fig. 9. Noise signal from (a) healthy vehicle (HV_Signal ), (b) faulty vehicle (DV_Signala0 ), and (c) faulty vehicle (DV_Signalb0 ).

0

Fig. 10. PSD spectra of noise signal from (a) healthy vehicle (HV_Signal ), (b) faulty vehicle (DV_Signala0 ), (c) faulty vehicle (DV_Signalb0 ). Table 5 Mean and coefficients of variation (CV) (in %) of statistical parameters of engine noise signal samples. 0

Parameter name (1–6)

DV_Signala0

HV_Signal

Crest factor SNR Skewness Kurtosis Peak frequency (Hz) Peak frequency amplitude

DV_Signalb0

Mean

CV %

Mean

CV %

Mean

CV %

3.59 0.82 0.07 2.5 69.9 0.048

4.16 2.29 206.8 12.08 41.6 24.95

4.25 0.8  0.11 2.87 156.63 0.044

5.88 2.3 198.6 12.24 62.2 25.1

4.2 0.78  0.21 3.3 184.2 0.04

3.23 2.43 105.2 12.19 45.6 14.5

Table 6 Evaluation of acoustic signal (after passive filter) statistical parameters with respect to kernel functions of SVM using least square method in (%). Parameters

Linear

Quadratic

Polynomial

RBF

MLP

Parametr1 Parametr2 Parametr3 Parametr4 Parametr5 Parametr6 Parameters1–6

92.8 84.36 69.9 78.64 80.02 78.96 94.62

93.7 89.47 72.8 89.96 76.58 82.64 97.6

92 85.46 69.28 85.6 80.2 78.8 95.48

94.68 88.6 71.29 88.46 83.64 79.84 97.7

92 66.57 71.28 84.56 92.68 78.59 92.25

for testing and validation of the trained model. The six statistical parameters, identified earlier, have been used to train the machine learning algorithm and have been tabulated in Table 5 (mean and coefficients of variation of all samples). Next, the previously designed multi-layer (with 10 hidden layers) feed-forward back-propagation neural network (FBNN) has been put to use. After successive training of the FBN network using 800 training samples, the model has been tested with 150 test samples (that were kept aside initially). The healthy and the defective vehicles can be observed to be classified with an average of 98.59% and 98.16% success rates, respectively. The observed recognition rate is 98.67% which is higher than earlier observation 91%. Afterwards, earlier discussed machine learning technique, i.e., SVM has been investigated using the same six statistical parameters. The computed average recognition rate using individual parameters and their combination are tabulated in Tables 6 and 7. From these tables, it is clear that the maximum recognition rate of 99.33% can be achieved using RBF kernel function along with sequential minimal optimization method using only six statistical parameters which is higher than

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

945

Table 7 Evaluation of acoustic signal (after passive filter) statistical parameters with respect to kernel functions of SVM using sequential minimal optimization in (%). Parameters

Linear

Quadratic

Polynomial

RBF

MLP

Parametr1 Parametr2 Parametr3 Parametr4 Parametr5 Parametr6 Parameters1–6

92.24 88.46 68.74 86.04 89.56 76.8 92.4

93.98 89.94 66.75 85.8 82.45 78.6 97.2

92.48 87.63 70.1 83.56 89.47 77.58 96.64

93.42 90.12 71.8 87.42 91.48 78.5 99.33

92.8 86.4 66.35 81.28 84.26 77.64 71.4

Fig. 11. Evaluation of piston-bore defect recognition-rate using ANN and SVM with different numbers of training samples.

earlier observation 94.4%. The dependency of the proposed methods on number of training samples has been shown in Fig. 11. From figure, it is well understood that SVM is more capable with lower number of training samples.

4. Conclusion Experimental investigations divulge that the engine noise carries motor bike piston-bore defect signature in a band of 50–300 Hz. No particular frequency peak can explain the presence of the defect. The machine learning mechanisms such as ANN and SVM are quite capable to produce fault model using identified standard six statistical parameters such as crest factor (CF), SNR, skewness, kurtosis, peak frequency (Hz), peak frequency amplitude. However, the achieved fault detection rate is less than 94%. In order to achieve high reliability of the expert system, a passive acoustic band-pass filter has been optimized and designed using genetic algorithm. On implementing the designed passive acoustic filter in noise path during the acoustic signal acquisition, the accuracy of FBNN and SVM have been improved and observed to be 98.67% and 99.33% using only six statistical parameters. The whole exercise has been carried out in TVS workshop environment and the adequacy of the proposed system has been demonstrated. SVM using RBF kernel function along with sequential minimal optimization method is to be recommended when less data samples are available for analysis.

Acknowledgment First author acknowledges the scholarship received from the Ministry of Human Resource Development, Government of India, for carrying out this work at Indian Institute of Technology, Bhubaneswar.

References [1] A.K.S. Jardine, D. Lin, D. Banjevic, A review on machinery diagnostics and prognostics implementing condition-based maintenance, Mech. Syst. Signal Process. 20 (2006) 1483–1510. [2] P. Konar, P. Chattopadhyay, Bearing fault detection of induction motor using wavelet and Support Vector Machines (SVMs), Appl. Soft Comput. J. 11 (2011) 4203–4211. [3] Y. Lei, J. Lin, M.J. Zuo, Z. He, Condition monitoring and fault diagnosis of planetary gearboxes: a review, Measurement 48 (2014) 292–305. [4] A. Albarbar, F. Gu, A.D. Ball, A. Starr, Acoustic monitoring of engine fuel injection based on adaptive filtering techniques, Appl. Acoust. 71 (2010) 1132–1141. [5] S.G. Barad, P.V. Ramaiah, R.K. Giridhar, G. Krishnaiah, Neural network approach for a combined performance and mechanical health monitoring of a gas turbine engine, Mech. Syst. Signal Process. 27 (2012) 729–742. [6] Y.S. Wang, Q.H. Ma, Q. Zhu, X.T. Liu, L.H. Zhao, An intelligent approach for engine fault diagnosis based on HilbertHuang transform and support vector machine, Appl. Acoust. 75 (2014) 1–9.

946

D.P. Jena, S.N. Panigrahi / Mechanical Systems and Signal Processing 70-71 (2016) 932–946

[7] G.-X. Jing, Z.-Y. Hao, A new technique based on traditional wavelet transform used in NVH application of internal combustion engine, Mech. Syst. Signal Process. 23 (2009) 979–985. [8] N. Tandon, B. Nakra, B. Sarkar, V. Adyanthaya, Noise control of two-wheeler scooter engine, Appl. Acoust. 51 (1997) 369–380. [9] A. Albarbar, F. Gu, A.D. Ball, Diesel engine fuel injection monitoring using acoustic measurements and independent component analysis, Measurement 43 (2010) 1376–1386. [10] J.-D. Wu, C.-H. Liu, An expert system for fault diagnosis in internal combustion engines using wavelet packet transform and neural network, Expert Syst. Appl. 36 (2009) 4278–4286. [11] J.D. Wu, P.H. Chiang, Y.W. Chang, Y. Jung Shiao, An expert system for fault diagnosis in internal combustion engines using probability neural network, Expert Syst. Appl. 34 (2008) 2704–2713. [12] K.C. Gryllias, I.A. Antoniadis, A Support Vector Machine approach based on physical model training for rolling element bearing fault detection in industrial environments, Eng. Appl. Artif. Intell. 25 (2012) 326–344. [13] H.F. Olson, Dynamical Analogies, 1st edition, . D. Van Nostrand Company Inc., New York, 1944. [14] M.L. Munjal, Acoustics of Ducts and Mufflers With Application to Exhaust and Ventilation System Design, 1st edition, . John Wiley and Sons, New York, 1987. [15] F. Mechel, Formulas of Acoustics, 1st edition, . Springer, Berlin, 2002. [16] M. Munjal, Recent advances in muffler acoustics, Int. J. Acoust. Vib. 18 (2013) 71–85. [17] J. Rafiee, F. Arvani, A. Harifi, M.H. Sadeghi, Intelligent condition monitoring of a gearbox using artificial neural network, Mech. Syst. Signal Process. 21 (2007) 1746–1754. [18] A. Saxena, A. Saad, Evolving an artificial neural network classifier for condition monitoring of rotating mechanical systems, Appl. Soft Comput. 7 (2007) 441–454. [19] V. Sugumaran, K.I. Ramachandran, Expert Systems with Applications Effect of number of features on classification of roller bearing faults using SVM and PSVM, Expert Syst. Appl. 38 (2011) 4088–4096. [20] P.K. Kankar, S.C. Sharma, S.P. Harsha, Fault diagnosis of ball bearings using machine learning methods, Expert Syst. Appl. 38 (2011) 1876–1886. [21] B.S. Bernhard Scholkopf, Alexander J. Smola, Learning With Kernels: Support Vector Machines, Regularization, Optimization and Beyond, 1st edition, . The MIT Press Ltd, London, 2001. [22] L. Bottou, C.-J. Lin, Support Vector Machine Solvers, in Large Scale Kernel Machines, 1st edition, . MIT Press, Cambridge, 2007. [23] G.W. Stewart, Acoustic waves filters, Phys. Rev. 20 (1922) 528–551. [24] R.F. Barron, Industrial Noise Control and Acoustics, 1st edition, . Marcel Dekker, Inc., New York, 2003. [25] S.-H. Seo, Y.-H. Kim, Silencer design by using array resonators for low-frequency band noise reduction, J. Acoust. Soc. Am. 118 (2005) 2332–2338. [26] A. Selamet, V. Easwaran, Modified Herschel–Quincke tube: attenuation and resonance for n-duct configuration, J. Acoust. Soc. Am. 102 (1997) 164–169. [27] J.D. Smith, Gear Noise and Vibration, 2nd edition, . Marcel Dekker, Inc., New York, 2003. [28] D. Goldberg, Genetic Algorithms in Search, Optimization and Machine Learning, 1st edition, . Addison-Wesley, New York, 1989. [29] N. Kaya, Machining fixture locating and clamping position optimization using genetic algorithms, Comput. Ind. 57 (2006) 112–120. [30] T. Airaksinen, E. Heikkola, Multiobjective muffler shape optimization with hybrid acoustics modeling, J. Acoust. Soc. Am. 130 (2011) 1359–1369. [31] D. Bunnag, M. Sun, Genetic algorithm for constrained global optimization in continuous variables, Appl. Math. Comput. 171 (2005) 604–636.