EXPERIMENTAL AND NUMERICAL STUDY OF PREMIXED FLAME PROPAGATION
IN A CONSTANT VOLUME COMBUSTION CHAMBER
By
Ashwin Hariharan

A THESIS
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
MASTER OF SCIENCE
Mechanical Engineering
2012

ABSTRACT
EXPERIMENTAL AND NUMERICAL STUDY OF PREMIXED FLAME PROPAGATION
IN A CONSTANT VOLUME COMBUSTION CHAMBER
By
Ashwin Hariharan

This work presents the experimental and numerical investigation of premixed flame propagation
in a constant volume rectangular channel that serves as a combustion chamber. Flames inside
closed vessels undergo distinct transformations in their physical structure; ignition event is
followed by an accelerating convex shaped flame-front, then a deceleration of the flame is
followed by the formation of a concave “tulip” shaped flame-front. Eventually, the flame is
dissipated either because of quenching by the cold wall or by collision with a flame that has
been ignited on the opposite end of the confined channel. Phenomena observed include both the
formation of concave flame-fronts and oscillating convex-concave flame-fronts. The flame
structure and speed of flame propagation is characterized for changes in the equivalence ratio,
initial pressure, timing and placement of the ignition sources. Numerical and analytical
calculations of the combustion phenomena are performed to understand the effect of the pressure
waves and various instabilities causing the changes to the flame structure. The transient
numerical results are validated using the experimental results.

Copyright by
ASHWIN HARIHARAN
2012

ACKNOWLEDGEMENTS

First and foremost, I would like to thank my advisor, Professor Indrek S. Wichman, who has
guided and supported me throughout my graduate studies. I have been very fortunate to have
access to his vast knowledge on combustion physics and mathematics. His constant
encouragement and ideas have helped me greatly in working towards my research goals.
A special thanks goes to Professor Norbert Mueller, who has given me the freedom to work on
both the numerical analysis of the wave disk engine and on experimental investigation of flame
propagation in combustion chambers. I have benefitted greatly in sitting through the hours-long
weekly meeting with the WDE group, and have been able to understand and appreciate the
challenges faced in transferring computational design to an experimental prototype. I would also
like to thank Professor Tonghun Lee for his input on the wave engine design and for helping me
with combustion diagnostics.
I would like to thank Mike Varney for helping me setup the electronic controls and diagnostics
for the experimental apparatus and teaching me LabView. I would also like to thank all my coworkers and friends for the endless hours of discussion on a variety of topics.
Finally, I would like to thank my parents and sister for their endless love and support during my
time spent in graduate school thousands of miles away from home.

iv

TABLE OF CONTENTS

LIST OF TABLES ........................................................................................................................vii
LIST OF FIGURES .....................................................................................................................viii
KEY TO SYMBOLS ....................................................................................................................xii
INTRODUCTION ..........................................................................................................................1
CHAPTER 1
WAVE DISK ENGINE .................................................................................................................2
Wave Disk Engine Introduction ..............................................................................................3
Combustion Chamber Design ................................................................................................10
CHAPTER 2
LAMINAR FLAMES....................................................................................................................17
Laminar and Turbulent Premixed Flames................................................................................20
Laminar Flame Speed..............................................................................................................21
Mallard and Le-Chatelier Theory......................................................................................24
Zeldovich, Frank-Kamenetskii and Semenov Theory.......................................................29
Flame Speed Measurements....................................................................................................31
Instabilities in Flame Propagation...........................................................................................34
Landau-Darrieus Instability...............................................................................................35
Thermo-diffusive Instability..............................................................................................39
Thermo-acoustic Instability...............................................................................................42
CHAPTER 3
EXPERIMENTAL SETUP, TEST METHODOLOGY AND DATA ACQUISITION...............46
Combustion Chamber Setup....................................................................................................47
Schlieren Setup........................................................................................................................54
Test Procedure.........................................................................................................................55
CHAPTER 4
EXPERIMENTAL RESULTS......................................................................................................57
Single Spark Ignition...............................................................................................................58
Dual Spark Ignition..................................................................................................................67
Dual Spark Ignition, 0 ms delay........................................................................................68
Dual Spark Ignition, 25 ms delay .....................................................................................73
Dual Spark Ignition, 50 ms delay .....................................................................................77

v

CHAPTER 5
NUMERICAL RESULTS.............................................................................................................82
CHAPTER 6
RESULTS AND DISCUSSION....................................................................................................96
CHAPTER 7
CONCLUSIONS AND FUTURE WORK..................................................................................104
APPENDICES.............................................................................................................................105
Appendix A ..........................................................................................................................106
Appendix B ..........................................................................................................................112
Appendix C ..........................................................................................................................116
REFERENCES ...........................................................................................................................119

vi

LIST OF TABLES

Table 1.1 Wave engine cycle timing, theoretical flame speed .......................................................8
Table 1.2 Wave engine cycle timing, experimental flame speed ...................................................9
Table 1.3 Wave engine cycle timing, third generation rotor ........................................................15

vii

LIST OF FIGURES

Figure 1.1 Wave disk engine schematic representation ..................................................................2
Figure 1.2 Computation domain of the straight channel .................................................................6
Figure 1.3 Pressure inside the wave engine channel .......................................................................7
Figure 1.4 Wave engine cycle time, in percentage .........................................................................9
Figure 1.5 First generation wave disk engine design ....................................................................10
Figure 1.6 Second generation wave disk engine design ...............................................................12
Figure 1.7 Third generation wave disk engine design ..................................................................13
Figure 1.8 Contour of Mach number during exhaust stage ..........................................................14
Figure 1.9 Contour of static pressure during exhaust stage ..........................................................15
Figure 1.10 Third generation wave engine cycle time ..................................................................16
Figure 2.1 Structure of a premixed flame in one dimension .........................................................20
Figure 2.2 Mallard and LeChatelier model ...................................................................................24
Figure 2.3 Balances across a differential element of laminar flame .............................................31
Figure 2.4 Bunsen burner flame ....................................................................................................34
Figure 2.5 Propagating flame front ...............................................................................................37
Figure 2.6 Inclined flame front .....................................................................................................38
Figure 2.7 Local streamline deviation in a wrinkled flame ..........................................................40
Figure 2.8 Global streamline deviation in a wrinkled flame ........................................................40
Figure 2.9 Wrinkled premixed flame structure .............................................................................43
Figure 3.1 Exploded view of combustion chamber apparatus ......................................................49
Figure 3.2 Isometric view of combustion chamber apparatus ......................................................50

viii

Figure 3.3 Manufactured combustion chamber apparatus ............................................................47
Figure 3.4 Isometric view of center plate .....................................................................................48
Figure 3.5 Orthographic view of center plate ...............................................................................48
Figure 3.6 Isometric view of side plate .........................................................................................49
Figure 3.7 Orthographic view of side plate ..................................................................................50
Figure 3.8 LabView interface .......................................................................................................51
Figure 3.9 Schlieren imaging setup ..............................................................................................52
Figure 4.1 Flame propagation along combustion tube, showing evolution of flame structure ....58
Figure 4.2 Schlieren image of ignition .........................................................................................59
Figure 4.3 Schlieren image of tulip formation ..............................................................................60
Figure 4.4 Flame position vs time, showing the four stages of flame propagation (I=ignition;
II=finger flame formation; III=transition to planar flame; IV=propagation of tulip flame) .........61
Figure 4.5 Flame speed vs time ....................................................................................................62
Figure 4.6 Flame speed vs normalized flame position .................................................................63
Figure 4.7 Pressure vs time ...........................................................................................................64
Figure 4.8 Flame propagation along combustion tube, dual spark with 0 ms delay .....................67
Figure 4.9 Flame position vs time .................................................................................................68
Figure 4.10 Flame speed vs time ..................................................................................................68
Figure 4.11 Flame speed vs normalized flame position ...............................................................69
Figure 4.12 Pressure vs time .........................................................................................................70
Figure 4.13 Flame propagation along combustion tube, dual spark with 25 ms delay .................72
Figure 4.14 Normalized flame position vs time ............................................................................73
Figure 4.15 Flame speed vs time ..................................................................................................74

ix

Figure 4.16 Flame speed vs normalized flame position ...............................................................74
Figure 4.17 Flame propagation along combustion tube, dual spark with 50 ms delay ................76
Figure 4.18 Normalized flame position vs time ...........................................................................77
Figure 4.19 Flame speed vs time ..................................................................................................78
Figure 4.20 Flame speed vs normalized flame position ...............................................................79
Figure 5.1 Combustion chamber computational domain .............................................................82
Figure 5.2 Combustion chamber mesh, spark region ...................................................................83
Figure 5.3 Flame propagation along combustion tube, showing finger formation........................85
Figure 5.4 Flame propagation along combustion tube, showing tulip formation .........................87
Figure 5.5 Velocity vector of flow inside combustion chamber, finger flame .............................90
Figure 5.6 Velocity vector of flow inside combustion chamber, tulip flame ...............................91
Figure 5.7 Unburnt mixture data point .........................................................................................92
Figure 5.8 Pressure vs time in the unburnt mixture ......................................................................93
Figure 5.9 Temperature vs time in the unburnt mixture ...............................................................94
Figure 5.10 Pressure comparison in the unburnt mixture .............................................................94
Figure 5.11 Lewis number of unburnt mixture .............................................................................95
Figure 6.1 Analytical representation of spherical flame ...............................................................97
Figure 6.2 Analytical representation of finger flame ....................................................................98
Figure 6.3 Analytical representation of planar flame ...................................................................99
Figure 6.4 Analytical representation of tulip flame ....................................................................100
Figure 6.5 Analytical representation of tulip flame collision with cold wall .............................101
Figure 6.6 Flame surface area vs time ........................................................................................103
Figure 6.7 Burnt volume vs time ................................................................................................103

x

Figure A.1 Z-type schlieren schematic .......................................................................................115
Figure A.2 Spherical mirror .......................................................................................................116

xi

KEY TO SYMBOLS

m - Mass flow rate
Cp – Specific heat, at constant pressure
λ – Thermal conductivity
δ – Flame thickness
A – Area of cross-section
ρ – Density
U – Velocity
SL – Laminar flame speed
k – Rate of reaction
Ti – Ignition temperature
Tu – Unburnt temperature
Tf – Flame temperature
w – Reaction rate
α – Thermal diffusivity
δL – Laminar flame thickness
Ea – Activation energy

xii

INTRODUCTION
This thesis is devoted to the investigation of the structure and propagation of premixed flames in
constant volume combustion chambers and their application to wave rotor devices. The structure
and speed of flame propagation in a constant pressure combustor is nearly constant as the
combusted gases are allowed to expand under isobaric conditions. In a constant volume
combustor, the temperature and pressure of the unburnt premixed mixture continuously changes
throughout the combustion process, which leads to a varying flame speed. Experimental
investigation performed by Ellis [1] showed that flames in closed tubes undergo distinct changes
in their structure and speed when ignited at one end of the tube. The flame propagation was split
into 3 stages: (1) ignition by spark followed by the growth of a spherical flame front; (2)
quenching of the flame-skirt by the cold wall and flame shape inversion; and finally, (3) collision
with the cold wall on the opposite end of the tube. The inward pointing cusp formed at stage two
was referred to as the “tulip flame” by Salamandra [2].
The primary motivation for this research is the study of the influence of combustion chamber
shape and initial conditions on the propagation of premixed flames in a wave disk engine. The
wave disk engine research, funded by the US Department of Energy [DoE] and the Advanced
Research Projects Agency-Energy [ARPA-E] involves the design of a rotary engine prototype
which potentially utilizes shock-wave compression and constant volume combustion to generate
power at a high thermal efficiency. The design of the engine prototype involved detailed analysis
of the flame propagation and combustion time for the determination of the intake and exhaust
port timing as well as the rotational speed of the engine (RPM).

1

CHAPTER 1
WAVE DISK ENGINE
A wave rotor is an unsteady flow device which uses shock waves to transfer energy between
fluids at different pressure levels. The wave engine uses either an external turbine or curved
channels to extract shaft work and uses the unsteady shock wave compression as a part of the
compression stage. The initial prototype of an axial wave disk engine was constructed by the
Ruston-Hornsby turbine company in the mid-1950s, which produced a maximum power of 26
kW at its design point. Further development of the axial wave engine took place at ABB and
ETH-Zurich which designed an engine to generate up to 100 kW. The wave disk engine concept
as proposed by M. Vagani et. al. [3] at Michigan State University consists of a radial disk with
curved channels which compress premixed fuel-air mixture and expand combusted gases,
producing shaft work. The radial design enables a compact layout of the engine which reduces
cost and results in a high power density. The compression and expansion of the premixed
mixture and burnt gases are achieved through shockwaves and expansion waves. The momentum
of the flow is used to extract shaft work. The wave disk engine developed at MSU is discussed in
greater detail in the following subsection.

