17th Australasian Fluid Mechanics Conference

Auckland, New Zealand

5-9 December 2010

Statistical Modeling to Determine Resonant Frequency Information Present in Combustion

Chamber Pressure Signals
T. A. Bodisco and R. J. Brown
School of Engineering Systems, Faculty of Built Environment and Engineering
Queensland University of Technology, Queensland 4001, Australia

Abstract

A statistical modeling method to accurately determine combustion chamber resonance is proposed and demonstrated. This
method utilises Markov-chain Monte Carlo (MCMC) through
the use of the Metropolis-Hastings (MH) algorithm to yield a
probability density function for the combustion chamber frequency and find the best estimate of the resonant frequency,
along with uncertainty. The accurate determination of combustion chamber resonance is then used to investigate various
engine phenomena, with appropriate uncertainty, for a range
of engine cycles. It is shown that, when operating on various
ethanol/diesel fuel combinations, a 20% substitution yields the
least amount of inter-cycle variability, in relation to combustion
chamber resonance.
Introduction

Fluid phenomena frequently involve pressure fluctuations

which are analysed using a variety of techniques. In-cylinder
pressure has been a focal point of engine research since
its introduction[8, 13]. Routinely during engine research
in-cylinder pressure is used to investigate the indicated power,
indicated mean effective pressure, peak pressure, maximum
rate of change of pressure, heat release, and thermal efficiency
of engines[10, 16]. However, as the greatest rate of change of
pressure occurs during combustion near top dead centre (TDC)
of the piston motion, where the velocity of the piston is low, it
is necessary to examine the Pressure vs Crank-angle (p-), or a
Pressure vs Time, plot to obtain the detailed characteristics of
combustion[10]. It is this region that is of interest.
Some interesting studies have been performed by other
researchers in relation to resonant frequencies in combustion
chambers from in-cylinder pressure[11, 15, 17, 18]. Isolation
of the resonant frequency is important as it is related to the
speed of sound and hence temperature[11, 15, 17, 19]. Since,
many assumptions are made when analysing data, it is of
great importance that the implications of these assumptions are
understood[3]. Typically, researchers interested in combustion
resonance have used Fast Fourier Transforms (FFT) to find the
desired frequency information[11, 17, 19].
Fast Fourier transforms are common practice in basic spectral
analysis; mostly because of their ease of use and computational
efficiency[4, 6, 14]. However, assumptions which underpin
them, for example the assumption of periodic stationary
frequencies, and their limited resolution can make them an undesirable tool for use in spectral analysis[3, 12]. Further, FFTs,
and many other traditional techniques, also have problems
when there is noise or incomplete data[12]. Limited data is also
a major draw back with the FFT approach to spectral analysis,
to produce any useful results multiple periods of data are
required[9]in some applications, such as the one proposed,
analysing multiple periods simultaneously is counter to the aim.
Bretthorst explains that if there are complex phenomena

or evidence of more than a single stationary harmonic frequency the Fourier method may yield incorrect or misleading
results[3]. Such results obtained with FFTs are not incorrect
because the FFT is wrong; they are simply not correct because
the FFT is attempting to answer a different question[5]. Hence,
other methods of spectral analysis need to be investigated to
ensure that the information found is accurate and therefore
usable for making decisions from or determining frequency
information present in the data. Emphasis should be given more
to scientific interest and less to mathematical convenience[2].
For the purpose of this investigation spectral analysis
was performed using Markov-chain Monte Carlo (MCMC) in a
Bayesian framework; specifically, using an implementation of
the Metropolis-Hastings Algorithm[7]. In situations where the
user can approach the problem from an informed perspective,
using this methodology has the advantage that it requires the
user to explicitly state any assumptions being made in the
calculation[3]. Consequently, the user always knows exactly
what problem is being solved. A further advantage is that this
methodology can allow for noise in the datawhich is always
present for in-cylinder pressure data.
Experimental Configuration

Experiments were conducted on a modern turbo-charged 6cylinder Cummins diesel engine (ISBe220 31) at the QUT Biofuel Engine Research Facility (BERF). The engine has a capacity of 5.9 l, a bore of 102 mm, a stroke length of 120 mm, a compression ratio of 17.3:1 and maximum power of 162 kW at 2500
rpm. The engine was coupled to an electronically controlled hydraulic dynamometer with load applied by increasing the flow
rate of water inside the dynamometer housing. In-cylinder pressure was measured by a Kistler piezoelectric transducer with
a Data Translation simultaneous analogue-to-digital converter
connected to a desktop computer running National Instruments
LabView. Data was collected at a sample rate of 200,000 samples per second. During testing the engine was run on neat
diesel fuel and diesel/ethanol at 2000 rpm on full load.
Experimental Data