2

1.1 WAVE DISK ENGINE INTRODUCTION
The simplified representation of the wave disk engine is shown in Fig. 1.1. The primary design
emphasis at Michigan State University has been on the radial wave engines where power is
generated by the radial turbine action of the channel geometry.
Scavenging

Compression

Jet Propulsion

Expansion

Loading
Constant Volume
Combustion
Figure 1.1 - Wave disk engine schematic representation [3]. “For interpretation of the references
to color in this and all other figures, the reader is referred to the electronic version of this thesis”

The wave disk engine uses the Humphrey cycle, which involves isentropic compression followed
by constant volume combustion and isentropic expansion of the combusted gases. The basic
design element of the radial wave engine consists of a circular rotor which has curved channels
machined inside, which act as combustion chambers. The rotor is surrounded on the outside by
the exhaust assembly, on the inside by the intake port assembly. The intake and exhaust timing is
controlled by the relative angle of opening of the intake and exhaust ports.
3

The initial calculation, based on 1-dimensional straight channel geometry and laminar
combustion speed enabled the design and manufacture of the first generation wave engine at
MSU [3].
The four distinct stages involved in the wave disk engine operation are as follows:
Stage I – Injection of fuel-air mixture
Once the hot exhaust gases at high pressure have been expelled from the previous cycle, the
intake port is opened in conjunction with the exhaust port. This allows the premixed fuel-air
mixture to flow into the combustion chamber of the wave rotor. The exhaust gas expulsion
creates a temporary low pressure regime throughout the combustion chamber due to the action of
expansion fans. This permits the injection of the fuel-air mixture through the mechanism of selfaspiration.
Stage II – Hammer shock compression
While the intake port is open, the exhaust port is abruptly closed. This creates a hammer shock
which travels upstream of the flow of the fuel-air mixture. The hammer shock phenomenon
occurs due to the sudden change in the velocity brought about by the closing of the exhaust port
by the rotation of the channel. The hammer shock also increases the pressure and temperature of
the air-fuel mixture and increases the effective mass inside the combustion chamber. The intake
port is closed once the hammer shock reaches it. The compressed fuel-air mixture is essentially
at zero velocity with respect to a rotating frame of reference attached to the combustion chamber.

4

Stage III – Constant volume combustion
The combustion chamber is maintained at a constant volume by the closure of the intake and
exhaust ports. The fuel-air mixture is ignited by a spark plug placed either on the outer radius of
the combustion chamber or axially above the channel. This enables the fuel-air mixture to
undergo combustion. Once the combustion process is complete, the combustion chamber will
contain gases at high pressure and temperature.
Stage IV – Isentropic expansion
The exhaust port of the channel is suddenly opened once combustion is complete. This enables
the expulsion of the combusted gases through the formation of expansion waves. The pressure
and temperature of the channel gradually decreases. After a certain time, the pressure inside the
channel decreases below the ambient pressure. At this stage, the intake port is opened and fresh
fuel-air mixture is let in as described in stage I.
The four stages described above repeat continuously as the rotor spins to generate power. The
curved shape of the channel enables the extraction of work through jet propulsion principles. The
amount of work extraction can be further increased by the addition of multiple turbine stages.
The emphasis of this thesis was to optimize both the extraction of power from the wave engine
through a single stage turbine expansion process and to develop an accurate numerical model to
estimate the time required for combustion. Extensive numerical analyses were carried out to
optimize the structure of the channels and to study the structure and speed of flame propagation

5

inside the wave rotor. The simplified numerical calculation for the various stages for a straight
channel is described below.
The numerical calculation for the wave disk engine is performed using commercially available
computation fluid dynamics solvers to obtain the combustion and port timings and the theoretical
power production. The explicit Navier-Stokes solver module in Ansys Fluent was used to obtain
the initial results for the simplified straight channel geometry. The numerical scheme used in the
following simulation results is discussed in greater detail in Appendix C.
The computational domain consists of a scaled straight channel which measures 50 mm in length
and 10 mm wide. The two dimensional channel is meshed with a structured grid which consists
of square cells measuring 0.1 mm. The schematic representation of the computational domain is
shown in Fig 1.2.

Exhaust
port

Inlet
port

Figure 1.2 - Computational domain of straight channel.
The inlet and outlet ports of the straight channel are located at the left and right ends of the
computational domain respectively. The domain is initialized with a pressure of 7 bar and a
temperature of 2100 K to obtain the conditions after complete combustion. The exhaust port is
opened and the combusted gases are allowed to escape. This creates expansion waves which
reduce the pressure and temperature of the channel representing stage IV of the wave engine
operation.

6

The inlet port is opened by setting it as a pressure inlet boundary when the average total pressure
inside the channel achieves its minimum value. The inlet boundary injects stoichiometric
methane-air mixture at ambient temperature and 10 kPa overpressure. This represents stage I of
the wave engine operation.
Once the channel is completely filled with methane, the exhaust port is closed. This generates a
hammer shock, which propagates towards the inlet. The hammer shock increases the pressure of
the premixed gas by a factor of 1.75 representing stage II of the wave engine operation. The inlet
port is closed once the hammer shock reaches the opposite end of the chamber. The premixed
mixture is then ignited with a spark plug to produce constant volume combustion process,
representing stage III of the wave engine cycle. Once the combustion is complete, the burnt
gases are exhausted in stage IV described earlier.
The average static pressure inside the channel is shown in Fig. 1.3.

Figure 1.3 - Pressure inside the wave engine channel.

7

The time taken for the different stages depends on the length of the channel and its aspect ratio.
For the present study, the time taken for combustion is calculated by utilizing both the theoretical
laminar flame speed of premixed methane-air stoichiometric mixture at 1.75 bar and the
experimentally observed combustion time. The time required for the different stages for the
present channel design using the theoretical laminar flame speed is shown in table 1.
Wave Engine Cycle Required Time (μs)
Stage I

148.3

Stage II

51.1

Stage III

30,000.0

Stage IV

66.0

Table 1.1 – Wave engine cycle timing, theoretical flame speed
An experimental analysis was performed using a channel of the same geometric dimensions to
obtain an estimate of the total time needed for combustion. The combustion chamber was
initially filled with premixed methane-air mixture at 75 kPa of overpressure to obtain the
condition similar to the end of stage II. The total combustion time was observed to be 6,300
microseconds. The shorter time required was due to changing flow-field and temperature of the
unburnt mixture during the combustion process. The data for the experimental time for
combustion are shown in table 1.2.

8

Wave Engine Cycle Required Time (μs)
Stage I

148.3

Stage II

51.1

Stage III

6,300.0

Stage IV

66.0

Table 1.2 – Wave engine cycle timing, experimental flame speed
The comparative time taken for the four stages is shown graphically in Fig. 1.4.

1% 2% 1%

Stage I, Fuel-air mixture injection
Stage II, Hammer shock compression
Stage III, Constant volume combustion
Stage IV, Isentropic expansion
96%

Figure 1.4 - Wave engine cycle time, in percentage.
The time required for combustion introduces restrictions in the speed of rotation of the wave
rotor as well as inlet and exhaust port timings and opening angles.

9

1.2 COMBUSTION CHAMBER DESIGN
The first generation wave rotor was constructed using the research data obtained by a
computational analysis of the one dimensional channel geometry [3]. This enabled the
construction of a wave engine with “C” shaped channels with inlet and exhaust ports optimized
for 20,000 RPM. The schematic representation of the first generation wave rotor is shown in Fig.
1.5. The wave rotor design was tested extensively for various operating conditions and port
timings. The corresponding power production and pressure results were used in the setup and
calibration of the numerical model.

Figure 1.5 - First generation wave disk engine design.

10

The results from the numerical simulation enabled the design of the second generation of the
channel geometry. The turbine efficiency of the wave rotor was numerically computed to be 16%
which prompted modification to the channel geometry. The second generation rotor was
designed to increase the turbine power extraction by the use of a converging channel section near
the exhaust port. This allowed the high pressure combusted gases to increase the speed at exhaust
allowing a greater torque. The second generation rotor design also featured a two cycle operation
which enabled a slower rotational speed by decreasing the time allowed for the various stages.
The design of the channel geometry, intake and exhaust port openings was optimized for a
design rotational speed of 10,000 RPM. The engine assembly was modified to accommodate
direct injection of fuel into the combustion chamber prior to spark ignition in order to enable
stratified mixture combustion. The lean combustion process was used to reduce NOx and SOx
emission while increasing the thermal efficiency of the engine. A schematic representation of the
second generation wave disk engine is shown in Fig. 1.6.

11

Figure 1.6 - Second generation wave disk engine design.
The second generation design incorporated additional features such as a secondary inlet, return
channels for multistage turbine extraction and replaceable exhaust port geometries. The inlet
geometry was modified to enable either the separate injection of air and fuel or the injection of
scavenging air to create a buffer zone between hot exhaust gases and cold fuel-air mixture.

12

The third generation wave engine design featured a modification of the channel geometry to
increase the turbine work extraction of the wave engine. The channel design was modified to
accommodate a converging-diverging geometry to act as a nozzle to accelerate the speed of
ejection of hot gases. The schematic representation of the wave engine is shown in Fig. 1.7.

Figure 1.7 – Third generation wave disk engine design.
The nozzle geometry was optimized to enable the maximum torque generation for an initial
pressure of 3 bar. The size of the inlet chamber was reduced to enable the design of a longer
channel and increase the volume of the combustion chamber. The intake and exhaust port
timings were optimized for a rotation speed of 5000 RPM and a two cycle design feature was
utilized. The construction of the nozzle geometry allowed the exhaust gas to be expelled at
supersonic velocities, further increasing the turbine power generated by the nozzle.
13

The use of the nozzle also allowed an increase in the time taken for the exhaust gases to be
expelled. The contour plot of the Mach number during the exhaust stage of the channel is shown
in Fig. 1.8. The channel is initially filled with burnt gases at high pressure and temperature of 7
bar and 2100K and the exhaust port is opened. This creates a chocked flow at the nozzle and
generates power due to the channel curvature.

Figure 1.8 - Contour of Mach number during exhaust stage.
The nozzle reduces the static pressure of the exhaust gases, which increases the efficiency of the
turbine. The contour plot of the static pressure is shown in Fig 1.9. The efficiency of the turbine
is 86%.

14

Figure 1.9 - Contour of static pressure during exhaust stage, units in kPa-g.
The time required for the four stages is shown in Table 1.3 and Fig 1.10.
Wave Engine Cycle Required Time (μs)
Stage I

341.1

Stage II

117.5

Stage III

3,638.5

Stage IV

1580.6

Table 1.3 – Wave engine cycle timing, third generation rotor

15

6%

2%

28%
Stage I, Fuel-air mixture injection
Stage II, Hammer shock compression
Stage III, Constant volume combustion
Stage IV, Isentropic expansion
64%

Figure 1.10 – Third generation wave engine cycle time, in percentage.
The total time taken for the four stages is equal to 5.7 milliseconds. For a two cycle wave disk
engine, the maximum rotational speed for peak efficiency is computed to be 5300 RPM.
The two biggest challenges encountered in the operation of the wave disk engine were the
efficiency of the turbine work extraction and the time required for constant volume combustion.
The efficiency and work output of the radial turbine is directly proportional to the rotational
speed of the turbine. The combustion process acts as a rate-limiting step which decreases the
work and thermal efficiency of the engine. This prompted the need for a detailed study of
combustion in a constant volume combustion chamber. Improvement in the flame speed would
lead to faster rotational speed of the wave engine, potentially leading to an increase in the
efficiency and power.

16

CHAPTER 2
LAMINAR FLAMES
Combustion in premixed fuel-air mixtures usually consists of a rapidly moving flame-front
which travels into the unburnt mixture with a finite reaction rate. This type of reaction is
typically called deflagration. The flame front propagates towards the unburnt mixture through
the action of species and heat diffusion into the unburnt mixture.
The two main classifications of combustion depending on their speed of propagation are
deflagration and detonation. When the speed of flame propagation is slower than the speed of
sound, it is called deflagration. When the speed of flame propagation is faster than the speed of
sound, it is known as detonation. The propagation of deflagration usually occurs due to the
diffusion of species and thermal energy into the unburnt mixture. Detonation is sustained by the
action of a shock wave travelling just ahead of the flame-front, which causes the reactants to
dissociate and recombine into products, with the heat release of the chemical reaction further
reinforcing the shock.
The part of the reaction zone which releases heat and light is typically called the flame. The
location of the flame-front is dependent on the stoichiometric conditions of the reactants. In
diffusion flames, the flame-front is sustainable only in those places where the fuel-air ratio is
suitable for combustion. However, in premixed flames, the fuel and air are intimately mixed and
the flame-front tends to propagate into the unburnt mixture at a speed proportional to the reaction
rate.

17

Premixed flames undergo numerous instabilities caused by to the shape of the combustion
structure, the Lewis number of the fuel-air mixture and intrinsic effects like the Darrius-Landau
instability mechanism.
The premixed flame propagation process can be achieved in open or closed channels and
burners. The flame structure can be described by establishing a reference frame and defining a
profile for the variables under consideration (temperature, species mass fraction). The coordinate
system is fixed to the moving flame front. In this coordinate system, the unburnt premixed
mixture will be seen as approaching the flame front with a speed equal to the flame speed SL.
The spatial variation of the reactant concentration, temperature and reaction rate in a premixed,
laminar flame front is shown in Fig. 2.1. The expansion ratio of the burnt gas is usually defined
by σ, which is given here by ρu/ρb.

Reaction Rate

Y=1

Tb

Fraction of
Reactants

σS

SL

Temperature
L

Reactants

Products

Tu

Y=0

Preheat Zone
Reaction Zone
Figure 2.1 – Structure of a premixed flame in one dimension [4].
18

The two main zones are defined as preheat and reaction zones, dividing the flame front into two
regions in which diffusion of radicals and heat occurs at different intensity. In the preheat zone,
the amount of heat released is very small and the reaction rate is practically negligible. As the
observer moves closer to the reaction zone, the heat transported from the flame front increases
leading to higher temperature and faster reaction rates. The reaction occurs in a very thin zone
typically about 1mm thick at atmospheric pressure. This region undergoes fast chemical reaction,
and due to high reaction rate and heat release, the concentration gradients are very. The large
thermal and concentration gradients drive the transport of energy and species into the upstream
region know as the preheat zone, shown in Fig. 2.1 and render the flame propagation process
self-sustaining.

19

2.1 LAMINAR AND TURBULENT PREMIXED FLAMES
There are two types of premixed flame depending on the nature of flow of the fluid: laminar and
turbulent flames. The flow regime can be determined by the dimensionless number which is a
function of flow velocity, kinematic viscosity, density and characteristic geometric dimension of
the channel called the Reynolds number.

In laminar flames, the combustion process is

dominated by the chemical kinetics of the mixture and flame propagation occurs through
diffusion of heat and species. In turbulent flames, the velocity and scalar properties are subjected
to rapid, random changes at a given point in space; therefore the combustion process is strongly
affected by it. The fluctuating temperature and species concentration due to turbulence is the
dominant force for flame propagation in turbulent flames.
Combustion processes are usually characterized by highly exothermic chemical reactions, in
which the reaction rate is a strong function of temperature. As mentioned previously, the
combustion of a premixed mixture is strongly dependent on the chemical-kinetics of the mixture,
thus the modeling of a detailed and accurate multi-step and multi-species reaction system
represents a major challenge in the analysis of laminar flames.

20

2.2 LAMINAR FLAME SPEED
A detailed study of the propagation of a flame through an unburnt mixture is made in this
section. First, a definition of laminar flame speed and laminar burning velocity is established. It
is common to find ambiguities in literature in which laminar flame speed and laminar burning
velocities have the same meaning, which is the velocity at which unburnt gas moves through the
combustion wave in the direction normal to the wave surface.
The two flame speeds are defined below:
(i)

The laminar burning velocity is defined as the velocity at which unburnt gases move
through the combustion wave in the direction normal to the wave surface.

(ii)

The laminar flame velocity is the velocity of the flame with respect to a fixed outside
observer.

In the case of a standing or stationary flame, like the flame in a Bunsen burner, the laminar flame
velocity is equal to the laminar burning velocity. For freely propagating flames, the laminar
burning velocity is equal to the sum of the laminar flame speed and the flow velocity of the
unburnt mixture. For example, the laminar flame velocity of methane-air and stoichiometric
mixture at STP is 36.2 cm/s.
The three main flame propagation theories are broadly classified as thermal, diffusion and
comprehensive theories.
The first theory describing flame propagation was the thermal theory of Mallard and LeChatlier
[5], who stated that the controlling mechanism of flame propagation is heat diffusion from the
burnt region to the unburnt region. This formulation defined the flame as a composite region

21

comprising two zones, separated by a layer where ignition occurs. Ignition was postulated to
occur when the local temperature reached the ignition temperature and δ was described as the
reaction zone thickness. This model is shown in Fig. 2.2.
Zone 2 (Preheat Zone)
T
Tb
Tb
Ti
Zone 1 (Reaction Zone)

Tu

Figure 2.2 - Mallard and LeChatelier model [4].
Improvements were made to the thermal energy model by Zeldovich and Frank-Kamenetskii [6].
Their formulation included the diffusion of molecules and heat, but did not involve dissociation
and the diffusion of free radicals. Tanford and Pease [7] developed this theory further, stating
that the mechanism for flame propagation depends only on the diffusion of radicals and not on
heat propagation as proposed by earlier theories.
These two theories of heat and radical diffusion drove the experimental research on the
determination of flame velocity on initial parameters, pressure and temperature in order to
examine the validity of the theoretical models of flame propagation. As is often the case, neither
model was entirely correct.

22

The following observations follow from the two models of flame propagation –


According to the thermal theory, increasing the initial pressure leads to a higher
temperature. This leads to a faster rate of reaction and a higher flame speed in achieved.



According to the diffusion theory, increasing the initial temperature leads to higher
dissociation. The higher concentration of radicals causes greater species diffusion leading
to a higher flame speed.

A comprehensive theory combining the species and thermal diffusion theories was developed by
Hirshfelder et al. [8]. A complex set of non-linear equations describing the diffusion and energy
balance were developed that could be solved by numerical methods. The Arrhenius equation for
the reaction rate, by definition, has a very small non-zero value at regions far upstream of the
flame, where the local temperature is equal to the ambient temperature. This was modified to
enable a zero reaction rate model at ambient temperature and equating the wall surface to act as a
heat sink. The exact solution of laminar flame propagation requires the use of fluid dynamics
equations with modifications done to accommodate heat release and species/thermal diffusion.
The two theories are comprehensively discussed by Irvin Glassman and Richard A. Yetter in
their book “Combustion” [4]. A short description of the two theories is provided in the following
section.

23

2.2.1 Mallard and Le Chatelier Theory
The theory developed by Mallard and Le Chatelier states that the heat conducted from zone 1 in
Fig. 2.2 is equal to that required to raise the temperature of the unburnt mixture to the ignition
temperature. The rate of temperature increase is assumed to be linear and is approximated to
equal [(Tf – Ti)/δ]. Here, Tf is the flame temperature, Ti is the ignition temperature and δ is the
thickness of the reaction zone.
The energy balance between the downstream inwardly convecting flow and the upstream
outwardly conducting thermal energy is given by the simple relation -


mCp (Ti Tu )  (Tf  Ti ) / A

(1)

where λ is the thermal conductivity of the mixture, m is the mass flow rate of the unburned gases
into the combustion wave, Tu is the temperature of the unburnt gases and A is the area of cross
section. The combustion wave is considered to be one dimensional.
From the continuity equation, we have, for the steady mass flow rate, the expression -


m  UA  SLA

(2)

where ρ is the density and U is the velocity of the unburnt gas, which is equal to the laminar
flame speed, SL.

Upon using the mass balance equation with the energy balance, we obtain -

SLCp (Ti  Tu )  (Tf  Ti ) / 

24

(3)

Solving for the flame speed, we obtain:

SL 

(Tf  Ti ) 1
Cp (Ti  Tu ) 

(4)

For a steady flow condition, the total mass rate entering the reaction zone must be equal to the
rate of consumption of the reactants, thus:

m

 U  SL  w
A

(5)

where, w represents the reaction rate in terms of concentration. Using this expression in the
flame speed equation yields:
1

1

  (T f  Ti ) w  2


 w 2
SL  
   
 
 C p (Ti  Tu )  



(6)

where α is the thermal diffusivity. From the expression above, we obtain an expression for the
mass consumed by the unburnt gas is equal to 1

  2

SL   w 
 Cp 



(7)

Thus, the Mallard-Le Chatelier approach yields the following result for the laminar flame speed 1


 w2
SL   
 

25

(8)

Equation (8) enables the estimation of laminar flame propagation speed for changes in the
physical and chemical parameters. We note that in this analysis, the chemical variables appear
implicitly in w.
A similar approach was used by Linan and Williams [9] based on Mikhelson’s studies [10]
where the heat release in the reaction zone is equated to the energy conduction from the hot
products to the cool unburnt gases.
The total energy per unit mass that is conducted to the unburned mixture is:
(9)
And the energy balance is calculated as:


hwL  (Tf  Tu ) / L

(10)

where, δL is the thickness of the reaction zone, which includes both region 1 and region 2 from
the Mallard-Le Chatelier approach.
Combining the energy balance equations gives -


Cp (Tf  Tu )wL  (Tf  Tu ) / L

(11)

Solving the above equation for the reaction zone thickness leads to
1
  1 2

L  


 Cp w 


26

(12)

Using Equation (5) derived earlier, we obtain:
1

1

  w2  w2


SL  
   
 
 Cp  



(13)

The expression derived above is the same as that derived by Mallard-Le Chatelier, but the
derivation here does not employ the temperature ratio between the ignition temperature and
initial/final temperature.
Additionally, following the Zeldovich, Frank-Kamenetskii and Semenov theory, in which the
temperature term is evaluated by assuming asymptotically narrow reaction rate region, we obtain
the following expression for the temperature term:

RTf2
E(Tu  Tf )

(14)

This term was identified as the Zeldovich number -

Ze 

RTf2

 L
E(Tu  Tf )


(15)

The expression for the laminar flame speed can now be derived as:
1


  w2
SL  

 Ze  

(16)

This expression allows a reasonably accurate global analysis of the flame speed with respect to
changes in the physical and chemical parameters of the problem. We obtain that the reaction
term is proportional to the pressure of the reactants as follows:
27

SL

1
n 2 2
 (P
)

(17)

Here, n is the overall order of the reaction.
The laminar flame speed’s dependence on the temperature is mainly dictated by the exponential
in the rate expression w. This gives the expression:
1

SL ~ exp(  E
) 2
RTf 




(18)

We note that the maximum reaction rate and heat release takes place close to the highest
temperature if Arrhenius kinetics as the controlling factor for the flame chemistry.

28

2.2.2 Zeldovich, Frank-Kamenetskii and Semenov Theory
The analytical model of flame propagation suggested by Zeldovich, Frank-Kamenetskii and
Semenov, serves as a more detailed method to determine the flame propagation speed. This
model combines both the effect of thermal diffusion and the diffusion of species into the unburnt
mixture on the flame speed. This model ignores the effect of diffusion of radicals and their effect
on reaction rate.
The main underlying assumption in this model is that the ignition temperature is very close to the
flame temperature. The ignition temperature is defined as the threshold temperature above which
most of the reaction occurs. In addition, the specific heat and thermal diffusivity are both
assumed to be constant over the entire flame region. This enables the use of the non-dimensional
Lewis number, which is the ratio between the thermal and mass diffusion. When the rate of
thermal and mass diffusion is identical, the Lewis number is unity.

Figure 2.3 – Balances across a differential element of laminar flame [4].
29

The flame is assumed to be a one-dimensional wave of thickness Δx as shown in Fig. 2.3. Here,
‘a’ is the mass of reactant per cubic meter, w is the rate of reaction and Q is the heat of reaction
per unit volume. Since this is a steady propagation process, there is no accumulation of species
or heat within the control volume.
The equations of mass balance, mass diffusion and energy balance equation are formulated
within the limits of the unburnt and burned gases from X=-∞ to X=+∞. The combustion wave is
considered to consist of two parts. The unburnt part in which there is no reaction and the reaction
part in which the reaction and diffusion processes dominate and the convective term is neglected.
The energy equation is centered at x=0, which allows the local temperature at the origin to be
specified as the ignition temperature.
Upon integration of the system of equation, the solution for the flame speed is obtained as
follows:
1
 2

 2     E  RT 2
f
 Ae RT 

SL   
 ET  T  
a 0  oCp 
 

f
0 




(19)

Here, A is the pre-exponential factor in the Arrhenius equation and the subscript 0 refers to the
initial conditions.

30