Data was collected with the engine running on neat diesel

fuel, and also with ethanol fumigation substituting 10%,
20%, 30%, 40% and 50% of the diesel on an energy basis.
Specific data that was collected included: in-cylinder pressure,
band-pass filtered in-cylinder pressure (180 kHz), crank-angle
information (updated every degree) and injection information.
In-cylinder pressure was collected as a differential voltage signal and for this application is not converted to absolute pressure.
As this investigation has been focused on combustion
resonance the injection information was used to define the
region of interest and the band-pass filtered in-cylinder pressure
data was used for the calculation. The frequencies of interest
ranged from approximately 5 kHz to 10 kHz, the band window
was set to be very wide to ensure that the frequencies of interest

were unaffected by the filteringthe methodology employed

allows for ignoring or modeling any additional phenomena
present in the signal making filter selection less important.
In this scenario it was only important to remove the low
frequencies dealing with the speed of the engine. An example
of this data is shown in Figure 1.

In this implementation each parameter is updated 5 times

per iteration. Further, before any values are saved there is
a burn-in period of 10,000 iterations to ensure that each
parameter has converged to a value. At the end of this period
1000 iterations are then performed where each parameter of
interest is saved in an array at the end of each iteration. These
saved values allow a probability density function (pdf) of
each parameter to be produced where the mode is taken to
be the best estimate of the parameter and the spread gives an
indication of the uncertainty.
The Model

Figure 1: Combustion chamber pressure fluctuation

after application of a band-pass filter (180 kHz):
2000 rpm, full load, neat diesel fuel
Mathematical Overview

The Metropolis-Hastings algorithm is an adaption of a random

walk that uses an acceptance/rejection rule to converge to the
specified target distribution[7]. For this implementation the random walk is performed by sampling from a Normal Distribution
with the current value as the mean and some adaptable standard
deviation. Each parameter has its own adaptable standard deviation that is updated every 300 trials, with the aim of controlling
the acceptance to rejection ratio at 40%, using the following formula:
= eAcc0.4 ,
(1)
where is the updated standard deviation, the current
standard deviation and Acc is the ratio of acceptance/rejection
for the previous 300 trials. This ensures that the results
obtained are useful; a 100% acceptance rate would indicate that
we are not exploring the space properly and a 0% acceptance
rate would indicate that we are not exploring the space at all.
Acceptance or rejection of each candidate value is driven
by the data and by prior information. If the signal is defined as
being normally distributed with a time varying mean:
y(n) N((n), ),

(2)

where y(n) is the signal, (n) is the mathematical model that

we define and is the standard deviationmodeling this way
accounts for noise present in the datathen the probability of
acceptance, R, can be defined as:


p(yi |i , ) p( )
R = min 1,

,
(3)
p(yi |i , ) p()
where is the parameter in question and denotes the use of
the candidate parameter. The more probable that the candidate
value is the true value, compared to the current value, the closer
to 1 R becomes. For the case that R < 1 a random number is
sampled between 0 and 1 from a uniform distribution and if R
is greater than this value the candidate value is accepted.

Examination of Figure 1 reveals four important features. Firstly,

a period of very high frequency fluctuation associated with
diesel injection; secondly, a peak associated with ignition;
thirdly, a rapidly decaying low frequency shape; finally, the resonant frequency information which can be seen superimposed
on the basic signal. Ignoring all but the resonant frequency and
the low frequency shape, a model can be created that fits this
data and finds the required information. More complex models
can be employed if further information is desired, for example a term could be included to model the injection frequency
or other frequency modes (radial, circumferential, axial); however, in this investigation, of interest is the first circumferential
mode resonant frequency only. In this instance, to eliminate the
need to model the injection, or model anything before the initial peak, a step function is used. Modeling just these features
fits the data sufficiently to produce usable results whilst still being general enough to work across a range of cycles and engine
operating modes, allowing cyclic behaviour to be examined as
well as a comparison between engine operating conditions. The
model is, therefore, as follows:



21 ea1 n
n + 2
(n) = H(n ) A1 e1 n sin
samplerate


a
n
2
22 e
+ A2 e2 n sin
n + 2
,
(4)
samplerate
where H(n ) is a step function that is equal to 0 when n <
and 1 when n , and 2 ea2 n represents the resonant frequency.
The use of prior information for each parameter ensures
that there are no issues with label switching and that we fit the
data in the intended way. For instance, a1 is given a uniform
prior to contain it between 0.001 and 0.01. This ensures that
it is positive, large enough to decay at a rate sufficient to
model the lower frequency signal and that it does not cause
complications with the modeling of the higher frequency
signal. Similarly, a2 is also given a uniform prior; however,
it is contained between 0.0001 and 0.001 to ensure that it
models the specific frequency of interest. In this framework
it is acceptable to give a parameter an uninformative prior
which has little effect on whether a candidate value is accepted
or rejected; for example A1 and A2 were given normally
distributed priors with relatively large standard deviations
making them effectively entirely data driven. An example of
the model fit, from equation (4), is shown in Figure 2. All
results were obtained using NetBeans as the development
environment utilising the standard gcc c++ compiler.
Results

In order to determine the value, and uncertainty, of any model

parameter we typically plot its multiple values as determined by
the Monte Carlo as a pdf. Figure 3 shows the pdf of parameter
2 and Figure 4 shows the pdf of a2 . Results obtained from
figures 3 and 4 allow characterisation of the resonant frequency
as a function of time. However, for the purpose of this

investigation the instantaneous frequency located at the top of

the first peakseen in Figure 1will be termed the reference
resonant frequency and henceforth denoted as r . An example
of the results running the engine at 2000 rpm at full load for
15 s is shown in Figure 5. This figure depicts the significant
inter-cycle variability present in a diesel engine running under
normal conditions. r ranges from 5200 to 6400 Hz with
the modal reference resonant frequency being 5450 Hz with
cycle-to-cycle deviations as great as 1000 Hz.

The resonant frequency is directly

temperature[11]:
T f 2,

related

to

the
(5)

therefore we can use the resonant frequency as a way of

commenting on the repeatability of combustion. Figure 5
clearly shows that even whilst running on neat diesel fuel there
is significant inter-cycle variability. Further, an observable
feature of this plot is that typically prior to a high frequency
peak is a cycle with a slightly lower frequencyperhaps an
indication of unburned fuel present from one cycle having an
effect on the next. Figure 6 shows the pdfs of the r for 250
consecutive cycles with the engine running at 2000 rpm on full
load on neat diesel fuel (D100E00) and with 10% (D090E010),
20% (D080E020), 30% (D070E030), 40% (D060E040) and
50% (D050E050) ethanol substitutions.

Figure 2: Comparison of expected resonant frequency signal from the second term of the model,
using equation (4) on the signal in Figure 1.

Figure 5: Reference resonant frequency, r , for 250

sequential cycles with the engine at 2000 rpm and
full load at selected ethanol substitutions.

Figure 3: Pdf of 2 from equation (4) and the signal

in Figure 1.

Figure 6: Pdfs of r of various engine operating conditions.

Figure 4: Pdf of a2 from equation (4) and the signal

in Figure 1.

It can be seen in Figure 6 that small ethanol substitutions lead

to a slight decrease in r . A 10% ethanol substitution gave
a decrease of approximately 2%. Also, from a combustion
frequency perspective Figure 6 indicates that at 20% ethanol
substutition there is the least variation from cycle-to-cycle and
when compared to neat diesel fuel a slight increase in frequency. A 20% ethanol subsitution leads to improved thermal
efficiency[1]. This result validates the improved efficiency,
indicated by the slightly higher combustion temperature and
the improved inter-cycle consistency.

At high ethanol substitutions, 30%, 40% and 50%, the

repeatability of each combustion becomes increasingly low,
compared to neat diesel and lower ethanol substitutions. It can
be seen in Figure 6 that the pdf of r becomes significantly
flatter, with a marked increase in standard deviation, as
the ethanol substitution increases, indicating the increasing
inter-cycle variability as the ethanol substitution is increased.
We conclude that up to 20% ethanol substitution, from a
combustion chamber resonance perspective, does not significantly alter inter-cycle variability in a negative way, with some
evidence to suggest improved inter-cycle variabilty at the 20%
ethanol substitution. However, at ethanol substitutions of 30%
and higher, inter-cycle variability is severely influenced in a
negative manner.
Conclusions

This paper demonstrated a statistical modeling method