2.3 Flame Speed Measurements
Depending on the method of measurement and the accuracy requirement, the measurement of
flame speed can either be simple or complicated. Flat, adiabatic flames are difficult to obtain in
laboratories due to heat loses and gravity, hence the most common flame shapes are either
conical or spherical. The most common methods used in the calculation of flame speeds are by
Bunsen burner, cylindrical tube, spherical flame and the flat flame burner methods. These are
discussed in detail below.
Bunsen burner method –
In a Bunsen burner with long tubes, the velocity profile at exit is parabolic, leading to variable
burning velocity over the flame surface. For this case, an average method is used to determine
the normal velocity component, Un over the “flame cone” that the Bunsen burner flame
produces. The normal velocity Un is equal to the flame speed SL. It is assumed that Un is
constant over the flame area Af. If m is the mass flow rate of the gas, and ρ is the unburnt gas
density, from mass continuity, we obtain:

Un 


m
u A f

(20)

When suitably shaped nozzles are used, the flow is uniform and the flame has straight edges. The
flame speed can then be computed as -

Un  Uu sin u
The structure of the flame in a Bunsen burner is shown in Fig. 2.4.
31

(21)

Figure 2.4 - Bunsen burner flame [4].
Cylindrical tube method –
This method utilizes a cylindrical open tube, typically made of glass, where the premixed fuel-air
mixture is ignited at the open end. The speed of flame propagation is measured optically. Due to
buoyancy effects, the flame shape is not planar and the area of the flame is larger than the crosssectional area of the tube. The combustion induces pressure waves which travel faster than the
flame and induce localized velocity in the unburnt mixture which can affect the speed of
propagation of the flame with respect to the laboratory coordinates. Additionally, heat loses to
the wall can influence the flame speed as the non-quenched region mixes with the quenched
region of the flame.

32

Spherical flame method –
In this method, the premixed mixture is contained in a soap bubble and is ignited at the center, so
that the flame will propagate radially through the unburnt mixture. Since the gas is enclosed by
the soap film, the pressure remains constant. There is no heat loss to the outside as the soap
bubble contains no heat sinks. The bubble grows as the premixed mixture burns and this growth
is measured optically. The flame speed can be measured as follows:

 
SL  U n  b 
 
 u

(22)

Flat flame burner method –
This method of flame speed measurement is performed by supplying combustible gas mixture to
a series of thin, parallel tubes of diameter not larger than 1 mm. The mixture flow rate is adjusted
to obtain a thin flat flame resting on top of the burner surface. The flame speed is obtained by
dividing the volumetric flow rate with the flame area. Since the flame is perfectly flat, it is the
most accurate method of measuring the flame speed of a fuel-air mixture.

33

2.4 Instabilities in Flame Propagation

Combustion in both open and closed channels is subjected to a number of instabilities which
have been studied theoretically, computationally and experimentally. The three main types of
instabilities which have the greatest impact on flame structure and propagation are listed below:
i.

The Landau-Darrieus instability, an intrinsic instability which occurs due to
hydrodynamic effects.

ii.

Instability caused by the coupling of combustion and acoustics of the system. These
instabilities are characterized by low-frequency oscillations in the longitudinal
dimension of the channel.

iii.

Instability caused due to the coupling between combustion and acoustic modes. These
instabilities are characterized by high frequency oscillations in the transverse
dimension of the channel.

The last two types of instabilities are commonly classified as thermo-acoustic instabilities and
the first one is referred to as a hydrodynamic instability. The three instabilities are discussed in
great detail by G. Searby in the paper “Instability Phenomena during Flame Propagation” [11]
and in Combustion Physics by C. K. Law [13]. A short description is provided in the following
section.

34

2.4.1 Landau-Darrieus Hydrodynamic Instability
All premixed flames are intrinsically unstable due to the hydrodynamic effects, which is caused
due to the expansion of hot, combusted gas into the denser, unburnt mixture. A one-dimensional
representation of an infinitesimally thin laminar flame is shown in Fig. 2.5. The planar flame
front is essentially a thin interface (a sheet) which separates the unburnt and burnt gases. The
cold reactant gases at temperature and density, T0 and ρ0, are separated from the hot combusted
gases of temperature and density, Tb and ρb. The flame front propagates at velocity SL into the
unburnt gas mixture. Upon combustion, the gases undergo thermal expansion and they exit the
flame front at Ub=SL.(σ). The density ratio here is almost equal to the temperature ratio for
standard hydrocarbon fuels, since the pressure is nearly constant.

ρb
ρu

Figure 2.5 – Propagating flame front [11].

35

To conserve the momentum, the velocity jump must be accompanied by a small increase in
pressure which is typically around 1 Pascal, given by:



  



1
1
2
P  . b U 2  o Uo  . 0S2 . o  1
b
L 

2
2
 b 

(23)

When the flame front is inclined with respect to the axis of propagation, as shown in Fig. 2.6, the
propagation of the exhaust gases changes to conserve momentum.

Figure 2.6 – Inclined flame front [11].

The incoming flow of the premixed fuel-air mixture can be decomposed into vector components
which are parallel and normal to the flame front as Ull and Un respectively. If the flame front is
the reference frame, then normal component of the incoming premixed mixture must be equal to
36

the flame velocity, Un=SL. The burnt gas leaving the flame front will have a normal component
equal to Un(ρo/ρb) (thermal expansion). The parallel component of the flow remains unchanged
as there is no physical mechanism to sustain the parallel pressure jump that would be necessary
to accelerate the parallel component of the flow. Upon vector addition of the parallel and normal
components, the resultant vector direction is now deviated to the outgoing normal. This causes
the flow to change its direction after passing through the combustion wave if it is not perfectly
normal to the flame-front sheet.
Consider a flame which is not planar, but is wrinkled at a wavelength λ as shown in Fig. 2.7. At
the localized places where the streamlines are exactly normal to the flame-front, they are
accelerated, but do not deviate, when they cross the flame. At places where the flame-front is
inclined to the incoming mixture, the streamlines are deviated to the rear normal of the flame
front. Although the streamlines in Fig 2.7 are locally correct, streamlines in a flow cannot cross
each other and they must curve to become parallel in the far downstream as shown in Fig. 2.8.

37

Figure 2.7 - Local streamline deviation in a wrinkled flame [11].

Figure 2.8 - Global streamline deviation in a wrinkled flame [11].

38

When the streamlines are curved, there are pressure gradients along the flow. This introduces
perturbations which are not local. This non-locality of the presence of flame-induced pressure
gradients will not only affect the downstream flow, but also the upstream flow. This effect
causes the gas expansion through a curved flame to converge at places where the flame-front is
concave and diverge where the flame-front is convex as shown in Fig. 2.8. Due to mass
conservation, the upstream flow is either accelerated or decelerated where the flame front lags
behind or is ahead of the mean position. Since we have assumed that the propagation velocity is
constant, the flame front becomes unconditionally unstable and the wrinkling will grow with
time.

2.4.2 Thermo-Diffusive Effects
The Zeldovich-Frank-Kamenetskii model of premixed flame propagation has high activation
energy and the chemical reactions are confined to a thin layer on the high-temperature side of the
flame front. The mechanism of flame propagation is primarily affected by the diffusion of the
heat and species within the flame thickness, δ. In the case of a curved or wrinkled flame-front,
the gradients of temperature and species concentration are not parallel to the average direction of
flame propagation and therefore the local flame velocity can change.
At a place where the flame-front is concave towards the unburnt mixture, the heat flux is locally
convergent. The local flame temperature and the local propagation velocity increase as a result.
At a place where the flame-front is convex towards the unburnt gas, the heat flux is locally
divergent. The local flame temperature and the local propagation velocity decrease as a result.
The effect of thermal diffusion is to stabilize the wrinkled flame.
39

The gradient of the thermal diffusion is shown in blue arrows and that of species concentration is
shown in green arrows in Fig. 2.8. The species gradient and the thermal gradient are in the
opposite directions. At a place where the flame-front is concave towards the unburnt mixture, the
species flux is locally divergent. The flux of reactive species into the reactive zone decreases as a
result, leading to a decrease in the local propagation velocity. The effect of species diffusion is to
destabilize a wrinkled flame. The net result of these two diffusive fluxes will depend on the
ration of the thermal and species diffusion coefficients. This ratio is called the Lewis number,

Le = Dth/Dmol

(24)

If the Lewis number is greater than unity, the effect of heat diffusion is prominent and the flame
is thermo-diffusively stable. There is however an additional stabilizing contribution that arises
from the inclination of the streamlines within the preheat zone. This internal inclination of the
streamlines creates an additional transport of heat and species that is convergent or divergent
with respect to the average direction of propagation. It has the effect of contributing an additional
term in the expression for the flame velocity. This term is always stabilizing, independent of the
Lewis number, and it increases with the gas expansion ratio.

40

Figure 2.9 - Wrinkled premixed flame structure [11]. In the concave region of convergent gas
flow, the streamlines move closer together, where the diffusion and mass flux is divergent (and
the upstream heat flow is convergent).

41

2.4.3 Thermo-Acoustic Instabilities
Combustion in constant volume generates acoustic waves that can interact with the flame front
and cause perturbations in the flames. These acoustic waves are usually reflected off of enclosing
walls and internal obstacles. Instability mechanisms that produce such perturbations have been
studied under different experimental conditions. Thermo-acoustic instabilities usually manifest
itself as flame distortion leading to increased flame surface area, increased flame propagation
speed and peak pressure enhancement. Lining of the combustion chamber walls with materials
that can absorb acoustic waves, or using obstacles to deflect the acoustic waves can help to
reduce the intensity of these acoustic instabilities.
Unsteady combustion is a strong source of acoustic noise, whose emission is governed by the
conservation equations given by:
Mass conservation
D
 .V  0
Dt

(25)

DV
 p
Dt

(26)

Momentum conservation



Energy Conservation

Cp

DT
D

q
 .(T)
Dt
Dt

42

(27)

Equation of state





P
C
 Cp  C v T  v a 2

Cp

(28)

where D()/Dt is the material derivative, V is the gas velocity vector and a is the local speed of
sound.
Combining the energy conservation and the equation of state, we obtain:

Dp
Dp  Cp  Cv 
.q

 a2

Dt
Dt  Cv 



(29)

When an acoustic wave interacts with a flame-front, the heat released due to reaction will be a
function of the acoustic pressure magnitude. This phenomenon was described by Rayleigh [12]
as follows – “If heat is periodically supplied and removed from a mass of vibrating air, the effect
produced will depend on the phase of vibration. If the heat is supplied at the moment of highest
condensation and removed at the moment of highest rarefaction, the application of heat would
lead to a support of the vibration.”
The Rayleigh criterion states that positive energy is transferred to the acoustic wave if pressure
and heat release fluctuations are in phase and negative energy is transferred to the acoustic wave
if pressure and heat release functions are out of phase.
The two major types of coupling between acoustic wave and the combustion process are pressure
coupling and acceleration coupling, which are discussed in greater detail below.

43

Pressure Coupling
The standard form of the chemical reaction rate in the Zeldovich-Frank-Kamenetskii model [6]
using Arrhenius law is described below:



Ea
AY RT

(30)

Where A is the pre-exponential factor, Y is the mass fraction of the limiting reactant, Ea is the
activation energy and ρ is the reactant density.
According to the above expression, the reaction rate is directly proportional to the density of the
reactants. The density of the reactant is proportional to its pressure if the reactants obey ideal gas
law, this leads to variation in the reaction rate with changes in density caused due to pressure
waves. The propagation of acoustic waves also causes an increase in the local temperature due to
its adiabatic nature. For reactions with high activation energies, the heat release is primarily
sensitive to the temperature oscillation produced by the acoustic waves.
The coupling of acoustic pressure waves and the reaction depends on the acoustic frequency,
angular velocity and the Lewis number of the reactant. The use of one-step Arrhenius equation to
model the reaction leads to a mismatch between the theoretical result and experimental
investigations performed on planar methane-air flame subjected to acoustic oscillation.

44

Acceleration Coupling
The local velocity field in an acoustic wave imparts a localized acceleration on the flame front,
which leads to changes in the reaction rate. The flame-front is also influenced by other
acceleration fields such as gravity due to the differences in the densities of the fluids separated
by a flame-front.
When the localized velocity of the pressure wave is oriented towards the unburnt mixture, the
size of the reaction cell increases and when the localized velocity of the pressure wave is
oriented towards the burnt mixture, the size of the reaction cell decreases. The total reaction rate
of the combustion process is directly proportional to the size of the reaction cell, hence the
acceleration coupling between the flame-front and pressure wave has a major effect of the
reaction rate of the flame.

45

CHAPTER 3
EXPERIMENTAL SETUP, TEST METHODOLOGY AND DATA ACQUISITION
The experimental setup consists of a constant volume combustion chamber with a square cross
section, high speed camera to capture the ignition and propagation of the flame, pressure
transducer and ignition control circuit to record the pressure and electronically control the
ignition time and a Schlieren diagnostic setup to visualize the density change during combustion.
The pressure transducer and the ignition circuit are connected to a computer through data
acquisition and control hardware. LabView program interface is used to initiate the spark and
record the pressure inside the combustion chamber.