for determining resonant frequency information present in
combustion chamber pressure signals. As such, accurate
determination of instantaneous frequency information was also
demonstrated. Results of this were shown for a large number
of consecutive cycles with the engine operating under various
ethanol substitution rates to demonstrate a novel utility of
this methodology. Using this cycle-to-cycle evidence evinces
support for improved inter-cycle variability with 20% ethanol
substitution. Higher ethanol substitutions were shown to
significantly increase inter-cycle variability.
A further implication of the results is that evidence was
found to support not using ad hoc methodology, such as cycle
averaging, when the desired output is frequency information.
This is particularly true when using frequency information
to make arguments during potentially unstable operating
conditions, such as the example used for this paper of alternative fuel substitutions, with a variation upwards of 2000 Hz.
Only through the use of more sophisticated techniques can
information such as this be accurately obtained, particularly in
light of the small number of data points and noise present.
Acknowledgements

We wish to thank Mr Tony Morris and Mr Glenn Geary for assisting with the design and running of the experimental campaigns. Further thanks also to technologist Mr Ken McIvor for
his assistance in setting up the data acquisition software and to
the HPC and Research Support department for their assistance
with debugging the c++ code, specifically Dr Neil Kelson and
Mr Timothy Gurnett. Also, thanks to Mr Nicholas Surawski
for his work in making the experiments a success. This work
was undertaken under an Australian Research Council Linkage
Grant (LP0775178) in association with Peak3 P/L.
References

[1] Abu-Qudais, M., Haddad, O. and Qudaisat, M., The effect of alcohol fumigation on diesel engine performance
and emissions, Energy Conversion and Management, 41,
2000, 389399.
[2] Box, G. and Tiao, G., Bayesian Inference in Statistical
Analysis, John Wiley and Sons, 1992.
[3] Bretthorst, G. L., Excerpts from Bayesian spectrum analysis and parameter estimation, Maximum Entropy and
Bayesian Methods in Science and Engineering, 1, 1988,
75145.

[4] Cooley, J. W., Lewis, P. A. W. and Welch, P. D., The Fast

Fourier transform and its applications, Education, IEEE
Transactions on, 12, 1969, 2734.
[5] Dou, L. and Hodgson, R., Bayesian inference and Gibbs
sampling in spectral analysis and parameter estimation. i,
Inverse Problems, 11, 1995, 10691085.
[6] Duhamel, P. and Vetterli, M., Fast Fourier transforms: a
tutorial review and a state of the art, Signal Processing,
19, 1990, 259299.
[7] Gelman, A., Carlin, J. B., Stern, H. S. and Rubin, D. B.,
Bayesian Data Analysis: Second Edition, Chapman &
Hall/CRC, 2003.
[8] Greene, A. and Lucas, G., The Testing of Internal Combustion Engines, The English Universities Press Limited,
1969.
[9] Gregory, P. C. and Loredo, T. J., A new method for the detection of a periodic signal of unknown shape and period,
Astrophysical Journal, 398, 1992, 146168.
[10] Heywood, J. B., Internal Combustion Engine Fundamentals, McGraw-Hill, Inc., 1988.
[11] Hickling, R., Feldmaier, D. A., Chen, F. H. K. and Morel,
J. S., Cavity resonances in engine combustion chambers
and some applications, Acoustical Society of America
Journal, 73, 1983, 11701178.
[12] Jaynes, E. T., Bayesian spectrum and chirp analysis,
in in Proceedings of the Third Workshop on MaximumEntropy and Bayesian Methods (1983), editors C. Smith
and D. Reidel, Boston, 1987.
[13] Lujn, J. M., Bermdez, V., Guardiola, C. and Abbad, A., A
methodology for combustion detection in diesel engines
through in-cylinder pressure derivative signal, Mechanical
Systems and Signal Processing, In Press, Corrected Proof.
[14] Papandreou-Suppappola, A., editor, Applications in TimeFrequency Signal Processing, CRC Press, 2003.
[15] Payri, F., Broatch, A., Tormos, B. and Marant, V., New
methodology for in-cylinder pressure analysis in direct injection diesel engines: application to combustion noise,
Measurement Science and Technology, 16, 2005, 540
547.
[16] Randolph, A. L., Methods of processing cylinder-pressure
transducer signals to maximize data accuracy, Society of
Automotive Engineers.
[17] Ren, Y., Randall, R. B. and Milton, B. E., Influence of the
resonant frequency on the control of knock in diesel engines, Proceedings of the Institution of Mechanical Engineers: Part D : Journal of Automobile Engineering. London, 213, 1999, 127133.
[18] Stankovic, L. and Bohme, J. F., Time-frequency analysis
of multiple resonances in combustion engine signals, Signal Processing, 79, 1999, 1528.
[19] Torregrosa, A. J., Broatch, A. and Margot, X., Combustion chamber resonances in direct injection automotive
diesel engines: a numerical approach, International Journal of Engine Research, 5, 2004, 8391.