46

3.1 COMBUSTION CHAMBER SETUP
The combustion chamber apparatus consists of three metal parts, the center plate and two side
plates which are bolted together to form a constant volume combustion chamber. The enclosed
combustion chamber is 12 inches long and has a square cross section with each side measuring 2
inches. The total volume of the combustion chamber is 48 cubic inches. An exploded view of the
apparatus is shown in Fig. 3.1. Its dimensions are specified in Fig. 3.5.

Center plate
Side plates
Fig 3.1 – Exploded view of combustion chamber apparatus.

47

The combustion chamber apparatus is designed to withstand a maximum pressure of 35 bar and a
pressure relief valve rated at 20 bar is mounted to ensure pressure remains within the safe limits.
The isometric view of the combustion chamber model is shown in Fig. 3.2. The manufactured
apparatus is shown in Fig. 3.3. Teflon tape is used to form an air-tight seal between the pressure
transducer, spark plugs, quick-disconnect valves and the combustion chamber apparatus. Hightemperature transparent silicone sealant is used to attach the glass window to the side plate of the
apparatus.

Fig 3.2 – Isometric view of combustion chamber apparatus.

48

Fig 3.3 – Manufactured combustion chamber apparatus.

The center plate made of aluminum measures 16.5x6.5x2 inches. A hollow channel measuring
12x2x2 inches is machined into the plate which forms the combustion channel. The center plate
is machined further to mount the flush-diaphragm pressure transducer, spark plugs, and inlet and
exhaust ports. The inlet and exhaust ports are used to remove the combusted gases and fill the
combustion chamber with fresh fuel-air mixture. The ports have a diameter of 0.125 inches,
located near both ends of the combustion chamber to avoid disrupting the flame structure. The
isometric and orthographic projection of the center plate is shown in Fig. 3.4 and Fig. 3.5.

49

Fig 3.4 – Isometric view of center plate.

Fig 3.5 – Orthographic view of center plate. Dimensions are in inches.

The center plate is sandwiched between two side-plates. The side plate holds a glass window that
forms the side wall of the combustion chamber. The side plate is made of steel and measures
16.5x6.5x1 inches. It has a recessed opening measuring 12x3x0.75 inches with a 0.25 inch lip
along the sides to hold the glass window as shown in Fig. 3.6 and Fig. 3.7. A borosilicate glass
plate measuring 12x3x0.75 inches is mounted on the recessed opening and forms the side-wall of
the combustion chamber. The glass plate is attached to the side plate using Loctite hightemperature silicone sealant which forms an airtight seal. The side plate and center plate are
bolted together with brass bolts of 0.5 inch diameter. An airtight seal is formed between the side
plate and center plate by lining the contact surface with Loctite silicone gasket. The isometric
and orthographic view of the side plate is shown in Fig. 3.6 and Fig. 3.7.

50

Recessed lip to hold glass plate
Fig 3.6 – Isometric view of the side plate.
16.5
12.0
0.25

Φ0.50

Recessed lip to hold glass plate

3.0

6.5

Fig 3.7 – Orthographic view of the side plate. Dimensions are in inches.

To maintain a constant volume during combustion and avoid leakage from the inlet and exhaust
ports during rapid pressure changes, Swagelok

TM

quick-connect couplings (BQC4D2PM) are

mounted on the inlet and exhaust ports. The quick-connect coupling has a shut-off valve with a
pressure rating of 25 bar. An Omegadyne PX43EO-100GI flush-diaphragm pressure transducer
51

is mounted on the center plate for pressure measurement. The transducer is positioned to ensure
the diaphragm is flush with the upper wall of the combustion chamber and Teflon tape is used to
seal the screw interface. Two MSD 10mm spark-plugs are flush-mounted on both ends of the
combustion chamber along its length. The spark plugs are powered by MSD Blaster 2 ignition
coils which supply 10mJ of energy for every spark discharge. An optional high-energy discharge
module can be plugged in to increase the energy deposited to 50mJ. The pressure transducer is
connected to a National Instruments NI9203 data acquisition module. The ignition coils are
attached to a National Instruments NI9841 digital relay module. The acquisition and relay
modules are mounted on a National Instruments NI cDAQ-9172 data controller module which is
connected to the computer. LabView is used to control the firing of the spark plugs and to record
the pressure inside the combustion chamber. The pressure reading from the transducer is sampled
at 10 kHz. A screenshot of the LabView interface is shown in Fig. 3.8. The flame propagation is
optically recorded using a Photron SA5 high speed camera which is focused onto the glass
window of the apparatus. The camera is configured to capture the combustion event at 20,000
FPS frame rate.

52

Fig 3.8 – LabView interface.

53

3.2 SCHLIEREN SETUP
To visualize the change in the density across the flame front and also study the structure of the
flamefront, a Schlieren imaging system was constructed. The system consists of two spherical
concave mirrors 6 inches in diameter with focal length 60 inches. A pin-hole light source forms a
collimated beam of light along the test area. The light beam is focused onto the lens of the high
speed camera after passing through a knife edge. A schematic arrangement of the Z-type
schlieren imaging setup is shown in Fig. 3.9. The theory of schlieren imaging is discussed in
Appendix B.

Fig 3.9 – Schlieren imaging setup.

54

3.3 TEST PROCEDURE
The test procedure involves filling the combustion chamber with a fuel-air mixture. Upon
equilibration, the mixture is ignited with a spark plug. The pressure data along with optical
imaging of the flame propagation is recorded electronically. The combustion chamber is initially
filled with compressed air and tested for leakages every time it is disassembled and reassembled.
The procedure used in performing a test is as follows –
1. The combustion chamber is assembled to the required configuration and tested for leaks
using compressed air.
2. The exhaust valve is closed and fuel is injected into the combustion chamber to the
required equivalence ratio.
3. The fuel-air mixture is allowed to settle inside the combustion chamber for a period of 30
minutes to obtain a fully premixed mixture.
4. The mixture is ignited using at least one spark plug. The pressure and a video of the
flame propagation is recorded and stored electronically.
5. The combusted gases are removed by flushing the system with compressed air for a
period of 10 minutes and the cycle is repeated for a new set of experimental conditions.
The waiting time for the diffusion of fuel in air to obtain a fully premixed mixture and the time
required to remove combusted gases from the combustion chamber by flushing compressed air
through the system was determined after performing tests for varying mixing and flushing times
and empirically ascertaining the minimum time required for repeatable test results.

55

The amount of fuel that needs to be injected to obtain the required equivalence ratio is obtained
using Dalton’s law of partial pressure, which states that the total pressure exerted by the mixture
of non-reactive gases is equal to the sum of the partial pressures of individual gases.
Pi
Ptotal

 xi

(31)

Using this equation, we obtain the final pressure of the combustion chamber after filling it with
fuel. Depending on the fuel used, the initial air pressure is either maintained at ambient pressure
or a partial vacuum is created in the combustion chamber to maintain a final total pressure that is
close to atmospheric pressure.
The optical measurements were conducted by focusing the camera on the side of the combustion
chamber. For optical measurement of flame propagation, one of the glass plates of the
combustion chamber is closed to produce a dark background for the flame. When using the
Schlieren technique, collimated light is allowed to pass through either side of the glass plates to
form an image of the flame structure. The camera records the experiments at either 5000 or
20,000 frames per second (FPS) depending on the flame speeds and test conditions. In order to
reduce noise and capture an accurate structure of the flame, background light is eliminated by
running the tests in a dark room. The pressure data are recorded at 10 kHz and stored along with
the ignition data. Once the experiment is completed, the recorded images are post processed in
MATLAB to obtain the instantaneous position of the flame along with the flame speed. The
MATLAB program used to post-process the recorded images is described in appendix A.

56

CHAPTER 4
EXPERIMENTAL RESULTS
A stoichiometric mixture of methane and air is used as the combustible mixture in all the
experimental tests performed. This is done to ensure that the results of different combustion
chamber configuration and ignition conditions can be compared with each other. Based on
Dalton’s law of partial pressure, for stoichiometric methane-air mixtures, the final pressure of the
combustion chamber is calculated to be 111.96 kPa. The initial condition is obtained by flushing
the system with compressed air. Methane is pumped into the combustion chamber at ambient
pressure until the final pressure is 111.96 kPa. The inlet and exhaust valves are closed. The
mixture is allowed to settle for 30 minutes to allow complete mixing to form a fully premixed
fuel-air mixture.
The structure and propagation of flames in a variety of test conditions was studied. A detailed
analysis for the structure of a flame in a long, straight channel is discussed in Chapter 6. In order
to maintain consistency between different test results, the pressure transducers are recalibrated
every 10 cycles of usage. To ensure an accurate mixture ratio inside the combustion chamber, an
Omega PX43EO-25GI pressure transducer, with an accuracy of 0.1% over the range of 25 psi, is
used to measure the amount of methane injected into the combustion chamber.

57

4.1 SINGLE SPARK IGNITION
In the single spark ignition setup, the combustion chamber is filled with stoichiometric methane
air mixture at an initial pressure of 1 bar. The premixed fuel air mixture is ignited by a single
spark source. A record of the pressure and flame propagation is made. The video is postprocessed in MATLAB to obtain the position of the flame front for every frame.
Upon ignition by a spark plug, a growing spherical flame kernel is initially formed. The flamefront grows faster along the length of the tube when compared to the width, which leads to the
elongation of the flame. The flame area grows rapidly to form a “finger” shaped flame front.
When the sides of the finger flame touch the cold walls, the side skirt is extinguished, the surface
area of the flame rapidly decreases and it progressively flattens and takes the shape of a planar
flat flame. The planar flame subsequently changes its shape to eventually form an inward
pointing cusp near its center; this shape is has been called the “tulip-flame”. The flame front
maintains its tulip shape as it propagates along the length of the tube until it is extinguished by
the cold wall on the opposite end of the tube. Small distortions to the tulip flame are caused by
the influence of pressure waves and gravity. The optical record of the flame propagation is
shown in Fig. 4.1 for every 10 milliseconds time intervals from ignition until extinction.

58

t=0 ms, Ignition.

t=10 ms, Spherical flame.

t=20ms, Finger-flame formation.

t=30 ms, Planar flame.

t=40 ms, Formation of tulip flame.

t=50 ms, Propagation of tulip flame.

t=60 ms

t=70 ms
Figure 4.1 – Flame propagation along the tube, showing the evolution of the flame structure

59

Figure 4.1 (cont’d)

t=80 ms

t=90 ms

t=100 ms

t=110 ms, Collision with cold wall

Upon ignition, the spherical flame-front grows rapidly and 10 milliseconds after ignition, a
finger flame develops. At t=30 milliseconds, the finger flame skirt is quenched by the cold wall
and the flame surface area rapidly decreases. This causes the flame to slow down and the flame
front collapses to form an inward pointing cusp called the “tulip flame” at t=40 milliseconds. The
flame front eventually collides with the cold wall at the opposite end of the tube at t=120
milliseconds while preserving its tulip shape. The schlieren imaging of the ignition event is
shown in Fig. 4.2.

60

t=0 ms

t=2.5 ms

t=5 ms

t=7.5 ms

t=10 ms

t=12.5 ms

t=15 ms

t=17.5 ms

t=20 ms

t=22.5 ms

t=25 ms

t=27.5 ms
Figure 4.2 – Schlieren image of ignition.

We observe that ignition by the spark plug is followed by a rapidly expanding spherical flame
front, which can be observed from 2-10 milliseconds. The tip of the flame travels faster than the
sides of the flame. This leads to the formation of an elongated “finger flame”, which can be seen
at 15 milliseconds in Fig. 4.2. The sides of the finger flame are quenched by the cold wall and
cause the flame surface area to decrease rapidly. A planar flame shape occurs at 27.5
milliseconds.
61

t=30 ms

t=32 ms

t=34 ms

t=36 ms

t=38 ms

t=40 ms

t=42 ms

t=44 ms

t=46 ms

t=48 ms

t=50 ms

t=52 ms

t=54 ms
t=56 ms
Figure 4.3 – Schlieren image of tulip formation.
The transition from a planar flame front to a tulip flame takes place in 10 milliseconds. At t=38
milliseconds, we observe the 3-D structure of the tulip, which is an inward pointing cusp at the
center from which four “petals” emerge. The two vertical lines near the outer edge of the lobe are
62

the flame close to the glass walls. At t=50 milliseconds, the tulip shape undergoes distortion to
its structure from thermo-acoustic instabilities and the influence of gravity. This causes the
formation of smaller cells to grow along the lobes of the tulip flame. The effect of gravity driven
buoyancy causes the upper lobe of the tulip flame to travel a little faster than the lower lobe.
The normalized position of the flame front is shown in Fig. 4.4. The various stages of flame font
structure are marked. Stage I indicates the ignition event followed by the growth of the spherical
flame. Stage II begins with the formation of the finger flame. Stage III indicates the transition
from a finger flame to a planar flame. The planar flame which subsequently changes to a tulip
flame is shown in stage IV.

I

II

III

IV

Figure 4.4 - Flame position vs time, showing the four stages of flame propagation (I=ignition;
II=finger flame; III=transition to planar flame; IV=propagating tulip flame).

63

The normalized flame position data shows the mean position of the flame front for each instant
of time. The flame position changes rapidly from 0 to 20 milliseconds, after which the rate of
change in the mean flame position is slower. The flame speed data are obtained by calculating
the time-derivative of the flame position as shown in Fig. 4.5.

I

II

III

IV

Figure 4.5 - Flame speed vs time.

The flame speed plot in Fig. 4.5 indicates that the fastest physical flame speed occurs shortly
after ignition. The spherical flame kernel accelerates towards the unburnt mixture. Its surface
area increases correspondingly. At t=15 ms, the surface area of the flame-front reaches its
maximum. The speed of the flamefront starts to decrease as the flame shape becomes planar. At
t=25 ms, the flame transforms to a tulip shape. The speed of propagation remains approximately
constant until it is quenched by the cold wall on the opposite end of the tube at t=118 ms.
64

The instantaneous flame speed for the normalized position on the tube is shown in Fig. 4.6. It is
observed that the transformation of the flame front from planar to tulip shape takes place at
roughly half the length of the tube. The speed of propagation of the flame front is approximately
constant after the tulip shaped flame front has formed.

I

II

III

Figure 4.6 - Flame speed vs normalized position.

65

IV

The corresponding pressure plot for the single spark combustion is shown in Fig. 4.7. The
pressure inside the combustion chamber reaches its maximum value about 15 ms before flame
extinction. This occurs due to the loss of heat from the burnt mixture to the cold aluminum wall
which causes a decrease in temperature in the burnt mixture with increase in time.

I

II

III

IV
Flame Extinction

Figure 4.7 - Pressures vs time plot.

66

4.2 Dual Spark Ignition
The dual spark ignition setup consists of two spark plugs on the opposite ends of the tube along
its length. This allows the effect of ignition from multiple sources to be studied in more detail.
Both spark plugs are flush mounted on the combustion chamber. The exact time of the spark
discharge can be electronically controlled through a delay circuit. Experiments are conducted for
spark delay timings of 0, 25 and 50 milliseconds. The flame front from the initial spark is
affected differently based on the stage of propagation. The ignition from the opposite end has a
bigger influence when the tulip shape is fully developed than when the initial flame is still a
spherical or planar flame. The experiments are discussed in greater detail in the following subsections.

67

4.2.1 Dual Spark, 0 millisecond delay
In this experiment, the two spark plugs are fired at the same time, and two separate flame fronts
are formed which travel towards each other. The amount of energy deposited by each spark is the
same as that released by a single spark. The optical record of the flame propagation is shown in
Fig. 4.8. We observe that the different stages of flame propagation are preserved when compared
to a single spark ignition. Shortly after ignition, a spherical flame front is formed at each end of
the combustion chamber near the spark plug. The accelerating spherical flame fronts change to
planar flame fronts after which the tulip shape is formed. The time taken for complete
combustion is approximately equal to half the time taken by a single spark.

t=10 ms, spherical flame

t=20 ms, finger flame formation

t=30 ms

t=40 ms, planar flame
Figure 4.8 – Flame propagation along tube, dual spark with 0 ms delay.

68

Figure 4.8 (cont’d)

t=50 ms, tulip flame formation

t=60 ms

t=70 ms

The normalized position of the flame fronts is shown in Fig. 4.9. The flame travelling from the
left end of the tube to the right end is called flame L-R and the flame travelling in the opposite
direction is called flame R-L. We observe from Fig. 4.9 that the traces for both the flame fronts
follow a similar profile as they propagate from each end of the tube to the middle. The flame
fronts collide with each other at half the length of the tube at the normalized position of 0.5.
Flame R-L ignition takes place at X=1 in the normalized scale and flame L-R ignition takes place
at X=0 in the normalized scale.

69

I

II

III

IV

Figure 4.9 – Flame position vs. time.

The instantaneous speed of flame propagation can be calculated by calculating the time
derivative of the flame position and is shown in Fig 4.9.
I

II

III

IV

Figure 4.10 – Flame speed vs time.
70

It is observed that the flame speed profile is similar to the single spark combustion process. The
ignition event is followed by an accelerating spherical flame on both ends of the combustion
chamber. The spherical flame undergoes elongation along its length to form a finger shaped
flame front. The side skirts of the finger flame are then quenched by the cold walls to create a
planar flame front. This causes a decrease in the flame surface area and also reduces the speed of
propagation. The planar flame then transforms to a tulip structure, which preserves its shape and
speed until the two flames collide with each other. An oscillating behavior is observed in the
speed of propagation of the tulip flame caused by small changes in the shape of the tulip as it
propagates along the combustion chamber. The speed of flame propagation with the normalized
tube length is shown in Fig 4.11.

Figure 4.11 – Flame speed vs normalized flame position.

71

Fig. 4.8 shows that the planar flame changes its shape to the tulip flame at approximately onethird of the distance from either ends of the tube. The tulip flames then travel at a constant
velocity until they collide with each other at the geometric center of the tube.
The pressure trace of the combustion process is shown in Fig. 4.12.
I

II

III

IV

Figure 4.12 – Pressure vs time.
We observe that the peak pressure precedes flame collision. The rate of heat loses to the cold
walls is maximum before flame collision which reduces the peak pressure. The pressure increase
becomes dominant after the formation of the tulip flame indicating an increase in the burn rate of
the premixed mixture once the tulip flame is formed.

72

4.2.2 Dual Spark, 25 milliseconds delay
In this experiment, the two spark plugs are fired with a delay of 25 milliseconds. The spark plug
on the left side of the combustion chamber is fired first followed by the one on the right side. The
optical record of the flame propagation is shown in Fig. 4.13.

t=0 ms

t=10 ms

t=20 ms

t=30 ms

t=40 ms
Figure 4.13 – Flame propagation along tube, dual spark with 25 ms delay.

73

Figure 4.13 (cont’d)

t=50 ms

t=60 ms

From the above figure, we observe that the ignition even at the left side of the combustion
chamber is followed by a growing spherical flame front. While the flame front is in the process
of transforming from planar to tulip shape, at 25 milliseconds, the second spark plug is ignited.
The flame ignited by the second spark follows the structural changes undergone by the first
flame front.
Ignition I occurs at the normalized position X=0 and ignition II occurs at the normalized position
X=1, 25 milliseconds after the initial ignition. The delay in the firing of the spark plug is set and
controlled through a digital delay circuit.
The instantaneous position of the flame front is shown in Fig. 4.14. The total time required for
complete combustion is observed to be 70 milliseconds.

74

Figure 4.14 – Normalized flame position vs time.

The speed of propagation of the flame front can be determined by calculating the time derivative
of the flame position data. The flame speed is shown in Fig. 4.15. Upon ignition of the first
spark, the flame front follows the structure and speed profile of a single spark combustion flame.
When the flame front is in the process of changing from planar to tulip shape, the second spark
plug is ignited. This causes the initial flame front to attain a negative flame speed in the
laboratory coordinate. A sinusoidal profile of the flame speed curves is obtained with a phase
difference as the flames undergo transformation in their shape and their physical flame speed
changes.

75

Figure 4.15 – Flame speed vs time.
The flame speed profile for the normalized position along the tube length is shown in Fig. 4.16.

Figure 4.16 – Flame speed vs normalized flame position.
76

4.2.3 Dual Spark, 50 milliseconds delay
In this experimental test, the delay between the firing of the two sparks is set at 50 milliseconds.
The spark on the left side of the tube is ignited first followed by the spark on the right side. The
firing of the two spark plugs is controlled by a digital delay circuit. The optical record of the
flame propagation is shown in Fig. 4.17.

t=0 ms, Ignition of spark I

t=10 ms, Finger flame formation

t=20 ms

t=30 ms, Planar flame

t=40 ms, Tulip flame formation
Figure 4.17 – Flame propagation along tube, dual spark with 25 ms delay.

77

Figure 4.17 (cont’d)

t=50 ms, Ignition of spark II

t=60 ms, Finger flame formation

t=70 ms, Planar flame

t=80 ms, Tulip flame formation

t=90 ms, Flame collision

From Fig. 4.17, we observe that the flame formed due to the initial spark ignition follows the
profile of the single spark ignition test condition. At t=50 milliseconds, the second spark plug is
ignited. This causes distortion in the structure of the tulip shape in the initial flame front
travelling from left to right. The R-L flame front undergoes the 4 stages of flame structure and
the two flame fronts collide with each other at 90 milliseconds from ignition.

78

The normalized position of the flame front is determined by post-processing the optical record of
the flame propagation and is shown in Fig. 4.18.

Figure 4.18 – Normalized flame position vs time.

The speed of propagation of the flame front can be obtained by calculating the time derivative of
the flame position. The instantaneous speed of flame propagation is shown in Fig 4.19.
We observe that the second ignition event at 50 ms causes the flame L-R to momentarily move
in reverse and causes a negative flame speed to occur in the laboratory coordinate. The speed of
flame L-R becomes negative again when the flame R-L transforms from a spherical to a finger
flame. A reversal in the speed of propagation of flame R-L occurs when it transforms from a
planar flame to a tulip flame. The R-L flame becomes stationary and L-R flame speed increases
once again until the two flame fronts collide. This alternating flame speed profile is similar to
that observed in ignition delay of 25 ms.

79

Figure 4.19 – Flame speed vs time.

The speed of flame propagation with the normalized flame position is shown in Fig. 4.20. It is
observed from the figure that the L-R flame undergoes reversal in its propagation direction
which is seen as a cyclic flame speed curve near the middle of the combustion tube.

80

Figure 4.20 – Flame speed vs normalized flame position.

81

CHAPTER 5
NUMERICAL RESULTS
A numerical simulation of the constant volume combustion process was performed to obtain the
flow field inside the combustion chamber for different stages in the combustion process and the
influence of pressure waves on the flame front. A detailed study of the formation of the tulip
flame will be performed.
A simplified 2 dimensional representation of the combustion chamber is constructed measuring
12 inches long and 2 inches wide. The computational domain consists of structured quadrilateral
cells, each with an edge of 0.1 millimeters as shown in Fig. 5.1 and Fig. 5.2. The domain
contains a total of 1.54 million elements.

Spark regions

Figure 5.1 - Combustion chamber computational domain.

82

Spark region magnified

Figure 5.2 - Combustion chamber mesh, spark region; square cell measuring 0.1 mm.
The commercial computational fluid dynamics software “Ansys Fluent” is used to perform the
numerical calculation. A Navier-Stokes solver is used to calculate the transient results using an
explicit scheme with a fourth-order Runge-Kutta solver. The Courant-Friedrich-Levy number is
set at 0.9 for the entire simulation domain. The solver methodology, boundary and operating
conditions are discussed in appendix C.
The control volume is initialized with a stoichiometric mixture of methane and air. A single-step
Arrhenius reaction rate mechanism is used. The walls surrounding the control volume are

83

initialized to be a no-slip boundary with constant temperature of 293 K. The combustion process
is initialized by patching a localized spot near the physical spark location to 2100 K.
The propagation of the flame from ignition until the formation of the finger flame is shown in
Fig. 5.3 as a contour plot of the temperature.

t=0 ms, Ignition

t=2 ms

t=4 ms

t=6 ms

t=8 ms, Spherical flame formation
Figure 5.3 – Flame propagation along the combustion tube, showing the evolution of the flame
structure from ignition to finger flame formation.
84

Figure 5.3 (cont’d)

t=10 ms

t=10 ms

t=12 ms, Finger flame formation

t=14 ms

t=16 ms, Maximum elongation of finger flame

t=18 ms, Collision of flame skirt with cold walls
The contour plot of the temperature representing the propagation of flame from the elongated
finger flame structure to the tulip shape is shown in Fig. 5.4.

85

t=20 ms

t=22 ms

t=24 ms, Note the flattening of the flame.

t=26 ms, Planar flamefront.

t=28 ms

t=30 ms, Tulip flame formation.
Figure 5.4 – Flame propagation along the combustion tube, showing the evolution of the flame
structure from planar flame to tulip flame propagation.

86

Figure 5.4 (cont’d)

t=32 ms

t=34 ms

t=36 ms

t=40 ms

t=45 ms

t=50 ms, Final tulip shape.

87

The formation and propagation of the tulip flame in both the spatial and temporal dimensions has
been captured with a numerical scheme using a single-step Arrhenius chemical reaction process
and laminar fluid flow model. The result obtained by uncoupling the pressure waves from the
Navier-Stokes (N-S) equations is comparable to that obtained by solving the coupled N-S
equations. We infer that the trigger for the formation of tulip flame is not dependent on the
influence of pressure waves. A detailed analysis of the flow field of the burnt gas during the
expansion and subsequent quenching of the finger flame skirts by the cold wall shows a reversal
in the direction of the local velocity field of the burnt gases as the flamefront nears the cold wall.
This induces a recirculation to be produced near the finger flame skirts. A moving stagnation
point along the centerline is observed, travelling through the burnt mixture towards the
flamefront as marked in Fig. 5.5. The velocity components are shown in the laboratory reference
frame. The local stagnation points are indicated by a yellow circle.

t=16 ms
Figure 5.5 – Velocity vector of flow inside combustion tube, finger flame extinction.

88

Figure 5.5 (cont’d)

t=17 ms

t=18 ms

t=19 ms

89

Figure 5.5 (cont’d)

t=20 ms

t=21 ms

As the stagnation point moves closer to the flamefront, a “stagnation line” develops which
eventually merges with the reaction front. As observed in laboratory coordinates, the burnt gases
move towards and the unburnt mixture moves away from the flamefront during the initial stages
of flame propagation (spherical and finger flame). The speed of the flamefront during the first
two stages is higher than the local velocity field, when observed from the flamefront coordinates,
the unburnt gas moves towards the flamefront and the burnt gas moves away from the
flamefront. As the speed of the flame decreases, there is a reversal in the direction of the burnt
gas and a subsequent reduction in the magnitude of velocity in the unburnt mixture as shown in
Fig. 5.6.

90

t=22 ms

t=23 ms

t=24 ms

t=25 ms

t=26 ms

t=27 ms

t=28 ms

t=29 ms

t=30 ms

t=31 ms
Figure 5.6 – Velocity vector of flow inside combustion tube, tulip formation.

91

It is observed that a recirculation of the burnt gas is produced when the flame skirt of the finger
flame collides with the cold wall. This recirculation moves towards the leading flamefront as its
speed decreases. For tulip flames in short tubes, it has been theorized by Dunn-Rankin that “the
recirculation produces the initial trigger for a Darrieus-Landau instability which subsequently
grows to the entire tulip flame” [5].
The properties of the unburnt mixture vary continuously during the combustion process. The
flamefront typically acts as a “piston” by compressing the unburnt mixture as the flame front
travels towards the other end of the tube. The data for the unburnt mixture properties is collected
at a point on the combustion chamber which is located 1 inch away from the cold wall of the
opposite end of the ignition point, along the centerline. There is a separation of an inch between
the data point and the cold wall to eliminate the effect of conductive heat loss. The data point on
the computational domain is marked with a yellow circle as shown in Fig. 5.7.

Figure 5.7 - Unburnt mixture data point.
The change in the pressure of the unburnt mixture during the combustion process is shown in
Fig. 5.8. It is observed that the effect of pressure waves produced by combustion is negligible
and pressure changes can be assumed to be uniform throughout the control volume.

92

Figure 5.8 - Pressure vs. time in the unburnt mixture.
The increase in pressure is accompanied with an increase in the temperature of the unburnt
mixture which increases by 180 K from ignition start to the point where the flame surface
reaches the data point. The change in temperature of the unburnt mixture is shown in Fig. 5.9 and
the comparison of the change in pressure for a constant specific heat ratio and variable heat ratio
is shown in Fig. 5.10.

93

Figure 5.9 - Temperature vs. time in the unburnt mixture.

Figure 5.10 - Pressure comparison in unburnt mixture.

94

It is observed that the adiabatic pressure and numerical pressure are similar for temperatures
close to ambient conditions and the two pressures diverge at higher temperatures due to the
change in the value of the specific heat ratio. The Lewis number remains unchanged when the
flamefront is far from the unburnt mixture. The change in temperature and pressure has little
effect on the Lewis number and the reaction proceeds with a small dominance of species
diffusivity as shown in Fig. 5.11. This leads to the formation and propagation of thermodiffusive instabilities in the flame front, which is observed both numerically and experimentally.

Figure 5.11 - Lewis number of unburnt mixture.

95

CHAPTER 6
RESULTS AND DISCUSSION
The results from the numerical and experimental and schlieren tests are used to develop a better
understanding behind each stage of flame propagation in a constant volume combustion
chamber.
An analytical model of the flame propagation is constructed where the flame front is
approximated to basic geometric shapes with which approximate value of flame surface area and
volume of burnt and unburnt regions can be calculated.
For a flame ignited by a single spark, the flame propagation process is divided into 4 stages. In
the first stage, spark ignition is followed by a growing hemispherical flame front. The maximum
radius attained by the flame front is visually identified to be three-fourths of the side of
combustion chamber. In the second stage, the hemispherical flame front undergoes an elongation
process along the length of the combustion chamber, which resembles a finger flame. The finger
flame is approximated to be a combination of a cylinder which expands along its radius and
length with a hemispherical cap. Upon reaching the normalized length L=0.45, the finger flame’s
side skirts are extinguished by the cold wall and a planar flame front is formed, which forms the
third stage. In the final stage of flame propagation, the planar flame begins to develop an inward
pointing cusp, which is approximated to be a pyramid with base equal to the area of cross section
and height visually determined to be half the length of its side. This shape is maintained until the
flame front collides with the opposite wall.

96

For the first 10 milliseconds after spark plug ignition, the flame front assumes the shape of a
uniformly growing hemisphere. The surface area of the flame front grows at a square rate of its
radius growth. In this stage, the effect of lateral walls and recirculation of unburnt gas is
negligible and the flame front expansion takes the shape of a physically unconfined premixed
flame. The analytical representation is shown in Fig. 6.1 which consists of a hemisphere of
radius ‘r’ which grows from 0-0.75 R, where R is the hydraulic radius of cross section of the
tube.

Figure 6.1 – Analytical representation of spherical flame front.

97

When the skirt of the flame is close to the tube wall, the axial velocity of the flame becomes
higher than the radial velocity, which leads to the formation of an elongated flame front called
the finger flame shown in Fig.6.2. The proximity of the wall prevents the flame front from
pushing the unburnt gas radially outward, which reduces the radial expansion of the flame front.
To develop the analytical representation, the flame front is modeled as a cylindrical surface of
radius ‘r’ which ranges from 0.75-1 R, and the length range from 0.75-5.5 R, with a
hemispherical cap. The radius and length of the cylinder increase over a period of 10
milliseconds.

Figure 6.2 - Analytical representation of finger flame front.

98

The side skirt of the finger flame is eventually quenched by the cold wall, which leads to the
formation of a planar flame front as shown in Fig. 6.3. The surface area of the flame front
decreases drastically, and this leads to the decrease in the rate of pressure increase in the
combustion chamber as shown in Fig. 4.7. The planar flame area is analytically modeled as a
cuboid with a square cross section measuring 0.5 L in length, where L is the total length of the
combustion chamber.

Figure 6.3 - Analytical representation of planar frame.

99

The planar flame rapidly undergoes a change in it structure to form an inward pointing cusp
called a tulip flame. The total time taken by the flame to change its shape from an outward
pointing finger flame to an inward pointing tulip flame takes 15 milliseconds. The dominant
factor influencing the transition appears to be the hydrodynamical effects. The recirculation of
the unburnt gas causes the mass burning rate to be the highest close to the center of the cross
section when compared to the periphery. The analytical model of the tulip flame is represented as
a pyramid with a four sided base whose height is R. The height of the pyramid remains constant
throughout its propagation from L=0.5 to L=1. The analytical shape of the tulip flame is shown
in Fig. 6.4 and Fig. 6.5.

Figure 6.4 - Analytical representation of tulip flame.

100

Figure 6.5 - Analytical representation of tulip flame collision with cold wall.

The approximations used in the calculation of the flame area and volume of burnt gases are
shown below:
Stage I – Spherical Flame Front

Volume of burnt gas:

Flame surface area:
Where 0<r<0.75R.

101

Stage II – Finger Flame Front

Volume of burnt gas:

Flame surface area:
where 0.75R<r<1R and 0<h<5.4R
Stage III – Planar Flame Front
Volume of burnt gas:

Flame surface area:
where 5.4R<h<6R
Stage IV – Tulip Flame Front

Volume of burnt gas:

Flame surface area:

where 6R<h<12R

The flame surface area and burnt gas volume is computed using a MATLAB program. The
corresponding plots are shown in Fig. 6.6 and Fig. 6.7.

102

Figure 6.6 - Flame surface area vs. time.

Figure 6.7 - Burnt volume vs. time.

103

CHAPTER 7
CONCLUSIONS AND FUTURE WORK

The results obtained through experimentation and numerical analysis correspond to those obtain
by C. Clanet [14], S. T. Lee [15] and R. Starke [16]. Additional experimental results of the
velocity field, using laser doppler anemometry, shown by D. Dunn-Rankin [20] enable the
comparison to the velocity field obtained numerically. Analytical calculation performed on the
formation of the tulip flame is shown by M. Matalon [17], M. Gonzales [18] and J. W. Dold [19].
The progression of the various structures formed during the propagation of a premixed flame
front in a closed tube is verified with previous research work. Previous experimental research
work involved the propagation of a premixed flame either in semi-open or closed tubes with
circular cross-section. It is observed that an inward pointing cusp occurs in tubes with both
rectangular and square cross-section, indicating that the aspect ratio is the primary factor in the
formation of a tulip flame and not the shape of cross-section.
The flame speed profile obtained from the straight channel combustor can be utilized to obtain
port timing and rotational speed of the wave disk engine. The use of obstacles along the length of
the channel and the use of other turbulence generation mechanisms to increase the flame speed
further can be investigated in the future. The use of obstacles in long channels has been
investigated by A. R. Marsi [21], G. Cicarelli [22] and V. N. Gamezo [23].

104

APPENDICES

105

APPENDIX A

FLAME ANALYSIS MATLAB PROGRAM
The flame analysis program is written using the MATLAB programming language. The high
speed camera generates a series of images recorded at the desired frame rate. The image is
recorded in the TIFF format. The focus area of the image is equal to the camera focus area.
Step 1 –
The first step involves cropping the image from the camera output to the required combustion
chamber window. This step involves batch processing the images using IrfanView image
processing toolkit.

Step 2 –
Once the image has been cropped, the MATLAB program is used to convert the grayscale image
to a binary image based on a threshold light intensity. The value of intensity in each pixel ranges
from 0 to 255 in the grayscale image. This value is converted to either a 0 or 1 based on the user
input threshold parameter. This step is performed to eliminate noise in the image from external
sources or reflection from the glass plate.

106

Step 3 –
The program then uses progressive scans to capture the exact position of the flame front for each
position along the height of the combustion chamber. This is done by capturing the position of
the image where the value of intensity changes from 0 to 1. This ensures that the only the
position of tip of the flame front is recorded by the program. The position of the flame front
along the entire height of the channel is captured. If the program is unable to find a flame for a
specific height, the position is reset to 0.
Step 4 –
The average position of the flame front is calculated by taking the mean position of the non-zero
flame position for each image. The average position is converted from pixel units to physical
distance in centimeters using a scaling ratio. The average position of the flame front for a
sequence of images is now used to find the flame speed and flame acceleration. Additional data
calculated by the program includes calculating the position of the lobes of the tulip flame,
centerline flame speed, flame front area and difference between centerline and average flame
front position.

107

The MATLAB program source code used to post-process the cropped images is given below –
clc
clear all
filedir=uigetdir('','Select image folder');
oldfolder=cd(filedir);
imagelist=[ls('*.tif') ls('*.bmp') ls('*.jpg')];
numimages=length(imagelist);
%% Initial Parameters define the focus area, tube length, frame rate and variables to compute
left2right=0;
right2left=1;
tubefocuswidth=4;
tubefocuslength=28;
framerate=20000;
rawimage=imread(imagelist(1,:));
[d1,d2]=size(rawimage);
time=zeros(numimages,1);
binconvert=zeros(d1,d2);
flamepos=zeros(d1,numimages);
avflamepos=zeros(1,numimages);
centerflamepos=zeros(1,numimages);
flamediv=round(d1/3);
upperlobeflame=zeros(1,numimages);
lowerlobeflame=zeros(1,numimages);
userinput=75;
userflame=zeros(1,numimages);
%% Image Conversion from grayscale to binary based on threshold intensity value
for imagenumber=1:numimages
rawimage=imread(imagelist(imagenumber,:));
[d1,d2]=size(rawimage);
maxint=max(max(rawimage));
minint=min(min(rawimage));
thresh=(maxint-minint)/24;
for i=1:d1
for j=1:d2
if rawimage(i,j)>thresh
binconvert(i,j)=1;
else
binconvert(i,j)=0;
end
end
end
108

dispimage=255*binconvert;
image(dispimage)
colormap gray
axis image
M(imagenumber)=getframe;
if right2left == 1
binconvert=flipdim(binconvert,2);
end
for i=1:d1
for j=d2:-1:1
if binconvert(i,j)==1
flamepos(i,imagenumber)=(j);
break
end
end
end
end
%% Flame Position Calculation based on intensity along the height
%% Average Flame Position
for i=1:numimages
avflamepos(1,i)=(sum(flamepos(:,i))/nnz(flamepos(:,i)));
end
%% Centerline Flame Position
for i=1:numimages
centerflamepos(1,i)=flamepos((d1/2),i);
end
%% Lower and Upper Lobe Flame Position
for i=1:numimages
lowerlobeflame(1,i)=flamepos((flamediv),i);
upperlobeflame(1,i)=flamepos((2*flamediv),i);
end
%% User Input Flame Position
for i=1:numimages
userflame(1,i)=flamepos(userinput,i);
end

109

%% Flipping position matrix
avflamepos=avflamepos';
centerflamepos=centerflamepos';
lowerlobeflame=lowerlobeflame';
upperlobeflame=upperlobeflame';
userflame=userflame';
%% Computing deviation from center flame
deviation1=avflamepos-centerflamepos;
deviation2=lowerlobeflame-centerflamepos;
deviation3=upperlobeflame-centerflamepos;
deviation4=userflame-centerflamepos;
%% Converting from pixel to cm
mulfac=tubefocuslength/d2;
avflamepos=avflamepos*mulfac;
centerflamepos=centerflamepos*mulfac;
lowerlobeflame=lowerlobeflame*mulfac;
upperlobeflame=upperlobeflame*mulfac;
userflame=userflame*mulfac;
deviation1=deviation1*mulfac;
deviation2=deviation2*mulfac;
deviation3=deviation3*mulfac;
deviation4=deviation4*mulfac;
%% Flame Speed Calculations
avflamespeed = diff(avflamepos)*framerate;
centerflamespeed = diff(centerflamepos)*framerate;
lowerlobeflamespeed = diff(lowerlobeflame)*framerate;
upperlobeflamespeed = diff(upperlobeflame)*framerate;
userflamespeed = diff(userflame)*framerate;
%% Flame Surface Area Calculation
flamearea=diff(flamepos);
flamearea=flamearea*mulfac;
flamearea=flamearea.*flamearea;
flamearea=flamearea+(1/tubefocuswidth);
flamearea=sqrt(flamearea);
for i=1:(d1-1)
for j=1:numimages
if flamearea(i,j)==(2/tubefocuswidth)
110

flamearea(i,j)=0;
end
end
end
totalflamearea=sum(flamearea);
totalflamearea=totalflamearea';
clear d1 d2 binconvert dispimage filedir flamearea flamedir flamepos framerate i imagenumber j
clear left2right right2left maxint minint numimages oldfolder rawimage thresh time
tubefocuslength tubefocuswidth flamediv mulfac
%% Writing Output to MS Excel files
heading1 = {'Average Flame Speed','Center Flame Speed','Lower Lobe Flame Speed', 'Upper
Lobe Flame Speed','User Flame Speed'};
datawrite=[avflamespeed,centerflamespeed,lowerlobeflamespeed,upperlobeflamespeed,userflam
espeed];
xlswrite('flamespeed.xls',datawrite);
heading2 = {'Average Flame Position', 'Center Flame Position', 'Lower Lobe Flame Position',
'Upper Lobe Flame Position', 'User Flame Position', 'Average Deviation', 'Lower Lobe
Deviation', 'Upper Lobe Deviation', 'User Flame Deviation', 'Flame Area'};
datawrite=[avflamepos,centerflamepos,lowerlobeflame,upperlobeflame,userflame,deviation1,dev
iation2,deviation3,deviation4,totalflamearea];
xlswrite('flamepos.xls',datawrite);
xlswrite('flamespeedheading.xls', heading1);
xlswrite('flameposheading.xls',heading2);

111

APPENDIX B

SCHLIEREN DIAGNOSTIC SETUP
The Schlieren imaging technique is the method of visualization of density changes in a
transparent media. The refractivity of a gas depends on gas composition, temperature and
density. Compressible flows and combustion lead to sharp changes in the refractive index of the
media due to the change in density. Optical inhomogeneities refract light rays in proportion to
their gradients of refractive index. The resulting ray curvature is given by

The angular ray deflection is obtained by integrating the ray curvature equation described above.

These expressions provide the mathematical basis for the schlieren and shadowgraph techniques.

112

Z-Type 2-Mirror Schlieren System
The z-type Herschellian system using two oppositely tilted, on-axis spherical mirrors is the most
popular schlieren setup. The combination of a diverging illuminator beam, an opposite
converging analyzer beam and a parallel beam between the two mirrors suggest the letter “z”,
hence the name. The mirrors used in this arrangement are usually symmetrical, on-axis
paraboloid or spherical mirrors. The advantage of using mirrors over lenses is primarily found in
a larger aperture for a given cost.
A minimum distance of 2f is required between the two field mirrors, where f is the focus length
of each mirror. While longer distance between the mirrors does not lead to any errors, shorter
distances lead to a reduced sensitivity in the schlieren image. The schematic representation of the
z-type schlieren setup is shown in Fig. 6.8.

Light source

Test area

Spherical Mirror

Spherical Mirror
Camera
Figure A.1 - Z-type schlieren schematic

The apparatus requires the mirrors, light source and camera to be on the same plane, which is
achieved by either mounting all the components on a breadboard or on a flat surface. First the
light source is setup and focused to the required angle with respect to the mirror. The collimating
spherical mirror is then set on the flat surface at a predetermined angle so as to form a collimated
113

beam of light. A white sheet of paper can be used to adjust the relative position of the light
source and collimating mirror in order to obtain a perfectly parallel beam of light. A circle is
drawn on the paper with diameter equal to the diameter of the mirror. The paper is then moved
downstream of the mirror to verify that the light beam remains within the bounds of the circle
and its size does not change. Once a parallel beam of light is obtained, the focusing mirror is set
up at a distance twice the focal length of the mirror on axis downstream of the collimating
mirror. The focusing mirror is set at an angle equal to that of the collimating mirror. The
focusing mirror collects the parallel beam of light and focuses it onto the lens of the camera after
it passes through a knife edge. The mirror used in the apparatus is a spherical mirror, 6 inches in
diameter and a focal length of 60 inches. The mirror used is shown in Fig. 34.

Figure A.2 - Spherical mirror
114

The two mirrors were placed 120 inches apart and a razor blade is used as a knife edge. The
combustion chamber is placed in the marked test area, 60 inches from either mirror. The focused
light is pointed to the high speed camera after it passes through the knife edge. As the diameter
of the mirror is 6 inches and the combustion chamber is 12 inches long, only half of the
combustion chamber setup can be viewed by the schlieren set for a single experiment. To view
the density change along the entire length of the combustion chamber, multiple experiments are
conducted focusing on various parts of the tube and a composite video is constructed after
verification is made on the consistency of the experiments. The video is compiled using
commercially available VideoMach software, which stitches the boundaries of the images to
form a continuous image.

115

APPENDIX C

NUMERICAL SCHEME
The numerical simulations are performed using the commercially available computational fluid
dynamics software, Ansys Fluent. The numerical code uses the density based explicit scheme to
discretize and solve the Navier-Stokes equation for a given control volume. The density-based
scheme solves the fully-coupled Navier-Stokes equation at every time step based on the CourantFriedrich-Levy conditions for stability.
This scheme solves the governing equations of continuity, momentum, energy and species
transport simultaneously as a set, or vector, of equations. Governing equations for additional
scalars are solved sequentially.
Governing Equations
The system of governing equations for a single component fluid, written to describe the mean
flow properties, is cast in the integral Cartesian form for an arbitrary control volume V with
differential surface area dA as follows:

Here, the vectors W, F and G are defined below,

116

The vector H contains source terms such as body forces and energy sources. Here, ρ, v, E and p
are the density, velocity, total energy per until mass and pressure of the fluid respectively. τ is
the viscous stress tensor and q is the heat flux.
The total energy E is related to the total enthalpy H by
E = H – p/ρ
Where
H = h + V2/2
The Navier-Stokes equations represented above can be solved to capture the pressure waves
produced due to combustion and the effect of pressure waves on the flame geometry when
solved by satisfying the CFL conditions.

117

A fourth order Runge-Kutta scheme is used for the calculation of unsteady flows in the explicit
solver. The solution is advanced from iteration n to iteration n+1 by the following steps –
0

Q = Qn
i

-1 i-1

ΔQ = -αiΔtΓ R
n+1

Q

m

=Q

i

i

n

Where ΔQ = Q – Q and i is the stage counter for the four stage scheme. α i is the multistage
i

i

coefficient for the i stage. The residual R is computed from the intermediate solution Q and is
given by –

The time step is computed from the CFL (Courant-Friedrichs-Lewy) condition –

Where V is the cell volume, A is the face area and λ is the maximum of the local eigen values
defined by the governing equations.
To enable successful capture of the pressure waves travelling through the control volume, the
CFL number should always be set to a value below or equal to 1. This enables the numerical
scheme to capture both the travel of pressure waves and calculate the thermal and species
diffusion between the flame front interfaces.
118

REFERENCES

119

REFERENCES

[1] Ellis, O.C De C., Fuel in Science and Practice, 7:502 (1928).
[2] Salamandra, G. D, Bazhenova, T. V, Naboko, I. M, “Formation of detonation wave during
combustion of gas in combustion tube”, Symposium (International) on Combustion, 7:851
(1958).
[3] Vagani, M. et. al., “A wave disk engine concept for micro power generation”, Proceeding of
PowerMEMS, 530 (2009).
[4] Glassman, I., Yetter, R. A, “Combustion”, Chap. 4, Academic Press (2007).

[5] Mallard, E. and Le Chatelier, H. L., Ann. Mines, 4:379 (1883).

[6] Semenov, N., NACA Technical Memo, 1282 (1951).

[7] Tanford, C. and Pease, R. N., Journal of Chem. Phys., 15:861 (1947).
[8] Hirschfelder, J. O., Curtiss, C. F. and Bird, R. B., “The molecular theory of gases and
liquids”, Chap. 11, Wiley Pub, New York (1954).
[9] Linan, A. and Williams, F. A., “Fundamental Aspects of Combustion”, Oxford University
Press (1994).

[10] Mikhelson,V. A., Moscow University, Moscow (1989).
[11] Searby, G., “Instability Phenomena during Flame Propagation”, Combustion Phenomena,
CRC Press (2009).
[12] Rayliegh, G. W. S, “The explanation of certain acoustical phenomena”, Nature, 18:319-321
(1878).
120

[13] Law, C. K., “Combustion Physics”, Chap. 10 pg. 456, Cambridge University Press (2006).
[14] Clanet, C. and Searby, G., “On the tulip flame phenomenon”, Combustion and Flame,
105:225-238 (1996).

[15] Lee, S. T and Tsai, C. H., Numerical investigation of steady laminar flame propagation in a
circular tube, Combustion and Flame, 99:484-490 (1994).
[16] Starke, R. and Roth, P., “An experimental investigation of flame behavior during cylindrical
vessel explosions”, Combustion and Flame, 66:249-259 (1986).
[17] Matalon, M., McGreevy, J. L., “The initial development of a tulip flame”, Symposium on
Combustion, 25:1407-1413 (1994).
[18] Gonzales, M., Borghi, R., Saouab, A., “Interaction of a flame front with its self-generated
flow in an enclosure: The tulip flame phenomenon”, Combustion and Flame, 88:201-220 (1992).
[19] Dold, J. W., Joulin, G., “An evolution equation modeling inversion of tulip flame”,
Combustion and Flame, 100:450-456 (1995).
[20] Dunn-Rankin, D., Sawyer, R. F., “Tulip flames: Changes in shapes of premixed flames
propagating in closed tubes”, Experiments in Fluids, 24:130-140 (1998).
[21] Marsi, A. R., Ibrahim, S. S, Nehzat, N., “Experimental study of premixed flame propagation
over various solid obstacles”, “Experimental Thermal and Fluid Science”, 21:109-116(2000).
[22] Ciccarelli, G., Fowler, C., G., Bardon, M., “Effect of obstacle size and spacing on the initial
stage of flame acceleration in a rough tube”, “Shock Waves”, 14:161-166 (2005).
[23] Gamezo, V. N., Oran, E. S., “Flame acceleration and DDT in channels with obstacles”,
“Combustion and Flame”, 155:302-315 (2008).

121