MODELING AND CONTROL OF PRE-CHAMBER INITIATED TURBULENT JET
IGNITION COMBUSTION SYSTEMS
By
Ruitao Song

A DISSERTATION
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
Mechanical Engineering – Doctor of Philosophy
2018

ABSTRACT
MODELING AND CONTROL OF PRE-CHAMBER INITIATED TURBULENT JET
IGNITION COMBUSTION SYSTEMS
By
Ruitao Song
Turbulent jet ignition (TJI) combustion is a promising concept for achieving high thermal efficiency and low NOx (nitrogen oxides) emissions. A control-oriented TJI combustion model with
satisfactory accuracy and low computational effort is usually a necessity for optimizing the TJI
combustion system and developing the associated model-based TJI control strategies. A controloriented TJI combustion model was first developed for a rapid compression machine (RCM) configured for TJI combustion. A one-zone gas exchange model is developed to simulate the gas
exchange process in both pre- and main-combustion chambers. The combustion process is modeled by a two-zone combustion model, where the ratio of the burned and unburned gases flowing
between the two combustion chambers is variable. To simulate the influence of the turbulent jets to
the rate of combustion in the main-combustion chamber, a new parameter-varying Wiebe function
is proposed and used for mass fraction burned (MFB) calculation in the main-combustion chamber.
The developed model is calibrated using the Least-Squares fitting and optimization procedure.
The RCM model was then extended to a TJI engine model. The combustion process is modeled
by a similar two-zone combustion model based on the newly proposed parameter-varying Wiebe
function. The gas exchange process is simulated by one-zone model considering piston movement
and intake and exhaust processes. Since the engine uses liquid fuel, a pre-chamber air-fuel mixing
and vaporization model is developed. And correspondingly, the pre-chamber uses a chemical kinetics based model for combustion rate calculation. The model was validated using the experimental
data from a single cylinder TJI engine under different operational conditions, and the simulation
results show a good agreement with the experimental data.
For control design, a nonlinear state-space engine model with cycle-to-cycle dynamics is developed based on the previous crank-angle-resolved (CAR) TJI engine model. The state-space

model successfully linked the combustion processes in the two chambers using the parametervarying Wiebe function. The validated CAR model is used to calibrate and validate the state-space
engine model. The simulation results of the two engine models show a good agreement with
each other. Thereafter, a linear-quadratic tracking controller is developed for combustion phasing
control. Simulation results are presented and a baseline controller has been implemented on the
research engine.
Combustion phasing control is very important for internal combustion engines to achieve high
thermal efficiency with low engine-out emissions. Traditional open-loop map-based control becomes less favorable in terms of calibration effort, robustness to engine aging, and especially control accuracy for TJI engines due to the increased number of control variables over conventional
spark-ignition engines. In this research, a model-based feedforward controller is developed for
the TJI engine, and a feedback controller is also designed based on the linear quadratic tracking
control with output covariance constraint. Since the TJI main-chamber combustion is influenced
by the pre-chamber one, the proposed controller optimizes the control variables in both combustion chambers. The proposed feedforward and feedback controllers show significant performance
improvement over a group of baseline controllers through a series of dynamometer engine tests.

Copyright by
RUITAO SONG
2018

ACKNOWLEDGEMENTS

This dissertation would not have been possible without the help of so many people. I would like to
express my sincere appreciation to those who provided me guidance, support, encouragement, and
also brought a lot of fun to my research and life.
I would like to express my deepest gratitude to my advisor, Dr. Guoming Zhu. It is my fortune
and great privilege to do my research work under his continued support and excellent guidance.
His comprehensive knowledge of control theory and automotive application inspired me to keep
pursuing and enriching my knowledge and experience. His systematic guidance encouraged me to
improve the way of thinking and exploring. Besides, he also gave me generous encouragement for
my future career and life, which is a lifetime of wealth!
I would like to express my appreciation to my committee members: Dr. Harold Schock, Dr.
Hassan Khalil, and Dr. Ranjan Mukherjee, for their discussions and insight comments during my
dissertation work. Their wonderful graduate courses also provided me valuable knowledge and
firm foundation for my research work, and brought me to the world of engine and control.
I would also like to say thanks to the graduate students and research assistants in our lab, who
helped me set up the engine bench and run the experiments, even day and night: Kevin Moran,
Tom Stuecken, Ravi Teja Vedula, Jie Yang, Tao Zeng, Yifan Men, Yingxu Wang, Chengsheng
Miao, Huan Li, Christine Li, Tianyi He, Shen Qu.
Lastly, I am overly grateful for my family’s continued support and warm hearted solicitude that
makes me full of courage to strive for the future.

v

TABLE OF CONTENTS

LIST OF TABLES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii
LIST OF FIGURES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

ix

CHAPTER 1 INTRODUCTION . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 Research Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2.1 Modeling of TJI combustion in a rapid compression machine . . . .
1.2.2 Combustion modeling for the Dual-Mode TJI engine . . . . . . . .
1.2.3 Optimal combustion phasing control for the Dual-Mode TJI engine .
1.3 Dissertation Contributions . . . . . . . . . . . . . . . . . . . . . . . . . .
CHAPTER 2
2.1
2.2
2.3

2.4
2.5
2.6

MODELING OF TURBULENT JET
RAPID COMPRESSION MACHINE
Introduction . . . . . . . . . . . . . . . . .
System Description . . . . . . . . . . . . .
Turbulent Jet Ignition Combustion Model .
2.3.1 Gas exchange model . . . . . . . .
2.3.2 Two-zone combustion model . . . .
2.3.3 Mass fraction burned model . . . .
Model Calibration . . . . . . . . . . . . . .
Model Validation and Simulation Results . .
Conclusion . . . . . . . . . . . . . . . . .

.
.
.
.
.
.
.

.
.
.
.
.
.
.

.
.
.
.
.
.
.

1
1
3
3
3
4
4

IGNITION COMBUSTION IN A
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .

.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.

6
6
8
10
10
12
16
18
21
28

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

30
30
31
31
32
33
36
38
42
43
45
54

CHAPTER 3
3.1
3.2

3.3
3.4
3.5

3.6
3.7

A CONTROL-ORIENTED COMBUSTION MODEL FOR
LENT JET IGNITION ENGINE USING LIQUID FUEL . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
TJI Combustion Model Architecture . . . . . . . . . . . . . . . .
3.2.1 Target engine description . . . . . . . . . . . . . . . . . .
3.2.2 Engine events and combustion model architecture . . . . .
Gas Exchange Model . . . . . . . . . . . . . . . . . . . . . . . .
Fuel Mixing and Vaporization Model . . . . . . . . . . . . . . . .
Two-zone Combustion Model . . . . . . . . . . . . . . . . . . . .
3.5.1 Pre-chamber combustion rate calculation . . . . . . . . . .
3.5.2 Main-chamber combustion rate calculation . . . . . . . .
Model Calibration and Validation . . . . . . . . . . . . . . . . . .
Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

CHAPTER 4
4.1
4.2

OPTIMAL COMBUSTION PHASING MODELING AND
OF A TURBULENT JET IGNITION ENGINE . . . . . . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Target Engine Description and Specifications . . . . . . . . . . .
vi

.
.
.
.
.
.
.

A TURBU. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .
. . . . . . .

CONTROL
. . . . . . . . . 55
. . . . . . . . . 55
. . . . . . . . . 57

4.3

4.4

4.5

4.6

Nonlinear State-Space TJI Engine Model . .
4.3.1 Main-chamber state-space model . .
4.3.2 Pre-chamber state-space model . . .
4.3.3 CA50 calculation . . . . . . . . . .
Controller Development . . . . . . . . . . .
4.4.1 Feedforward control . . . . . . . .
4.4.2 Output covariance constraint control
4.4.3 Baseline controllers . . . . . . . . .
Experimental Validation . . . . . . . . . . .
4.5.1 Model validation . . . . . . . . . .
4.5.2 Controller validation . . . . . . . .
Conclusions . . . . . . . . . . . . . . . . .

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

58
58
62
65
66
67
69
73
74
74
76
80

CHAPTER 5 CONCLUSIONS AND FUTURE WORK . . . . . . . . . . . . . . . . . . . 82
5.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
5.2 Recommendations for Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . 83
APPENDICES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
APPENDIX A NONLINEAR STATE-SPACE MODEL . . . . . . . . . . . . . . . . 86
APPENDIX B ENGINE MODEL PARAMETERS . . . . . . . . . . . . . . . . . . . 88
BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89

vii

LIST OF TABLES

Table 2.1: RCM specifications. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

8

Table 2.2: Experimental set-up. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
Table 2.3: Calibration results for cases 1 and 2. . . . . . . . . . . . . . . . . . . . . . . . . 21
Table 2.4: Simulation parameters. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
Table 2.5: Calibration results for cases 4-6. . . . . . . . . . . . . . . . . . . . . . . . . . . 28
Table 3.1: TJI Engine Specifications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
Table 3.2: Calibration Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
Table 4.1: PI controller parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
Table B.1: Parameters for Combustion Phase Model . . . . . . . . . . . . . . . . . . . . . 88

viii

LIST OF FIGURES

Figure 2.1: Rapid compression machine . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6

Figure 2.2: Working process of the RCM . . . . . . . . . . . . . . . . . . . . . . . . . . .

9

Figure 2.3: Two-zone combustion model. . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
Figure 2.4: The value of Îąb in the two cases. . . . . . . . . . . . . . . . . . . . . . . . . . 15
Figure 2.5: Inverse thermodynamic calculation results. . . . . . . . . . . . . . . . . . . . . 20
Figure 2.6: MFB obtained from inverse thermodynamic calculation and Least-Squares
fitting result. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
Figure 2.7: Experimental and calculated pressure traces for experimental case 1. . . . . . . 23
Figure 2.8: Experimental and calculated pressure traces for experimental case 2. . . . . . . 23
Figure 2.9: Experimental and calculated pressure traces for experimental case 3. . . . . . . 24
Figure 2.10: Experimental and calculated net heat release rates. . . . . . . . . . . . . . . . . 25
Figure 2.11: Chemical energy release rate calculation. . . . . . . . . . . . . . . . . . . . . . 25
Figure 2.12: The calculated value of Îąb during combustion. . . . . . . . . . . . . . . . . . . 26
Figure 2.13: Experimental and calculated pressure traces for experimental case 4. . . . . . . 27
Figure 2.14: Experimental and calculated pressure traces for experimental case 5. . . . . . . 27
Figure 2.15: Experimental and calculated pressure traces for experimental case 6. . . . . . . 28
Figure 2.16: Experimental and calculated Burn1050 in the main-combustion chamber. . . . . 29
Figure 2.17: Experimental and calculated Burn5090 in the main-combustion chamber. . . . . 29
Figure 3.1: TJI engine architecture . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
Figure 3.2: Key events over an engine cycle for TJI engine modeling. IVO: intake valve
opening; IVC: intake valve closing; EVO: exhaust valve opening; EVC: exhaust valve closing; SPK: spark timing; and INJ: pre-chamber fuel and air
injection. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
ix

Figure 3.3: Pre-chamber fuel vaporization model. . . . . . . . . . . . . . . . . . . . . . . 37
Figure 3.4: Two-zone combustion model. Îąb: area fraction for burned gas . . . . . . . . . 38
Figure 3.5: The value of Îąb in the two cases. . . . . . . . . . . . . . . . . . . . . . . . . . 41
Figure 3.6: TJI engine at MSU. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
Figure 3.7: Experiment matrix. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
Figure 3.8: Pressure traces obtained at 1500 rpm. . . . . . . . . . . . . . . . . . . . . . . . 49
Figure 3.9: Heat release rates in two chambers. . . . . . . . . . . . . . . . . . . . . . . . . 50
Figure 3.10: The value of Îąb and other related parameters during simulation. . . . . . . . . 51
Figure 3.11: Pre-chamber pressure traces and main-chamber heat release rates for different pre-chamber Îť . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
Figure 3.12: Main-chamber burn duration with the same main-chamber properties and different pre-chamber Îť . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
Figure 3.13: IMEP comparison between simulation and experimental results. . . . . . . . . . 53
Figure 3.14: Spk-CA50 comparison between simulation and experimental results. . . . . . . 53
Figure 3.15: Pre-chamber first peak pressures and their crank positions. . . . . . . . . . . . . 53
Figure 4.1: TJI engine architecture . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
Figure 4.2: State-space TJI engine model architecture. . . . . . . . . . . . . . . . . . . . . 59
Figure 4.3: Control system architecture. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
Figure 4.4: Experimental matrix. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
Figure 4.5: State-space engine model validation. . . . . . . . . . . . . . . . . . . . . . . . 75
Figure 4.6: Comparison between experimental and calculated CA50 under different operational conditions. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
Figure 4.7: Comparison between experimental and calculated CA50 when pre-chamber
fuel is changed. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
Figure 4.8: Test bench structure. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78
x

Figure 4.9: Experimental results with controller on and off. . . . . . . . . . . . . . . . . . 78
Figure 4.10: Experimental results of OCC and PI controllers. . . . . . . . . . . . . . . . . . 79
Figure 4.11: Experimental results of OCC controllers with and without using pre-chamber
fuel. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
Figure 4.12: Experimental results of OCC controllers under different conditions. . . . . . . . 81
Figure 4.13: Standard deviation of CA50 with different controllers. . . . . . . . . . . . . . . 81

xi

CHAPTER 1
INTRODUCTION

1.1

Motivation

A major focus in light-duty vehicle development is to improve the thermal efficiency of directinjection spark ignition (DISI) engines and reduce exhaust emissions using various strategies such
as high compression ratio, charge dilution, tumble enhancement, high ignition energy, and late
intake valve closure timing (Miller or Atkinson cycles). Most passenger cars utilize less than 10%
of their maximum engine power during daily commute, according to [1]. This fact highlights the
importance of improving engine part-load efficiency. According to a recent benchmarking study
conducted by the US EPA, the current production engines carry a part-load brake thermal efficiency
ranging from 30 to 35% [2]. An exception to this range are the Mazda’s 2.0L engine that showed a
peak brake efficiency of 37% under part-load conditions and Toyota’s 1.3L Atkinson cycle engine
of compression ratio 13.5 with a brake efficiency of 38%. Reese [3] made a “propulsion system
efficiency” analysis based on 2015 US EPA certification data to estimate the net improvements
required in the engine, transmission, and driveline efficiencies to meet US 2025 GHG regulations.
According to this analysis, a further 30% reduction in fuel consumption is required for gasoline
engines if all this improvement is to originate solely from engine development [3, 4]. If a thermal
efficiency of 35% is considered as the current industry standard, this analysis would indicate a
brake thermal efficiency of 50% (=35%/0.7) to meet the 2025 regulations.
New combustion technologies were proposed and studied during the past few decades. Some
of these technologies have already been brought to production, such as engine downsizing [5]
and Atkinson cycle based engines for many hybrid electric vehicles. Lean burn technologies, like
the homogeneous charge compression ignition (HCCI) combustion [6, 7], also attract significant
research interest due to their excellent fuel economy and extremely low emissions. The turbulent
jet ignition (TJI) combustion system is another promising combustion technology that has the

1

potential to be widely used in the next generation IC engines.
HCCI combustion has been widely investigated in past decades, and demonstrated the potential of providing higher fuel thermal efficiency and lower emissions than those of the conventional
spark ignition (SI) combustion. The un-throttled HCCI combustion significantly reduces the pumping loss, and the lean HCCI combustion results in relatively low in-cylinder flame temperature and
has the advantage of significantly reducing the level of NOx emissions. However, since the HCCI
combustion is not suitable for all engine operational conditions, especially under high load and
speed conditions, mode transition between SI for high load and speed operations and HCCI combustion for low or mediate load and speed conditions is required. The combustion mode transition
control, along with the combustion phase control, are two of the key challenges for the HCCI combustion technology [7, 8]. Whereas, the TJI combustion is able to cover the entire speed and load
range of a typical SI engine and the start of combustion can be easily controlled by adjusting the
pre-chamber spark timing in the TJI system. Note that as the engine load increases, the achievable
lean limit decreases.
Turbulent jet ignition is a pre-chamber initiated two-chamber combustion system [9]. The TJI
system mainly consists of three key components: a small pre-chamber, a multi-orifice nozzle, and
a large main chamber. A spark-ignited reactive mixture in the pre-chamber flows through the
nozzle orifices and results in multiple, chemically active, turbulent jets that emerge into the main
chamber. TJI combustion system is able to greatly reduce the NOx emissions while maintains a
comparatively low HC and CO emissions, especially when the relative air-fuel ratio (AFR), Îť , is
greater than 1.4. According to the early research[10], stable combustion can be achieved for the
TJI system when Îť is up to 1.8 with extremely low NOx emissions. Furthermore, the lean burn
TJI combustion has been recorded to have an 18% improvement in fuel consumption, comparing
to the conventional stoichiometric SI combustion [11].
The TJI combustion system introduces additional freedoms for combustion control. For instance, pre- and main-chamber AFRs, pre-chamber fuel injection and ignition timing. This makes
it difficult to utilize the traditional single variable feedback and look-up table based control, and it

2

often requires real-time optimization, especially under the transient operations. Therefore, modelbased control is required to handle multiple-input and multiple-output TJI combustion system.

1.2
1.2.1

Research Overview
Modeling of TJI combustion in a rapid compression machine

A control-oriented TJI combustion model with satisfactory accuracy and low computational effort is usually a necessity for optimizing the TJI combustion system and developing the associated model-based TJI control strategies. This research presents a control-oriented TJI combustion
model developed for a rapid compression machine (RCM) configured for TJI combustion. An
one-zone gas exchange model is developed to simulate the gas exchange process in both preand main-combustion chambers. The combustion process is modeled by a two-zone combustion
model, where the ratio of the burned and unburned gases flowing between the two combustion
chambers is variable. To simulate the influence of the turbulent jets to the rate of combustion in
the main-combustion chamber, a new parameter-varying Wiebe function is proposed and used for
mass fraction burned (MFB) calculation in the main-combustion chamber. The developed model
is calibrated using the Least-Squares fitting and optimization procedures. Experimental data sets
with different AFR in both combustion chambers and different pre-combustion chamber orifice
areas are used to calibrate and validate the model. The simulation results show a good agreement
with the experimental data for all the experimental data sets. This indicates that the developed
combustion model is accurate for developing and validating TJI combustion control strategies.

1.2.2

Combustion modeling for the Dual-Mode TJI engine

For engines equipped with the turbulent jet ignition (TJI) system, the interaction between the preand main-combustion chambers should be considered in the control-oriented model for developing
control strategies that optimize the overall thermal efficiency in real-time. Therefore, a two-zone
combustion model based on the newly proposed parameter-varying Wiebe function is proposed.
Since the engine uses the liquid fuel, a pre-chamber air-fuel mixing and vaporization model is also
3

developed. The model was validated using the experimental data from a single cylinder TJI engine
under different operational conditions, and the simulation results show a good agreement with the
experimental data.

1.2.3

Optimal combustion phasing control for the Dual-Mode TJI engine

Combustion phasing control is very important for internal combustion engines to achieve high thermal efficiency with low engine-out emissions. Traditional open-loop map-based control becomes
less favorable in terms of calibration effort, robustness to engine aging, and especially control accuracy for turbulent jet ignition (TJI) engines due to the increased number of control variables over
conventional spark-ignition engines. In this research, a model-based feedforward controller is developed for a TJI engine. A feedback controller is designed based on the linear quadratic tracking
control with output covariance constraint. Since the TJI main-chamber combustion is influenced
by the pre-chamber one, the proposed controller optimizes the control variables in both combustion chambers. The proposed feedforward and feedback controllers show a significant performance
improvement over a group of baseline controllers through a series of dynamometer engine tests.

1.3

Dissertation Contributions

The dissertation has the following major contributions:
• A newly proposed parameter-varying Wiebe combustion model is used to link the combustion processes in both pre- and main-combustion chambers. The developed model can be
calibrated using a simple and systematic calibration procedure based on the experimental
data.
• The air-fuel mixture in the two combustion chambers are modeled using two zones (burned
and unburned) to simulate the complicated gas exchange process between the two chambers.
• Different from many other TJI engines, liquid gasoline is used as the fuel for both the preand main-chambers. Since liquid gasoline is hard to be fully vaporized inside the small
4

pre-chamber before combustion starts, an air-fuel mixing and fuel vaporization model is
proposed for the pre-chamber in this resreach.
• The development of a state-space TJI engine model for model-based control and design of a
combustion phasing controller, consisting of feedforward and feedback control, for the TJI
engine based on the unique TJI combustion characteristics.

5

CHAPTER 2
MODELING OF TURBULENT JET IGNITION COMBUSTION IN A RAPID
COMPRESSION MACHINE

2.1

Introduction

The TJI combustion system was proposed almost a century ago. In 1918, Harry R. Ricardo first
developed and patented the engine using a TJI system[12]. In 1970s, more research efforts were
devoted to the development of new TJI systems. Honda developed the compound vortex controlled
combustion (CVCC) system [13] that is considered the most significant development in Ottocycle
engines with the TJI system. It was able to meet the 1975 emission standards without a catalytic
converter.
A typical TJI system consists of a main-combustion chamber and a small pre-combustion
chamber. Its volume is a few percent of that of the main-combustion chamber. The two combustion chambers are connected through a few small orifices. The air-fuel mixture is lean in the
main-combustion chamber and relatively rich (or close to stoichiometric) in the pre-combustion
chamber to make spark ignition (SI) easy. Consequently, the TJI system usually needs two fuel delivering systems for the two combustion chambers. The combustion process is initiated by a spark
inside the pre-combustion chamber. Then, the turbulent jets of the reacting products from the pre-

Figure 2.1: Rapid compression machine

6

combustion chamber flow into the main-combustion chamber and ignite the air-fuel mixture in the
main-combustion chamber.
The TJI combustion possesses many advantages over other combustion technologies. One of
the approaches to reduce the NOx (nitrogen oxides) emissions is to operate the engine at very lean
conditions with its relative air-to-fuel ratio (AFR) greater than 1 since the resulting relatively low
temperature combustion leads to a significant reduction of NOx formation. Note that significant
NOx emission reduction can only be achieved at the extremely lean condition for conventional SI
engines[14]. Extremely lean operation of conventional SI engines will lead to poor combustion
stability with high occurrence of misfire due to the narrow fuel flammability limits. Lean mixture
also has very slow laminar flame speed that often leads to incomplete combustion. As a result,
lean operation in conventional SI engines significantly increases the HC (hydrocarbon) and CO
(carbon monoxide) emissions. However, in the TJI combustion system, the mixture in the maincombustion chamber is ignited by the hot turbulent jet that contains much higher energy than what
a spark plug can provide[15]. As a result, lean air-fuel mixture can be ignited and burned at a very
fast rate with high combustion stability. Therefore, TJI combustion system is able to greatly reduce
the NOx emissions while maintain a comparatively low HC and CO emissions, especially when the
relative AFR, Îť , is greater than 1.4. According to the previous research, stable combustion can be
achieved for the TJI system when Îť is up to 1.8[10], approaching elimination of NOx emissions.
The spark timing, rate of combustion and other combustion parameters in the TJI system need
to be optimized by the control strategies to achieve the best fuel efficiency with reduced emissions.
To develop and validate the TJI combustion control strategies, a control-oriented TJI combustion
model is also required. Toulson [16] modeled a TJI engine using the computational fluid dynamics
(CFD) method. Ghorbani [17] modeled a transient turbulent jet by the probability density function (PDF) method. These investigations provide insight to better understand the TJI combustion.
However, these models are too detailed to be used for model-based control. Model-based combustion control requires a simple combustion model capable of capturing the TJI system dynamics
with good accuracy, low computational and calibration efforts [18].

7

In this chapter, a control-oriented TJI combustion model is developed for a rapid compression machine (RCM) equipped with a TJI system. The gas exchange process in the combustion
chambers before combustion is simulated by a one-zone gas exchange model. It is based on the
assumption that the air and fuel are uniformly mixed in both combustion chambers. After ignition,
both combustion chambers are divided into burned and unburned zones. The ratio between the
burned and unburned gases flowing through the orifices connecting the two combustion chambers
are adjusted due to the tiny pre-combustion chamber to improve the model accuracy. To link the
two combustion processes in both combustion chambers, a new parameter-varying Wiebe function
is proposed and used for the main-combustion chamber mass fraction burned (MFB) calculation.
The newly proposed Wiebe function allows the combustion rate in the main-combustion chamber
to vary based on the characteristics of turbulent jets from the pre-combustion chamber, which is
one of the key features of the TJI combustion.

2.2

System Description

Fig. 2.1 shows the basic architecture of the RCM equipped with a TJI system modeled in this
chapter. There is an auxiliary fuel system injecting the methane into the pre-combustion chamber.
The two combustion chambers are connected by a small orifice. The detailed parameters of the
system are listed in Table 2.1.
Table 2.1: RCM specifications.
Parameter

Value

bore
stroke
compression ratio
pre-combustion chamber volume
pre-combustion chamber orifice diameter

50.8 mm
20.2 cm
8.5:1
2.3 cm3
1.5-3.0 mm

Here, methane was used as the fuel for all the experiments. Two Kistler piezoelectric pressure
sensors were installed into the two combustion chambers for pressure measurements. During the
experiment, the combustion chamber wall was heated to 80◦ C. The combustion chambers were

8

Pre−chamber
Main−chamber
Fuel Injection Control
Dwell Control

30

40
30

20

20

10

10

0

0

5

10

15

20
25
Time (ms)

30

35

40

Voltage (V)

Pressure (bar)

40

0

Figure 2.2: Working process of the RCM

firstly evacuated by a vacuum pump and then filled with air-fuel mixture with a known AFR. Then,
the piston rapidly compressed the mixture in both combustion chambers. At the same time, a
charge of fuel was injected into the pre-combustion chamber; see Fig. 2.2. At the end of compression, the piston kept still; and therefore, the volume in the main-combustion chamber remained
constant. At the falling edge of the dwell control signal, the spark is initiated through the spark
plug inside the pre-combustion chamber and then the reacting products from the pre-combustion
chamber were injected into the main-combustion chamber and ignited the air-fuel mixture. Fig. 2.2
shows the control signals and the typical pressure traces measured during the experiment.
In order to have a good mixing, the fuel was injected during the compression. The fuel flow
mixes with the gas flowing through the orifice from the main-combustion chamber during the
injection, leading to a better mixing. Moreover, to allow enough time for the mixing process after
fuel injection, the injection timing was set at the beginning of the compression. In this way, we
are able to make sure that the air-fuel mixture in the pre-combustion chamber is close to uniformly
mixed.

9

2.3
2.3.1

Turbulent Jet Ignition Combustion Model
Gas exchange model

During compression, methane is injected into the pre-combustion chamber. The mass flow rate of
the injected methane can be calculated by the one-dimensional compressible flow equation [1].
Pin j
ṁin j = Cd1 Av1 p
ψ
RTin j
where




Ppre
, Pin j > Ppre
Pin j

v
u


 (Îş+1)
 Îş
u 


t
2
2
(κ−1)
(κ−1)


x<
 κ
Îş +1
Îş +1
ψ (x) =
s




 Îş


(κ−1)
1
2Îş
2
(κ−1)


1−x κ
x≥
x κ
κ −1
Îş +1

(2.1)

(2.2)

Note that the coefficient Cd1 is experimentally determined; Av1 is the orifice area of the fuel
injector; Îş is the ratio of specific heats; R is the gas constant; Pin j and Tin j are the upstream pressure
and temperature, respectively; and Ppre is the pressure in the pre-combustion chamber.
The gas exchange process between the two combustion chambers is modeled similarly. However, the pressure in the pre-combustion chamber can be either greater or less than that in the
main-combustion chamber. The mass flow rate between the two combustion chambers is calculated by the following equation.

ṁtur =



Ppre
Pmain


Ppre ≥ Pmain

 Cd2 Av2 pRT ψ P
pre
pre



Ppre
P

main

ψ
Ppre < Pmain
−Cd2 Av2 √
Pmain
RTmain

(2.3)

where Cd2 and Av2 are the discharge coefficient and the area of the orifice connecting the two combustion chambers. The subscripts pre and main denote the pre-combustion and main-combustion
chamber properties, respectively.

10

Before ignition, the pre-combustion chamber is considered as a control volume with mass and
energy exchange. The mass and energy conservation equations are used to describe such a control
volume.

dm pre
= ṁin j − ṁtur
dt
dU pre
= Ḣin j − Ḣtur − Q̇ht
dt

(2.4)

where m pre and U pre are the mass and internal energy of the gas in the pre-combustion chamber,
respectively; Q̇ht is the heat transfer rate through the chamber wall; and H is the enthalpy flow.
The subscript in j and tur represent the properties of the gas from the fuel injector and through the
orifice connecting the two combustion chambers, respectively.
Assuming that the gas can be considered as an ideal gas, the two equations can be coupled by
the ideal gas law below.

Ppre ¡Vpre = m pre ¡ R ¡ Tpre

(2.5)

where Vpre is the pre-combustion chamber volume. Substituting Eqn. (2.5) into Eqn. (2.4), the
following two equations are obtained to calculate the gas pressure and temperature.


dPpre
R 
=
c p ṁin j Tin j − c p ṁtur Ttur − Q̇ht
dt
Vpre cv
dTpre
Tpre R 
=
c p ṁin j Tin j − c p ṁtur Ttur
dt
PpreVpre cv



−cv ṁin j − ṁtur Tpre − Q̇ht

(2.6)

where c p and cv are the specific heat at constant pressure and constant volume, respectively, and

Ttur =



Tpre

ṁtur > 0

(2.7)


T
main ṁtur ≤ 0
The pre-combustion chamber volume is only around 2-4% of the main-combustion chamber
clearance volume. Therefore, the gas flowing between the two combustion chambers can be
11

Figure 2.3: Two-zone combustion model.

neglected for the main-combustion chamber model. The pressure and temperature in the maincombustion chamber can be solved using the energy and mass conservation equations. These
equations can be found in many other articles about engine modeling[6, 18] and thus will not be
shown here.

2.3.2

Two-zone combustion model

After ignition, the pre-combustion chamber is divided into two zones to improve the model accuracy. Both the burned and unburned zones can be regarded as control volumes. Besides the
mass, enthalpy, and work exchange between the two control volumes, there are also mass and enthalpy exchange through the orifice to the main-combustion chamber; see Fig. 2.3. The burned
and unburned gases in the pre-combustion (or main-combustion) chamber are assumed to enter
the burned zone and unburned zone in the main-combustion (or pre-combustion) chamber, respectively. Before the ignition of the main-combustion chamber, all the gas from the main-combustion
chamber is considered as unburned gas. The energy balance equation of the burned zone is shown
in Eqn. (2.8). To make the equations concise, the variables in the following equations in this
subsection are for the pre-combustion chamber, if not specified.

dV
d(mb Tb )
+ P b + xb Q̇ht =
dt
dt
mu dxb
Q̇ch +
c p Tu − ṁtur−b c p Tb
1 − xb dt
cv

12

(2.8)

The energy balance equation of the unburned zone is represented by

d(mu Tu )
dVu
+P
+ (1 − xb ) Q̇ht =
dt
dt
mu dxb
−
c p Tu − ṁtur−u c p Tu
1 − xb dt

cv

(2.9)

The masses of both burned and unburned zones are obtained based on the following mass
conservation law.

mu dxb
dmb
= −ṁtur−b +
dt
1 − xb dt
dmu
mu dxb
= −ṁtur−u −
dt
1 − xb dt

(2.10)

The subscripts b and u represent burned zone and unburned zone. Qht is the heat transfer to
the chamber wall. Qch is the chemical energy released by combustion. xb is the MFB. mtur−b and
mtur−u represent the burned gas and unburned gas flowing through the orifice. Correspondingly,
the area of the orifice is also divided into two parts. One for the burned gas and the other for the
unburned gas. Fig. 2.3 shows the basic idea of the two-zone combustion model. When the pressure
in the pre-combustion chamber is greater than that in the main-combustion chamber, the two mass
flow rates are calculated by


Ppre
Pmain
ψ
ṁtur−b = αbCd2 Av2 √
RTb
Ppre


Ppre
Pmain
ṁtur−u = (1 − αb )Cd2 Av2 √
ψ
Ppre
RTu

(2.11)

Similar result can be obtained when the pressure in the main-combustion chamber is greater than
that in the pre-combustion chamber.
The coefficient Îąb in Eqn. (2.11) is chosen as a function of the volume fraction of the burned
gas vb . Assuming that the burned and unburned gases were always well mixed, Îąb would be
always equal to vb . However, in reality, this is not the case. Îąb is combustion chamber structure
dependent. For our TJI system, the spark plug is located at the top of the pre-combustion chamber;
see Fig. 1. In this case, the combustion is initiated at the top of the pre-combustion chamber.
13

Since the orifice is at the bottom, it is hard for the burned gas to escape from the pre-combustion
chamber at the early stage of the combustion. As a result, the fraction of the burned gas flowing
through the orifice to the main-combustion chamber is much smaller than the burned gas fraction
inside the pre-combustion chamber. This is why Îąb is smaller than vb in the pre-combustion
chamber when Ppre > Pmain . Note that Îąb will be determined using experimental data. When
Ppre < Pmain , the gas in the main-combustion chamber flows through the orifice and Îąb will be
determined by the burned gas fraction in the main-combustion chamber. Since the combustion
in the main-combustion chamber is initiated by the turbulent jet (close to orifice), the orifice is
surrounded by the gas with high concentration of burned gas. Therefore, Îąb is larger than vb in
the main-combustion chamber. And again, the actual value will be determined by the experimental
data. This is the main reason why two-zone combustion model is used. The value of Îąb can be
expressed by Eqn. (2.12).

Îąb =





 f1 vb−pre , c pre

Ppre > Pmain

(2.12)


 f (v
2 b−main , cmain ) Ppre < Pmain
To simplify the calibration process, the two functions, f1 and f2 , are approximated by seconddegree Bézier curves [19]. Besides the control points (0,0) and (1,1), (c pre ,1 − c pre ) was added for
f1 and (cmain ,1 − cmain ) for f2 as the third control points; see Fig. 2.4. The parameters c pre and
cmain are experimentally determined. The BÊzier curve guarantees ιb ∈ [0, 1] as long as c pre ∈
[0, 1] and cmain ∈ [0, 1]. By changing c pre and cmain , the ratio of the burned and unburned gases
flowing through the orifice can be adjusted to better match the actual physical process and thus to
improve the model accuracy.
Applying the principle of mass conservation, the instant fuel mass in the pre-combustion chamber can be obtained by


m pre− f uel dxb
dm pre− f uel
1
=−
− ṁtur−u
dt
1 − xb dt
Îť (A/F)s + 1

14

(2.13)

f1 (Ppre>Pmain)

f2 (Ppre<Pmain)
1

0.8

0.8

0.6

0.6

0.4

0.4

(cmain,1−cmain)

Îąb

1

0.2
0

0.2
(cpre,1−cpre)
0

0.5

0

1

0

0.5
vb−main

vb−pre

1

Figure 2.4: The value of Îąb in the two cases.

where Îť is the relative AFR, and (A/F)s is the stoichiometric AFR. Note that only the fuel from
the unburned zone is considered.
From Eqn. (2.11) and Eqn. (2.13), it can be observed that the total amount of fuel burned inside
the pre-combustion chamber is highly influenced by Îąb .
The rate of chemical energy release (CER) is obtained by the following relationship.

Q̇ch = η pre QLHV

m pre− f uel dxb
1 − xb dt

(2.14)

where the combustion efficiency Ρ pre is experimentally determined, and QLHV is the lower heating
value of the fuel.
The rate of heat transfer to the combustion chamber wall can be modeled by the following
equation[20].

Q̇ht = A pre hc Tpre − Tw




(2.15)

where A pre is the pre-combustion chamber surface area; Tw is the mean wall temperature; and hc
is the heat-transfer coefficient calibrated by the experiment.
After the ignition in the pre-combustion chamber, the combustion in the main-combustion
chamber will not be initiated until the generation of the turbulent jet from the pre-combustion
chamber. Before the ignition of the main-combustion chamber, the mass flow from the burned

15

zone of the pre-combustion chamber to the main-combustion chamber is neglected. The amount
of the fuel in the main-combustion chamber is calculated by


dmmain− f uel
1
= ṁtur−u
dt
Îť pre (A/F)s + 1

(2.16)

where mmain− f uel is the fuel mass in the main-combustion chamber.
After ignition in the main-combustion chamber, the burned zone is created. Different from
the two-zone combustion model in a conventional SI engine, the combustion model of the maincombustion chamber needs to consider the gas flowing through the orifice into the pre-combustion
chamber. The mass and energy conservation equations for burned and unburned zones are very
similar to those of the pre-combustion chamber model presented in this subsection and is omitted
here. The major difference is that the total volume of the main-combustion chamber is varying.

2.3.3

Mass fraction burned model

The MFB in the pre-combustion chamber is obtained from the Wiebe function [21].
#


t − tign m+1
xb = 1 − exp −a
∆td
"

(2.17)

The coefficients, a and m, are chosen to be 6.908 and 2, respectively; tign is the start of ignition;
and ∆td is the burn duration that is calibrated by AFR before ignition.
At the early stage of the combustion in the main-combustion chamber, the rate of combustion
is determined by not only the gas properties in the main-combustion chamber but also the turbulent
jet from the pre-combustion chamber. This is due to the fact that the turbulent jets create multiple
and distributed ignition sites, which increases the overall flame front area in the main-combustion
chamber. Moreover, these turbulent jets increase the turbulence intensity in the main-combustion
chamber and thus the flame front propagation speed. After the turbulent jet disappears, the rate of
combustion reduces gradually to a relatively low level and its characteristics are mainly determined
by the gas properties only in the main-combustion chamber. Here, the term ’intensity’ of the turbulent jet is used to describe the resulting increment of the combustion rate in the main-combustion
16

chamber due to the turbulent jet. Since the intensity of the turbulent jet is determined by the combustion processes in both combustion chambers, estimating the rate of combustion before ignition
is difficult and requires significant calibration effort. Therefore, adjusting the rate of combustion
according to the turbulent jet intensity during the combustion process is preferred for the TJI combustion model. The conventional single-Wiebe function is not suitable for our combustion model.
Multi-Wiebe function is a possible approach for modeling the MFB. However, it requires to determine all associated parameters before ignition occurs. Therefore, a new parameter-varying Wiebe
function is proposed and used in this chapter; see Eqn. (2.18).

( 
m+1 )

t
−
t
(t)

ign_b
0


xb (t) = 1 − exp −a
∆t
d

Z t




tign_b (t) = t0 − [b (t) − 1] dt

(2.18)

t0

where t0 and ∆td are determined by the spark timing and the AFR in the main-combustion chamber.
The coefficient a and m are chosen to be 6.908 and 2, respectively.
If a, m, and ∆td are the same in Eqn. (2.17) and Eqn. (2.18), it can be proved that for any given
tign_b = tign
dxb0
dx
= b ¡ b (t) .
dt
dt

(2.19)

In other word, the combustion rate calculated by the new Wiebe function is b(t) times larger
than that calculated by the conventional Wiebe function. Therefore, the intensity of the turbulent
jet can be mathematically expressed by b (t). The combustion model is able to adjust the rate of
combustion by making b(t) as a function of some characteristics of the turbulent jet. Moreover,
b(t) can be changed at any time during the combustion process. As long as b (t) is greater than 0,
the combustion rate is greater than 0 and xb0 (t) equals to 1 as t goes to infinity. From the available
experimental results, it is found that the rate of combustion in the main-combustion chamber is
highly related to the mass flow rate of the turbulent jets from the pre-combustion chamber. As
a result, it is assumed that the intensity of the turbulent jet can be linked to its mass flow rate.

17

Although this assumption provides a good match between the modeled and available experimental
results, it is important to find an accurate method to calculate the intensity of the turbulent jet in
the future when more experimental data are available. Since the influence of the turbulent jet to the
combustion in the main-combustion chamber is also delayed, b (t) is modeled to be proportional
to the flow rate of the turbulent jet with a first order dynamics; see Eqn. (2.20). An offset is used
such that b(t) = 1 when the flow rate of the turbulent jet is zero.

+ ∗ f )(t)] + 1
b (t) = β · [(ṁtur
l
+ ∗ f )(t) =
(ṁtur
l

where ” ∗ ” is the convolution operator

Z t

+ (τ) f (t − τ)dτ
ṁtur
l

0
+
and ṁtur and fl

(2.20)

are defined as follows



 ṁtur ṁtur ≥ 0
+
ṁtur =

 0 ṁtur < 0

(2.21)

fl (t) = ωc e−ωct u(t)
+ with the exponential decay function f
Note that u(t) denotes unit step. The convolution of ṁtur
l

represents the first order dynamics and is used to emulate the time delay. To be more specific, fl
is the time response of a low-pass filter with cutoff frequency ωc . The parameters β and wc are
experimentally determined. When ṁtur ≤ 0, we have b (t) = 1 which means that the combustion
rate will not be altered if there is no turbulent jet from the pre-combustion chamber.

2.4

Model Calibration

The combustion model was calibrated using the experimental data collected from the RCM
at Michigan State University described in the System Description section. The model is firstly
calibrated using two experimental data sets and validated using another data set. Then, to further
validated the model, more data sets with different pre-combustion chamber orifice sizes are used.
The experimental set-ups for the first three cases can be found in Table 2.2, where cases 1 and 2
are used for model calibration and case 3 for model validation.
18

Table 2.2: Experimental set-up.
Parameter
orifice diameter (mm)
main-chamber
AFR
pre-chamber fuel addition (mg)
pre-chamber AFR (calculated)

Case
1

2

3

4

5

6

1.5

1.5

1.5

2.0

2.0

2.5

1.83

2.10

1.83

1.5

1.25

1.5

0.66

0.89

0.89

0

0

0

0.98

0.90

0.85

1.5

1.25

1.5

The first step is to calibrate the heat transfer model. To do this, the net heat release (NHR) rate
in the main-combustion chamber needs to be calculated from an inverse thermodynamic calculation [18, 22] based on the experimental pressure data. However, this calculation cannot be completed without knowing the mass flow rate between the two combustion chambers. Fortunately,
the mass flow between the two combustion chambers mainly influence the combustion process in
the pre-combustion chamber. Its effect to the main-combustion chamber is limited. For calibration
purpose, the mass flow can be neglected when calculating the NHR rate in the main-combustion
chamber. The result of the inverse calculation is shown in the upper plot of Fig. 2.5. The NHR
rate is the sum of the CER rate and the heat transfer rate from the cylinder wall. The heat transfer
model can be calibrated assuming that the heat transfer rate is dominant where the NHR rate is
negative. The calculated heat transfer rate is shown in the bottom plot of Fig. 2.5.
The next step is to calibrate ∆td for the main-combustion chamber. After the heat transfer
model is calibrated, the CER rate can be obtained by subtracting the heat transfer rate from the
NHR rate, shown as the dotted line in Fig. 2.5 . This allows one to calculate the MFB of the
main-combustion chamber, which is the solid line in Fig. 2.6. According to the later stage of the
MFB curve (shown in Fig. 2.6), the parameter, ∆td , in the Wiebe function can be determined by
a linear Least-Squares fitting procedure[23]. The dotted line in Fig. 2.6 is the curve fitting result.
The parameter b (t) is set to be 1 during this calibration procedure.
Because the heat transfer coefficient hc of the pre-combustion chamber is assumed to be equal
19

Net Heat
Chemical Energy

Release Rate (kW)

300
250
200
150
100
50
0
Heat Transfer Rate (kW)

−50

0

5

10

15

20

25

30

0

5

10

15
Time (ms)

20

25

30

15
10
5
0

Figure 2.5: Inverse thermodynamic calculation results.

Mass Fraction Burned

1
Calculated
Curve fitting

0.8
0.6
0.4
←Least−Squares fitting→
0.2
0

0

5

10

15
Time (ms)

20

25

30

Figure 2.6: MFB obtained from inverse thermodynamic calculation and Least-Squares fitting
result.

20

to that of the main-combustion chamber. The heat transfer rate to the pre-combustion chamber
wall is determined. This allows one to calibrate the other unknown parameters in Table 2.3 using a
nonlinear Least-Squares optimization procedure. Among these parameters, c pre , cmain , ωc and β
in Eqn. (2.12) and Eqn. (2.20) remain constant for all the first three experimental cases. The results
of the linear Least-Squares fitting and nonlinear Least-Squares optimization procedures are shown
in Table 2.3.
The nonlinear Least-Squares optimization problem is solved by the nonlinear Least-Squares
solver in Matlab using the Trust-Region-Reflective algorithm. This algorithm minimizes the following function.
n

n

∑ (Ppre−i − P̂pre−i)2 + ∑ (Pmain−i − P̂main−i)2

i=1

(2.22)

i=1

where n is the total number of the data points; Ppre−i and Pmain−i are the experimental pressure
points; P̂pre−i and P̂main−i represent the modeled pressure points; i is the data index.
Table 2.3: Calibration results for cases 1 and 2.
Parameter

2.5

Case 1

Case 2

c pre
cmain
ωc
β

0.655
0.01
3263
1.72

0.655
0.01
3263
1.72

pre-chamber

∆td (ms)
Ρ pre

4.17
0.941

3.94
0.884

main-chamber

∆td (ms)
Ρmain

19.3
0.912

28.6
0.975

Model Validation and Simulation Results

After the calibration procedure, the model is then validated by the third experimental data set
listed in Table 2.2. The simulation parameters for the third experimental case are determined based
on the following simple assumptions. Since the AFR in the main-combustion chamber of the third
case is the same as the first one, ∆td and ηmain of the main-combustion chamber are assumed
21

the same as the first case. In the pre-combustion chamber, the other combustion parameters are
assumed to vary linearly with the AFR for the three experimental cases, because their AFRs in
the pre-combustion chamber are within a relatively small range. Based on the calculated AFRs in
Table 2.2, the coefficients, ∆td and η pre of the pre-combustion chamber are calculated and listed
in Table 2.4. The parameters, c pre , cmain , ωc and β , remain the same. Figs. 2.7-2.9 show the
comparison between the modeled and experimental pressure traces in two combustion chambers
for the three experimental cases. The calculated pressure traces for the first two cases match
the experimental pressure traces very well, since their parameters are obtained by the calibration
procedure. Although the combustion parameters for the third case are calculated based on the very
simple assumptions discussed above, the agreement between the modeled and measured pressure
traces is satisfactory. The relative errors on the pressure traces are always below 10%. The errors
mainly occur after 31 ms in Fig. 2.9. These errors are mainly caused by the simple assumptions
used for determining the simulation parameters for case 3. In reality, the combustion process in
the main-combustion chamber does not only depend on the parameters in the main-combustion
chamber. The relationship between the combustion parameters and the AFR in the pre-combustion
chamber is also not exactly linear. Once more experimental data are available, further calibration
can be done and model accuracy can be improved.
Table 2.4: Simulation parameters.
Parameter

Case 3
c pre
cmain
ωc
β

0.655
0.01
3263
1.72

pre-chamber

∆td (ms)
Ρ pre

3.80
0.848

main-chamber

∆td (ms)
Ρmain

19.3
0.912

One of the main differences between TJI combustion and conventional SI combustion is that the
turbulent jet increases the rate of combustion in TJI combustion. The experimental and calculated
22

Pressure (bar)

Pre−chamber

40
30
20
10
0

Main−chamber

40
30
20
10
0

Calculated
Experimental
0

5

10

15

20

25
Time (ms)

30

35

40

45

Figure 2.7: Experimental and calculated pressure traces for experimental case 1.

Pressure (bar)

Pre−chamber

40
30
20
10
0

Main−chamber

40
30
20
10
0

Calculated
Experimental
0

5

10

15

20

25
Time (ms)

30

35

40

45

Figure 2.8: Experimental and calculated pressure traces for experimental case 2.

23

Pressure (bar)

Pre−chamber

40
30
20
10
0

Main−chamber

40
30
20
10
0

Calculated
Experimental
0

5

10

15

20

25
Time (ms)

30

35

40

45

Figure 2.9: Experimental and calculated pressure traces for experimental case 3.

NHR rates are compared in Fig. 2.10. We can find a pick on the NHR rate curve at the early stage
of the combustion, which is caused by the turbulent jet. As an example, Fig. 2.11 shows how the
model simulates the rate of combustion according to the mass flow rate of the turbulent jet. The
top plot of Fig. 2.11 shows the mass flow rate through the orifice connecting the two combustion
chambers. The curve in the middle plot of Fig. 2.11 is b (t) in Eqn. (2.20). The bottom plot of
Fig. 2.11 is the calculated CER rate. According to the simulation results, the MFB model with the
parameter-varying Wiebe function successfully links the combustion processes in two combustion
chambers.
The value of the coefficient Îąb during the simulation is shown in Fig. 2.12, where the dashed
line is the value of Îąb calculated by assuming Ppre is greater than Pmain and the dash-dotted line
is the value of Îąb calculated by assuming Ppre is less than Pmain . The actual value of Îąb used
for combustion calculation , that is shown by the solid line, is on the dashed line at beginning
because the pre-combustion chamber is ignited first. After the switch line (see Fig. 2.12), the
main-combustion chamber pressure becomes larger than that of the pre-combustion chamber, Îąb
jumps from the dashed line to the dash-dotted line.
To further validate the model, twelve more experimental data sets were used. The orifice

24

Net Heat Release Rate (kW)
300
Experimental
Calculated

Case 1

200
100
0
200

Case 2

150
100
50
0
−50

Case 3

300
200
100
0
26

28

30

32

34
36
Time (ms)

38

40

42

44

Mass Flow Rate (g/s)

Figure 2.10: Experimental and calculated net heat release rates.

4
2
0
−2

b(t)

6
4

Chemical Energy
Release Rate (kW)

2
0
300
200
100
0
20

22

24

26

28

30
32
Time (ms)

34

36

38

Figure 2.11: Chemical energy release rate calculation.

25

40

1
← Ppre<Pmain
0.8

Îąb

0.6
← switch line
0.4
0.2
← Ppre>Pmain
0
24

25

26

27
28
Time (ms)

29

30

31

Figure 2.12: The calculated value of Îąb during combustion.

diameter varied from 2.0 mm to 3.0 mm. The relative AFRs in the two combustion chambers varied
from 0.9 to 1.5. To simplify the presentation, the pressure traces were plotted only for the first three
cases; see Figs. 13-15. The experimental set-ups for the three cases are shown in Table. 2, as cases
4-6. The associated calibration results can be found in Table. 5. For the other cases, the calculated
10-50% burn duration (Burn1050) and 50-90% burn duration (Burn5090) of the main-combustion
chamber were compared with the experimental values; see Fig. 16-17. Due to large variations of
orifice areas and AFRs of these data sets, the model parameters need to be re-calibrated. However,
β and wc were kept unchanged for the cases with the same orifice sizes, like cases 4 and 5. It was
also found that c pre and cmain were very similar for all the experimental cases. This indicates that
these two parameters are mainly associated with the combustion chamber structure. In Table. 5, the
pre-combustion chamber burn durations are quite different from the previous experimental cases
1-3. This is due to the difference in the orifice area.
To conclude, the proposed model is able to fit the experimental data sets with large ranges
of AFRs in both combustion chambers and different pre-combustion chamber orifice areas. This
indicates that the developed combustion model has the potential to be used for the development of
TJI engine model.

26

Pressure (bar)

Pre−chamber

50
40
30
20
10
0

Main−chamber

50
40
30
20
Calculated
Experimental

10
0

0

5

10

15

20
25
Time (ms)

30

35

40

Figure 2.13: Experimental and calculated pressure traces for experimental case 4.

Pressure (bar)

Pre−chamber

50
40
30
20
10
0

Main−chamber

50
40
30
20
Calculated
Experimental

10
0

0

5

10

15

20
25
Time (ms)

30

35

40

Figure 2.14: Experimental and calculated pressure traces for experimental case 5.

27

Table 2.5: Calibration results for cases 4-6.
Parameter

Case 4

Case 5

Case 6

c pre
cmain
ωc
β

0.651
0
3829
0.824

0.651
0
3829
0.824

0.651
0
5064
0.764

pre-chamber

∆td (ms)
Ρ pre

1.97
0.980

2.28
0.939

1.79
0.901

main-chamber

∆td (ms)
Ρmain

16.5
0.927

11.2
0.905

14.0
0.932

Pressure (bar)

Pre−chamber

50
40
30
20
10
0

Main−chamber

50
40
30
20
Calculated
Experimental

10
0

0

5

10

15

20
25
Time (ms)

30

35

40

Figure 2.15: Experimental and calculated pressure traces for experimental case 6.

2.6

Conclusion

This chapter presents a control-oriented TJI combustion model for the RCM (rapid compression machine) at Michigan State University. A newly proposed parameter-varying Wiebe combustion model is used to link the combustion processes in both pre- and main-combustion chambers.
The developed model can be calibrated using a simple and systematic calibration procedure based
on the experimental data. The model validation process shows a good agreement between the
modeled and experimental pressure traces, which indicates that the developed model is capable of
accurately capturing the TJI combustion dynamics. The validation results also indicate that the

28

Burn1050 (ms)
3

Calculated Result

2.5
2
1.5
1
0.5
0

0

0.5

1
1.5
2
Experimental Result

2.5

3

Figure 2.16: Experimental and calculated Burn1050 in the main-combustion chamber.
Burn5090 (ms)
3.5

Calculated Result

3
2.5
2
1.5
1
0.5
0

0

0.5

1
1.5
2
2.5
Experimental Result

3

3.5

Figure 2.17: Experimental and calculated Burn5090 in the main-combustion chamber.

model is able to predict the combustion process that is not used to calibrate the model parameters.
This shows that the developed model has the potential to be used for studying TJI combustion engines and developing the associated control strategies. Although only methane is used as the fuel
in this chapter, this model can be extended to other gaseous fuels with new calibrations. However,
for the liquid fuel, the model structure may need to be changed, especially a gas-fuel mixing model
will be required. The next step is to extend the modeling work for TJI engines using gaseous fuel.
Note that in this case the piston dynamics and gas exchange (intake and exhaust) models need to
be added. If the liquid fuel is used, an gas-fuel mixing model will also be required.

29

CHAPTER 3
A CONTROL-ORIENTED COMBUSTION MODEL FOR A TURBULENT JET
IGNITION ENGINE USING LIQUID FUEL

3.1

Introduction

The TJI combustion system introduces additional freedoms for combustion control. For instance, pre- and main-chamber AFRs, pre-chamber fuel injection and ignition timing. This makes
it difficult to utilize the traditional single variable feedback and look-up table based control, and it
often requires real-time optimization, especially during the transient operations. Therefore, modelbased control is required to handle multiple-input and multiple-output TJI combustion system.
Development of the control-oriented TJI combustion model serves two purposes. One is to quantitatively understand the influence of control parameters to TJI combustion performance and the
other is to develop the physics-based TJI combustion model for model-based control.
The engine model developed in this article is for a single cylinder engine equipped with a TJI
system. Different from many other TJI engines, liquid gasoline is used as the fuel for both the
pre- and main-chambers. Since liquid gasoline is hard to be fully vaporized inside the small prechamber before combustion starts, an air-fuel mixing and fuel vaporization model is proposed for
the pre-chamber in this article. In order to simulate the combustion performance variations due
to the interactions between the two combustion chambers, a parameter-varying Wiebe function is
used to model the main-chamber burn rate. The air-fuel mixture in the two combustion chambers
are modeled using two zones (burned and unburned) to simulate the complicated gas exchange
process between the two chambers.

30

Figure 3.1: TJI engine architecture

3.2
3.2.1

TJI Combustion Model Architecture
Target engine description

Fig. 3.1 shows the basic architecture of the single cylinder TJI engine modeled in this article.
The engine specifications are listed in Table 3.1. In order to have a good mixing inside the prechamber, most of TJI engines use gaseous fuel in the pre-chamber, which is not practical for
liquid fuel engines due to the required dual-fuel systems, leading to inconvenience, high system
complexity and cost. Therefore, liquid gasoline is used in both chambers for the TJI engine at
Michigan State University. Due to the small volume of the pre-chamber and the comparatively
long fuel penetration, a large amount of fuel will be sprayed directly onto the pre-chamber wall.
To control the AFR and improve the air-fuel mixing in the pre-chamber, an air injector, located in
the opposite side of the fuel injector, is used. When the fuel injector is turned on, the pressurized
air is also injected against fuel jets to reduce the amount of fuel impingement on the pre-chamber
wall. The spark plug is mounted at the top of the pre-chamber. The two combustion chambers are
connected through six orifices. The size of each orifice is listed in Table 3.1.

31

Table 3.1: TJI Engine Specifications
Parameter

Value

Bore

95 mm

Stroke

100 mm

Con-rod

190 mm

Compression ratio

11.5:1

Pre-chamber volume

2.7 cm3

Pre-chamber orifice diameter

1.5 mm

Intake/exhaust valve lifts

8.3 mm/8.3 mm

Figure 3.2: Key events over an engine cycle for TJI engine modeling. IVO: intake valve opening;
IVC: intake valve closing; EVO: exhaust valve opening; EVC: exhaust valve closing; SPK: spark
timing; and INJ: pre-chamber fuel and air injection.

3.2.2

Engine events and combustion model architecture

Fig. 3.2 illustrates the engine events and the corresponding modeling methods over one engine
cycle. For the current TJI engine configuration, the intake and exhaust valve timings are fixed
but the proposed model is also suitable for engines with variable valve timing. During the intake
stroke, fresh air and fuel flow (through port-fuel-injection) into the main-chamber and mix with the

32

residual gas inside the main-chamber. It is assumed that the air, fuel and residual gas are uniformly
mixed. Therefore, the gas inside the main-chamber is considered as one single zone, and one-zone
gas exchange model is used to calculate the gas properties in the main-chamber. After intake valve
closing (IVC), the piston compresses the gas inside the main-chamber and thus the pressure in the
main-chamber increases. The pressure difference between the two combustion chambers pushes
the gas from the main-chamber into the pre-chamber. The pressure in the pre-chamber increases
and follows the main-chamber pressure. During this process, the gas inside the pre-chamber is also
considered as a single zone. To control the AFR and residual gas contents inside the pre-chamber,
additional fuel and air are injected into the pre-chamber during the compression stroke, where air
injection is also used during the fuel injection to reduce the fuel wall-wetting effect in the prechamber. The fuel and air injection (INJ) timings are shown in Fig. 3.2. After INJ, the fuel mixing
and vaporization model calculates the vaporized fuel amount inside the pre-chamber based on the
amount of injected fuel and the fuel vaporization rate on the pre-chamber wall. When the piston
moves close to the combustion top dead center (CTDC), the spark plug discharges and ignites the
mixture inside the pre-chamber. Then, the reacting products in the pre-chamber are pushed out and
ignite the mixture in the main-chamber. After the ignition in the pre-chamber, the gas inside the
combustion chambers is divided into two zones, burned and unburned zones. In the pre-chamber,
the fuel vaporization model calculates the fuel masses vaporized from the pre-chamber wall for
both zones. A single step reaction based model is used to calculate the combustion rate inside the
burned zone in the pre-chamber. In the main-chamber, the combustion rate is calculated by the
parameter-varying Wiebe function. The modeling details are discussed in the next section.

3.3

Gas Exchange Model

During the exhaust, intake and compression strokes before ignition, the engine is modeled by
the one-zone model. The following assumptions are made:
1. Fuel and air are uniformly premixed in the intake manifold.
2. The gas inside the pre- or main-chamber is homogeneous, and thus can be considered as one
33

zone. However, the properties of the gas in one chamber can be different from that in the
other chamber.
During the intake stroke, the mass flow rate of the intake fresh charge can be calculated by the
one-dimensional compressible flow equation [1].
P
ṁint = Cd,int Av,int √ int ψ
RTint
where



Pmain
Pint



v
u

 (Îş+1)
 Îş

u 


t
2
2
(κ−1)
(κ−1)


x<
 κ
Îş +1
Îş +1
ψ (x) =
s


 Îş




(κ−1)
1
2Îş
2
(κ−1)


1−x κ
x≥
x κ
κ −1
Îş +1

(3.1)

(3.2)

Note that Av,int is the intake valve reference area; the discharge coefficient Cd,int is experimentally determined; Îş is the ratio of specific heats; R is the gas constant; Pint and Tint are the upstream
pressure and temperature in the intake manifold, respectively; and Pmain is the main-chamber pressure .
During the exhaust stroke, the mass flow rate of the gas flowing throw the exhaust valve, ṁexh ,
can be calculated similarly.
The gas exchange process between the two combustion chambers is also modeled in the similar
way. However, the pressure in the pre-chamber can be either greater or less than that in the mainchamber. Therefore, the following equation is used to calculate the mass flow rate between the two
combustion chambers.
ṁtur =



Ppre
Pmain


Ppre ≥ Pmain

 Cd,tur Av,tur pRT ψ P
pre
pre



Ppre
P

main

ψ
Ppre < Pmain
−Cd,tur Av,tur √
Pmain
RTmain

(3.3)

where Cd,tur and Av,tur are the discharge coefficient and the area of the orifice connecting the two
combustion chambers; respectively. Note that the subscripts pre and main denote the pre- and
main-chamber properties, respectively.
34

Before ignition, the main-chamber is considered as a control volume with mass and energy
exchange. The following mass and energy conservation equations describe such a control volume.
dmmain
= ṁint − ṁexh + ṁin j − ṁtur
dt
d(mmain umain )
(3.4)
= ṁint hint − ṁexh hexh
dt
dW
+ ṁtur htur − Q̇main
ht − dt
where mmain and umain are the mass and internal energy of the gas in the main-chamber, respecis the heat transfer rate through the chamber wall; and ṁh is the enthalpy flow. The
tively; Q̇main
ht
work of the pressure force, W , is calculated by:
dW
dV
= Pmain main
dt
dt

(3.5)

where Vmain is the main-chamber volume. The Woschni correlation model [20] is used to calculate
the quantity of heat transfer.
Q̇main
= Amain
w hc (Tmain − Tw )
ht

(3.6)

−0.55
hc = qB−0.2 P0.8 (C1V̄p )0.8 Tmain

where V̄p is the mean piston velocity and C1 = 2.28; B is the cylinder bore; Ac is the main-chamber
wall area; and coefficient q is a calibration parameter.
The mass and energy conservation equations for the pre-chamber can be expressed as:
dm pre
= ṁin j − ṁtur
dt


(3.7)
d m pre u pre
pre
= ṁin j hin j − ṁtur htur − Q̇ht
dt
where m pre and u pre are the mass and internal energy of the gas in the pre-chamber, respectively;
pre

Q̇ht is the heat transfer rate through the pre-chamber wall. The subscripts in j and tur represent
the properties of the gas from the fuel and air injectors and through the orifice connecting the two
combustion chambers, respectively.
Assuming that the gas can be considered to be ideal, the above equations can be coupled using
the ideal gas law as follows.
P ¡V = m ¡ R ¡ T
35

(3.8)

Finally, the pressures and temperatures in both of the two combustion chambers can be obtained
by substituting Eqn. (3.8) into Eqns. (3.4) and (3.7).

3.4

Fuel Mixing and Vaporization Model

As discussed in the last section, although an air injector is used in the pre-chamber, it is difficult
to completely prevent the injected liquid gasoline from shooting onto the chamber wall. Therefore,
a fuel mixing and vaporization model is developed based on the following assumptions:
1. The injected fuel can be divided into two parts: vaporized and liquid fuel. The vaporized
fuel uniformly mixes with the gas in the pre-chamber immediately after injection.
2. The injected liquid fuel forms a fuel film on the chamber wall, and the area fraction of the
fuel film contacting with the burned (or unburned) zone is the same as the volume fraction
of the burned (or unburned) zone.
3. The fuel film surface temperature contacting the gas mixture is assumed to be the same as
the gas temperature. The actual temperature should be lower, but this modeling error can be
absorbed by the coefficient hm in Eqn. (3.11).
The amounts of vaporized and liquid fuel are simply calculated by the following equations
based on the total fuel injected into the pre-chamber.
pre

pre

ṁ f ,wall = k · ṁ f ,in j
pre
ṁ f ,mix

pre
= (1 − k) · ṁ f ,in j

(3.9)

pre

where ṁ f ,in j is the pre-chamber fuel injection rate and the superscript pre represents the prepre

pre

chamber properties; ṁ f ,wall is the rate of the fuel deposited on the pre-chamber wall; and ṁ f ,mix is
the rate of the injected fuel that directly vaporises and mixes with the gas inside the pre-chamber.
The fuel deposited on the pre-chamber wall forms a thin fuel film and continue vaporizing, see
Fig. 3.3. The mass of the fuel film can be calculated by
pre

pre

pre

ṁ f , f ilm = ṁ f ,wall − ṁ f ,evap
36

(3.10)

Figure 3.3: Pre-chamber fuel vaporization model.
pre

pre

where m f , f ilm is the mass of the fuel film on the chamber wall; and ṁ f ,evap is the fuel vaporization
rate from the fuel film. It is calculated based on the net molar diffusion flux of the fuel, N f .
According to Fick’s law:

pre

Nf =

ṁ f ,evap
Mf Af

= hm (c f , f ilm − c f ,∞ )

(3.11)

where hm is the mass transfer coefficient; M f is the average molecular weight of the fuel; A f is
the surface area of fuel film exposed to the pre-chamber gas; and c f , f ilm and c f ,∞ are the fuel
vapor concentrations near the fuel film surface and in the pre-chamber gas, respectively. Since
c f , f ilm  c f ,∞ , c f ,∞ is set to zero and c f , f ilm can be derived according to the ideal gas equation
as follows
c f , f ilm =

Ps
Rm T f

(3.12)

where Ps is the fuel saturated vapor pressure; T f is the fuel film temperature which is assumed to
be the pre-chamber gas temperature; and Rm is the molar gas constant. The saturated fuel vapor
pressure can be predicted by the Clausius-Clapeyron equation as
Ps = Ps0 e

∆H
Rm T f

(3.13)

where Ps0 is the fuel saturated vapor pressure constant; and ∆H is the fuel molar enthalpy of
vaporization. According to Eqns. (3.11)-(3.13), the fuel film vaporization rate can be derived as
pre

pre
ṁ f ,evap

m f , f ilm

∆H

P
R T
=
M f hm s0 e m f
δf ρf
Rm T f
37

(3.14)

Figure 3.4: Two-zone combustion model. Îąb: area fraction for burned gas

where δ f is the fuel film thickness and ρ f is the density of the liquid fuel.
After ignition, the pre-chamber is divided into two zones. Correspondingly, the fuel vaporization flow is also divided into two parts: one to the burned zone and the other to the unburned
zone; see Fig. 3.3. It is assumed that the mass of the fuel film exposed to the burned zone is
pre pre

pre

vb m f , f ilm , where vb

is the volume fraction of the burned gas in the pre-chamber. Therefore, the
pre

pre pre

fuel vaporization rate to the burned zone can be calculated by replacing m f , f ilm with vb m f , f ilm
in Eqn. (3.14) and using the burned zone temperature for T f . The fuel vaporization rate to the
unburned zone is calculated similarly.

3.5

Two-zone Combustion Model

After ignition, the two combustion chambers are divided into burned and unburned zones to
improve the model accuracy. The following assumptions are made in the two-zone combustion
model:
1. The gas mixture in each zone is homogeneous.
2. After entering the burned zone, the unburned gas burns immediately and releases the chemical energy inside the burned zone.
3. The heat transfer between the burned and unburned zones can be negligible.

38

4. The pressure of the gas mixture inside each combustion chamber is assumed to be evenly
distributed through the burned and unburned zones.
The burned and unburned zones can be considered as control volumes. Besides the mass,
enthalpy, and work exchange between the two zones in each chamber, there are also mass and
enthalpy exchanges through the orifice connecting the two combustion chambers. When both the
unburned and burned zones are present in the two chambers, the gas in the burned (or unburned)
zone of one chamber are assumed to enter the burned (or unburned) zone of another chamber as
shown in Fig. 2.3. Before the ignition of the main-chamber, all the gas from the main-chamber is
considered as unburned gas. After all the unburned gas in the pre-chamber becomes burned gas,
all the gas coming from the main-chamber goes into the burned zone in the pre-chamber. The
energy balance equation of the burned zone is shown in Eqn. (3.15). The following equations in
this subsection considered the case when both burned and unburned zones are present, and the
variables in these equations are for the pre-chamber, if not specified.

dV
d(mb ub )
+ P b + vb Q̇ht =
dt
dt

(3.15)

Q̇ch + ṁu−b hu − ṁtur,b hb
The energy balance equation of the unburned zone can is expressed by
d(mu uu )
dVu
+P
+ (1 − vb ) Q̇ht =
dt
dt

(3.16)

− ṁu−b hu − ṁtur,u hu
The masses of burned and unburned gases are calculated based on the following mass conservation
law.
dmb
= −ṁtur,b + ṁu−b
dt
dmu
= −ṁtur,u − ṁu−b
dt

(3.17)

For the above three equations, Qht is the heat transfer to the chamber wall; Qch is the chemical
energy released by combustion; subscripts b and u represent properties in burned and unburned
39

zones, respectively; vb is the volume fraction of the burned zone; ṁu−b is the mass transfer rate
from the unburned zone to the burned zone and it will be discussed in the next section. The gas
flowing through the orifice contains both burned gas, mtur,b , and unburned gas, mtur,u . Correspondingly, the area of the orifice is also divided into two parts: one for the burned gas and the other for
the unburned gas, see Fig. 2.3. The two mass flow rates are calculated by the following equations
when the pressure in the pre-chamber is greater than that in the main-chamber.


Ppre
Pmain
ψ
ṁtur,b = αbCd,tur Av,tur √
RTb
Ppre


Ppre
Pmain
ψ
ṁtur,u = (1 − αb )Cd,tur Av,tur √
Ppre
RTu

(3.18)

Similar results can be obtained when the pressure in the main-chamber is greater than that in the
pre-chamber.
The coefficient Îąb in Eqn. (3.18) is calculated based on the volume fraction of the burned
gas vb . Assuming that all the burned and unburned gas was evenly distributed in the combustion
chambers, Îąb would be always equal to vb . However, in reality, Îąb depends on the combustion
chamber structure. For the current TJI engine setup, the spark plug is located at the top of the
pre-chamber; see Fig. 3.1. As a result, the combustion in the pre-chamber is initiated at the top of
the chamber. Therefore, the burned gas is mainly located at the top of the pre-chamber at the early
stage of the combustion, which makes it hard for the burned gas to escape through the orifice to the
main-chamber. Consequently, the fraction of the burned gas flowing to the main-chamber is much
smaller than the burned gas fraction inside the pre-chamber when the gas flow direction is from
the pre-chamber to the main-chamber. This is why Îąb is smaller than vb in the pre-chamber when
Ppre > Pmain . When Ppre < Pmain , the gas flow changes its direction and Îąb will be determined
by the burned gas fraction in the main-chamber. Since the combustion in the main-chamber is
initiated at the end of the turbulent jet (far away from orifice)[24], the orifice is surrounded by the
gas with high concentration of unburned gas. Therefore, Îąb is smaller than vb in the main-chamber.
This is the main reason why the two-zone combustion model is used. To simulate this relationship
between Îąb and vb , the value of Îąb is calculated by Eqn. (3.19).

40

f1 (Ppre>Pmain)

f (Ppre<Pmain)
2

1

0.8

0.8

0.6

0.6

0.4

0.4

Îąb

1

0.2

0.2
(c

0

,1−c

pre

0

0.5

)

pre

0

1

pre

(cmain,1−cmain)
0

0.5

1

main

vb

vb

Figure 3.5: The value of Îąb in the two cases.





 f1 v pre , c pre
Ppre > Pmain
b
Îąb =



 f vmain , cmain Ppre < Pmain
2 b

(3.19)

where the two functions, f1 and f2 , are approximated by the second-degree BĂŠzier curves [19],
to simplify the calibration process. Each curve is defined by three control points: (0,0) and (1,1)
in both curves; and (c pre ,1 − c pre ) in f1 , (cmain ,1 − cmain ) in f2 ; see Fig. 2.4. The Bézier curve
guarantees ιb ∈ [0, 1] as long as c pre ∈ [0, 1] and cmain ∈ [0, 1]. When c pre and cmain are larger
than 0.5, Îąb is smaller than vb during combustion; as shown in Fig. 2.4. The values of c pre and
cmain are experimentally determined.
The main-chamber is also divided into two zones after ignition. Different from the two-zone
combustion model in a conventional SI engine, the combustion model of the main-chamber needs
to consider the gas flow through the orifice, especially the unburned fuel from the pre-chamber.
The governing equations for the main-chamber are very similar to those of the pre-chamber model
presented in this subsection and thus are omitted here. The major difference is that the total volume
of the main-chamber is varying.

41

3.5.1

Pre-chamber combustion rate calculation

The fuel inside the pre-chamber is composed of two parts: vaporized fuel and liquid fuel on the
chamber wall. During combustion, the fuel film continues vaporizing into the gas mixture inside the chamber. Therefore, the combustion process is more complicated than that in the mainchamber. After the spark is discharged, the burned zone is created. As the flame front propagates
inside the pre-chamber, the fuel in the unburned zone burns and the unburned gas becomes burned
gas. It is assumed that all the fuel in the unburned zone burns completely during this process. As
a result of combustion, the temperature inside the pre-chamber increases and the fuel vaporization
rate also goes up. Because the AFR in the pre-chamber before combustion is close to stoichiometry, there is usually no or little oxygen left inside the burned zone. Therefore, the concentration
of the unburned fuel in the burned zone keeps increasing. As all the unburned gas is burned, the
pressure in the pre-chamber starts dropping. At this point, the combustion in the main-chamber
has already started and the pressure increases rapidly. When the pressure in the main-chamber is
greater than that in the pre-chamber, the gas with high concentration of air from the main-chamber
enters the pre-chamber burned zone. Then the unburned fuel reacts with the air and the pressure in the pre-chamber increases and becomes the same or even higher than the pressure in the
main-chamber.
The mass transfer rate of the gas from the unburned zone to the burned zone is obtained from
the Wiebe function[21]:
#

t − tign m+1
xb = 1 − exp −a
∆td
"



(3.20)

where xb is the mass fraction of the burned gas; the coefficients, a and m, are chosen to be 6.908
and 2, respectively; tign is the start of ignition; and ∆td is the burn duration that is a function of the
AFR right before ignition.
The mass transfer rate of the gas from unburned to burned zone, ṁu−b , can be calculated by:
ṁu−b =

mu dxb
1 − xb dt

42

(3.21)

where mu is the mass of the unburned gas. Because the total mass inside the pre-chamber is
changing during the combustion, ṁu−b cannot be calculated directly by the total mass.
The rate of chemical energy release in Eqn. (3.15) can be obtained by the following relationship.
pre

Q̇ch = η pre QLHV α f ,u ṁu−b

(3.22)

where the combustion efficiency Ρ pre is experimentally determined; QLHV is the lower heating
pre

value of the fuel; and Îą f ,u is the mass fraction of the fuel in the unburned zone of pre-chamber.
The fuel combustion rate in the burned zone is determined not only by the fuel amount but also
by the air amount and other in-cylinder properties, such as the gas temperature. It is not trustworthy
to use a conventional Wiebe function to predict the chemical energy release rate. Therefore, a
chemical kinetics based model is used to calculate the rate of reaction [25].
d[C8 H18 ]
= −Aexp(−EA /Ru T )[C8 H18 ]m [O2 ]n
dt

(3.23)

where [C8 H18 ] and [O2 ] are the molar concentrations of the fuel and oxygen; A is the pre-exponential
factor; and EA /Ru is the activation temperature. m and n are chosen to be 0.25 and 1.5, respectively.
Then the chemical energy release rate can be calculated.

3.5.2

Main-chamber combustion rate calculation

During the TJI engine combustion process, both gas properties in the main-chamber and the turbulent jets from pre-chamber determine the rate of combustion in the main-chamber. The reason
is that the turbulent jets are able to create multiple and distributed ignition sites that increase the
overall flame front area in the main-chamber. At the same time, the turbulence intensity in the
main-chamber is also increased by the turbulent jets and thus the flame front propagation speed
goes up. After the turbulent jets disappear, the rate of combustion reduces to a relatively low level,
which is mainly determined by the gas properties in the main-chamber. To describe the increment
of the combustion rate in the main-chamber caused by the turbulent jets, the term ’intensity’ of the
turbulent jet is introduced.
43

Because the turbulent jets are able to change the rate of combustion in the main-chamber, the
conventional single-Wiebe function is not suitable for our combustion model. Using multi-Wiebe
function is a possible approach for calculating the mass fraction burned. However, this approach
requires all the associated parameters to be determined before ignition occurs, which requires
significant calibration effort and is not able to model the interactions between the two combustion
chambers. Therefore, adjusting the rate of combustion according to the turbulent jet intensity
during the combustion process is preferred for the TJI combustion model. Correspondingly, a new
parameter-varying Wiebe function is proposed and used in this chapter; see Eqn. (3.24).

( 
m+1 )

t
−
t
(t)

ign_b
0


xb (t) = 1 − exp −a
∆t
d

Z t




tign_b (t) = t0 − [b (t) − 1] dt

(3.24)

t0

where t0 is the start of ignition; t0 and ∆td are determined by the spark timing and the AFR in the
main-chamber. The coefficient a and m are chosen to be 6.908 and 2, respectively.
If a, m, and ∆td are the same in Eqns. (2.17) and (2.18), it can be proved that for any given
tign_b = tign
dxb0
dx
= b ¡ b (t) .
dt
dt

(3.25)

In other words, the rate of combustion calculated by the parameter-varying Wiebe function is
b(t) times larger than that calculated by the conventional Wiebe function. Therefore, the parameter
b(t) is actually the mathematical expression of the intensity of the turbulent jets. Moreover, b(t)
can be changed at any time during the combustion process. As a result, the parameter-varying
Wiebe function makes it possible for the combustion model to adjust the rate of combustion by
making b(t) as a function of some characteristics of the turbulent jets. According to the available
experimental results, it is found that the rate of combustion in the main-chamber is highly related
to the mass flow rate of the turbulent jets from the pre-chamber. So, it is assumed that the intensity
of the turbulent jet can be linked to its mass flow rate. Although this is a simple assumption, it

44

provides a good match between the simulated and experimental results. Other more complex models can also be developed to calculate the intensity of the turbulent jets more accurately. However,
this method is still a preferred one due to its simplicity that is very important for a control-oriented
model. Since the influence of the turbulent jets to the combustion in the main-chamber is also
delayed, b (t) is modeled to be proportional to the mass flow rate of the turbulent jets with the first
order dynamics; see Eqn. (3.26). An offset is used such that b(t) = 1 when the flow rate of the
turbulent jets is zero.

+ ∗ f )(t)] + 1
b (t) = β · [(ṁtur
l
+ ∗ f )(t) =
(ṁtur
l

Z t
0

+ (τ) f (t − τ)dτ
ṁtur
l

(3.26)

+ and f are defined as follows
where ” ∗ ” denotes the convolution operator and ṁtur
l



 ṁtur ṁtur ≥ 0
+
ṁtur =

 0 ṁtur < 0

(3.27)

fl (t) = ωc e−ωct u(t)
+ with the exponential decay
where u(t) denotes the unit step function. The convolution of ṁtur

function fl represents the first order dynamics and is used to emulate the time delay. To be more
specific, fl is the time response of a low-pass filter with a cutoff frequency of ωc . The parameters
β and ωc are experimentally determined. b (t) equals to one when ṁtur ≤ 0, which means that the
combustion rate will not be altered if there is no turbulent jets from the pre-chamber.

3.6

Model Calibration and Validation

To calibrate and validate the developed model, experiments were conducted on the single cylinder optical TJI engine, shown in Fig. 3.6, at MSU. The intake and exhaust valve timings are fixed
and shown in Fig. 3.2. The crank shaft is connected to the dynamometer. Two Kistler pressure
sensors are mounted in the two combustion chambers for capturing both chamber pressures. The
45

Figure 3.6: TJI engine at MSU.

pressure and other signals coming from the sensors are conditioned and recorded by a Phoenix
Combustion Analysis System (CAS) provided by the A&D Technology.
To calibrate the model, experiments were conducted at three different engine speeds. At each
engine speed, four different engine load conditions were used during the experiment. Upon calibration, the model was then validated by a set of experiments with different operational conditions
shown in Fig. 3.7, where the y-axis is the indicated mean effective pressure (IMEP). For the current
TJI engine, the maximum achievable air-fuel equivalence ratio, Îť , is about 2 with stable combustion. Note that the criterion of stable combustion is the coefficient of variation of IMEP to be less
than 3%. The data points around the dash-dot line in Fig. 3.7 corresponds to the engine operational
condition with wide open throttle (WOT) and λ ≈ 2. Above the dash-dot line, the engine runs with
WOT and Îť < 2. Below the dash-dot line, the amount of air flowing into the cylinder needs to be
reduced by closing the throttle to obtain a stable combustion. Besides the validation data shown
in Fig. 3.7, other experimental data with different pre-chamber fuel injection amounts are used to
further validate the model. Since the engine is an optical single-cylinder engine, the maximum engine speed is limited by its rotational balancing. Therefore, the maximum engine speed is limited
to 2000 rpm in the experiments.
The model is calibrated using the following steps. First, the discharge coefficients, Cd,int ,

46

8
7.5

Calibration data
Validation data

IMEP (bar)

7
6.5
6
5.5
5
4.5
4
1000

1200

1400
1600
1800
Engine speed (rpm)

2000

2200

Figure 3.7: Experiment matrix.

Cd,exh and Cd,tur , in Eqns. (2.1) and (2.3) are calibrated by matching the simulated pressure
traces with several experimental pressure traces obtained by motoring the engine. Second, the
heat transfer coefficient, q, in Eqn. (3.6) is calibrated using the similar method stated in [18].
This method first calculates the net heat release (NHR) rate by an inverse thermodynamic calculation and then calibrates the heat transfer model by assuming that the heat transfer is dominant
where the NHR rate becomes negative from the end of combustion to EVO. Note that NHR includes both the chemical energy and heat transfer to the chamber wall. Third, parameters, c pre in
Eqn. (2.12), k in Eqn. (3.9), and ∆td in Eqn. (2.17), are calibrated by matching the pre-chamber
pressure during the early stage of combustion (-25 to 0 deg ACTDC in Fig. 3.8). Then, parameters, cmain in Eqn. (2.12), δ f in Eqn. (3.14), A and EA /Ru in Eqn. (3.23), are adjusted according
to the pre-chamber pressure in the later stage of the combustion. After the parameters in the prechamber model are calibrated, the main-chamber parameters can be calibrated by matching the
main-chamber simulated pressure curve to be the same as the experimental one. Table. 3.2 shows
all the model calibration parameters. Most of these parameters are constant within the entire engine operational range. Some of them are functions of other parameters, like engine speed and
intake manifold pressure, Pint . These functions are implemented by look-up tables.

47

Table 3.2: Calibration Parameters
Parameter

Symbol

Function of/Value(if constant)

Intake valve discharge coefficient

Cd,int

0.9

Exhaust valve discharge coefficient

Cd,exh

0.9

Pre-chamber orifice discharge coefficient

Cd,tur

0.85

Heat transfer coefficient

q

389

Thickness of the fuel film

δf

0.03

BĂŠzier curve coefficients

c pre , cmain

0.8, 0.9

k

0.4

Pre-chamber Wiebe function coefficient

∆td

engine speed, Pint , Îť pre

Pre-chamber combustion efficiency

Ρ pre

engine speed, Pint , Îť pre

Pre-exponential factor

A

1.7 ¡ 1013

Activation temperature

EA /Ru

engine speed, Pint

K

∆td

engine speed, Pint , Îťmain

ms

Coefficient of parameter-varying Wiebe function

β

engine speed, Pint , Îťmain

Coefficient of parameter-varying Wiebe function

ωc

10000

Portion of fuel injected on pre-chamber wall

Main-chamber Wiebe function coefficient

48

Unit

mm

ms
(mol/cm3 )1−m−n /s

rad/s

Main−chamber pressure
Pre−chamber pressure

Pressure (bar)

50
40
←IMEP=7.1bar
30
20

Error (bar)

10
0
−50
5
0
−5
−50

IMEP=5.4bar→

0

0

50

100

150

200

50
100
150
Crank position (deg ACTDC)

200

Figure 3.8: Pressure traces obtained at 1500 rpm.

Fig. 3.8 shows the simulation results at 1500 rpm with IMEP of 5.4 and 7.1 bar. These two
operational conditions are marked in Fig. 3.7, where one is used for calibration and the other for
model validation. The upper plot in Fig. 3.8 shows the pressure traces in the two combustion chambers obtained from the developed TJI engine model. The bottom plot shows the simulation errors.
The error is always less then 1.8 bar and the relative error is less than 8%. This demonstrates that
the developed TJI engine model is capable of capturing the TJI engine combustion characteristics.
Fig. 3.9 shows the heat release rates in the two combustion chambers during combustion. The
upper plot shows the main-chamber NHR rates obtained from simulation (dashed line) and experiment (solid line). The experimental NHR is calculated from an inverse thermodynamic calculation
[18, 22] based on the experimental pressure data. The dash-dot line is the coefficient b(t) in
Eqn. (3.24). The bottom plot shows the chemical energy release (CER) rates in the pre-chamber.
The solid line is the CER from the fuel that vaporized immediately after injection. The dashed line
is the CER coming from the fuel vaporized from the fuel film. This part of CER cannot increase
rapidly before the main-chamber pressure is greater than the pre-chamber pressure, because the
vaporized fuel in the burned zone in pre-chamber burns rapidly only after a large amount of air is

49

3
Experimental HRR
Calculated HRR
2
b(t)

4
3
2
1
0

b(t)

Heat release rate (W)

5

x 10

1

−40

−20

0

20

40

60

80

100

120

Heat release rate (W)

4

4

x 10

Vaporized fuel
Fuel film

3
2
1
0
−1
−40

−20

0

20
40
60
80
Crank position (deg ACTDC)

100

120

Figure 3.9: Heat release rates in two chambers.

pushed into the pre-chamber from the main-chamber. It can also be found that the main-chamber
NHR rate is increased by the turbulent jets during the early stage (-20 to -10 deg ACTDC) and later
stage (0 to 10 deg ACTDC) of the combustion, where b(t) > 1. And these turbulent jets are caused
by the high CER rate in the pre-chamber at the corresponding two stages. In other words, without
the pre-chamber fuel vaporization model, the main-chamber NHR rate in the simulation cannot
match the experiment results. This confirms the necessity of the pre-chamber fuel vaporization
model.
Fig. 3.10 presents the value of coefficient Îąb in Eqn. (2.12); see the solid line. There are two
switch points around -11 and 0 deg ACTDC on this line, corresponding to the two points when the
mass flow between the two chambers changes its direction. The mass flow rate, ṁtur , is shown as
the light colored dash-dot line. Îąb is always smaller than the corresponding burned zone volume
fraction,vb , in the pre-chamber or in the main-chamber. This is consistent with the discussion in
Section 3.5 that the burned gas is hard to flow through the orifice between the two chambers.
Fig. 3.11 shows the simulation results with different pre-chamber fuel injection amounts. In

50

40

0.8

20

0.6

0

vbpre

0.4

−20

vbmain
Îąb
ṁtur

0.2

−40

0
−30

Mass flow rate (g/s)

vb and Îą b

1

−60
−20

−10
0
10
20
Crank position (deg ACTDC)

30

40

Figure 3.10: The value of Îąb and other related parameters during simulation.

5

Îť=1.20
Îť=1.09
Îť=1.00
Îť=0.87

30
20
10

5
4
3
2

←NHR
0

1

−10

0

−40
5
0
−5
−10
−40

−20

−20

0

20

40

60

0
20
40
Crank position (deg ACTDC)

60

Heat release rate (W)

40 Pre−chamber pressure→

5

x 10
2
1
0
−1

Error (W)

Error (bar)

Pre−chamber pressure (bar)

x 10

Figure 3.11: Pre-chamber pressure traces and main-chamber heat release rates for different
pre-chamber Îť .

51

Error (deg)

Burn duration (deg)

35
30

Spk−CA50
Spk−CA10

25
20
15
10
0.85
2
0
−2
0.85

0.9

0.95

0.9

0.95

1

1.05

1.1

1.15

1.2

1.25

1
1.05
1.1
1.15
Pre−chamber calculated λ

1.2

1.25

Figure 3.12: Main-chamber burn duration with the same main-chamber properties and different
pre-chamber Îť .

the upper plot, it can be found that the pre-chamber pressure rise is high for the case with rich prechamber AFR. High pressure rise generates the turbulent jets with large energy. Correspondingly,
the NHR rate is also high for the case with small pre-chamber Îť . This means that the main-chamber
burn rate increases as the pre-chamber fuel injection quantity is increased. As the pre-chamber Îť
becomes smaller than 1, adding more fuel will not increase the main-chamber burn rate much since
the pre-chamber is already very rich. This phenomenon can be shown more clearly in Fig. 3.12.
The upper plot presents how the main-chamber burn duration from spark timing (Spk) to CA50,
the crank angle when 50% of fuel is burned, is affected by the pre-chamber Îť . The bottom plots
in Fig. 3.11 and 3.12 show the simulation errors for each operational conditions. It can be found
that the errors are within a fairly small range, which indicates the developed TJI engine model is
able to simulate this unique phenomenon for the TJI combustion. This is very important because
it means that the model has the potential to be used for developing the corresponding model-based
control strategies using pre-chamber fuel injection to optimize the main-chamber combustion.
Fig. 3.13 compares the IMEPs obtained from the TJI engine model and those from the experiments with the validation operational conditions shown in Fig. 3.7. The two dash-dot lines on each
side of the data points show where the 5% error is. All of the IMEPs calculated by the model match
52

IMEP

Main−chamber peak pressure
55
50

6

Simulation (bar)

Simulation (bar)

7

5

4

3

45
40
35

3

4

5
6
Experiment (bar)

30
30

7

35

40
45
Experiment (bar)

50

55

Figure 3.13: IMEP comparison between simulation and experimental results.
Spk−CA50

Burn1090
25
Simulation (deg)

Simulation (deg)

30
25
20
15

20

15

10
10
10

15

20
25
Experiment (deg)

30

10

15
20
Experiment (deg)

25

Figure 3.14: Spk-CA50 comparison between simulation and experimental results.
Pre−chamber first peak pressure

Crank position of pre−chamber first peak pressure
20
Simulation (deg BCTDC)

Simulation (bar)

35

30

25

20

15
15

20

25
30
Experiment (bar)

35

15

10

5

0

0

5
10
15
Experiment (deg BCTDC)

Figure 3.15: Pre-chamber first peak pressures and their crank positions.

53

20

the experimental results very well. To further validate the model, other pressure related variables of
the main-chamber calculated by the model are compared with experimental results, including the
burn duration from Spk to CA50, 10-90% burn duration (Burn1090) and peak pressure. Finally,
Fig. 3.15 compares the experimental and simulation results of the first peak pressures in the prechamber (around 12 deg BCTDC in Fig. 3.11) and their points of occurrence (in crank position).
It is seen that the TJI engine model is able to reproduce the experimental results very well.

3.7

Conclusions

A control-oriented combustion model is proposed for a single cylinder turbulent jet ignition
(TJI) engine using liquid fuel for both combustion chambers. Several techniques are used to keep
the developed model simple enough for real-time simulations with satisfactory accuracy. The
model was calibrated using experimental data at engine speeds between 1200 and 2000 rpm with
indicated mean effective pressure between 4.2 to 7.2 bar, and then, the calibrated model is validated
using a set of data different from that used for calibration. The simulation and experimental results
show a fairly good agreement. The simulation results also shows the effectiveness of the fuel
vaporization model and the combustion model based on the proposed parameter-varying Wiebe
function that models the coupling combustion effect between the pre- and main-chambers. In
conclusion, the developed model is able to capture the TJI combustion characteristics in both
chambers, especially the interaction between them. This indicates that the developed model can
be used for developing the TJI engine control strategies that optimize both pre- and main-chamber
control variables. Future work is to further validate and calibrate the model over the entire metal
engine operational range.

54

CHAPTER 4
OPTIMAL COMBUSTION PHASING MODELING AND CONTROL OF A
TURBULENT JET IGNITION ENGINE

4.1

Introduction

In recent years, improving the thermal efficiency of internal combustion (IC) engines becomes
increasingly important for powertrain researchers and engineers [26]. According to the analysis
based on 2015 EPA certification data, a further 30% reduction in fuel consumption is required for
gasoline engines to meet the US 2025 GHG regulations, if the improvement is originated solely
from IC engines [27].
TJI combustion is a promising combustion technology that is able to achieve low engine out
emissions and high thermal efficiency. The application of the TJI combustion technology ranges
from engines for passage cars to commercial vehicles, such as delivering trucks, and stationary
power generator, with a thermal efficiency improvement of more than 25 percent[27]. Moreover,
the TJI technology was even used in Ferrari’s Formula One racing car [28]. The history of TJI
engines can be traced back to the 1940s [14]. The TJI system mainly consists of three parts: a
large main-chamber, a small pre-chamber (a few percent of main-chamber volume), and a singleor multi-orifice nozzle connecting the two combustion chambers. Once the pre-chamber is ignited,
the pressure difference between the two sides of the nozzle drives the hot reacting products into
the main-chamber. The hot reacting products from the pre-chamber are called turbulent jets. Since
the turbulent jets contain significantly large ignition energy, the main-chamber can be ignited under ultra lean conditions. There are two major advantages for lean combustion. The first is the
improvement of engine thermal efficiency. It is achieved by increasing the specific heat ratio of
the air-fuel mixture inside the main-chamber with extra air [14]. Lean combustion can also reduce
the NOx (nitrogen oxides) formation due to the low combustion temperature [14]. Therefore, lean
or diluted combustion attracts more researchers’ attention as the regulations on NOx emissions

55

become increasingly stringent [29].
Another factor that directly affects the engine thermal efficiency is the combustion phase [30].
Among several critical combustion phasing locations of interest, CA50 (crank location when 50%
of fuel is burned) is of primary interest, and thus, is widely used for combustion phasing analysis and control [31]. The ignition timing associated with the most efficient combustion is called
maximum break torque (MBT) timing that can be achieved by regulating the CA50 location to its
optimum. Once the CA50 is far away from that location, the efficiency decreases and the combustion also becomes less stable, leading to high coefficient of variation (COV) of the indicated mean
effective pressure (IMEP). Therefore, combustion phasing control is important for both conventional spark ignition (SI) and TJI engines.
Different methods have been used for combustion phasing control. Generally, these methods
can be divided into two groups. The first focuses on finding the optimal combustion phase. For
example, [32–34] use the extreme seeking control method to find the optimal spark timing and reference [35] uses an adaptive optimal control method similar to extreme seeking. The second group
mainly focuses on combustion phase detection and use it for optimal combustion phase control.
In reference [31], a combustion phasing model is developed based on an artificial neural network
and used for closed loop control, reference [36] uses an adaptive control method based on a radial
basis function network, and reference [37] utilizes the hypothesis test for combustion phasing control. Reference [38] developed an ionization based closed-loop combustion phase detection and
control strategy. In this research, the optimal CA50 values under different engine operational conditions are derived from the developed engine model and is regulated at its optimal value through
feedforward and feedback control. Different combustion phasing control methods have been developed for a various IC engines, but not for TJI engines. For TJI engines, besides spark timing,
pre-chamber control variables can also affect the combustion phasing in the main-chamber. This
makes the open-loop map-based controller less accurate, because the main-chamber combustion
phasing is subject to the characteristics of turbulent jets from the pre-chamber. Therefore, modelbased feedforward and feedback control are preferable. The increased number of control variables

56

Figure 4.1: TJI engine architecture

makes the control problem more complicated, but they provide the additional degree of freedoms
in control to improve the closed-loop system performance, which will be shown in the controller
validation section of this chapter.

4.2

Target Engine Description and Specifications

Fig. 4.1 shows the basic architecture of the single cylinder Dual-Mode TJI (DM-TJI) engine
studied in this paper. Each cylinder has two combustion chambers, pre and main ones. The prechamber is located on top of the main-chamber. Two chambers are connected through six small
orifices between them. The engine has two fuel injectors, one for pre-chamber (see Fig. 3.1)
and one for main-chamber in the intake port. A distinctive feature of this DM-TJI engine is the
inclusion of an air injector in the pre-chamber. The air injector is located at the opposite side of
the pre-chamber fuel injector. When air injector is turned on, the pressurized air is injected against
the fuel jets to reduce the amount of fuel deposited on the pre-chamber wall. More importantly,
the air injector is also used to control the air-to-fuel ratio (AFR) in the pre-chamber to be close to
stoichiometry and to improve the pre-chamber combustion characteristics, especially under heavy
exhaust-gas-recirculation. For TJI engines, the combustion is firstly initiated in the pre-chamber
by the spark plug, and then the hot reacting products in the pre-chamber are pushed into the mainchamber through the connecting orifices between two chambers and ignite the lean mixture in the

57

main-chamber. The engine specifications are listed in Table 3.1. Both chambers use EPA LEV-II
liquid gasoline.

4.3

Nonlinear State-Space TJI Engine Model

During engine operation, the combustion in the current engine cycle is influenced by the residual gas from the previous cycle. For the TJI engine discussed in this article, the fuel film left from
the previous cycle on the pre-chamber also influences the combustion during the current cycle.
These are the sources of the engine cycle-to-cycle dynamics. The states of the state-space engine
model thus need to describe the properties of the residual gas in both chambers and the fuel film on
the pre-chamber wall. Fig. 4.2 shows how the state-space engine model is developed. The whole
engine cycle is divided into several different stages for each state. The state-space equations can
be obtained by combining the governing equations of each stage. The details of the governing
equations for each stage is presented in the next subsections.

4.3.1

Main-chamber state-space model

The major assumptions for developing the main-chamber state-space model are:
1. Because the combustion in the main-chamber is always lean, it is assumed that there is no
unburned fuel left after the combustion.
2. The intake and exhaust manifold pressures, PInt and PExh are assumed to be equal to the
main and Pmain , respectively.
main-chamber pressures at IVC and EVC, PEVC
IVC

3. The mass of the gas inside the main-chamber is assumed to be constant when both the intake
and exhaust valves are closed.

58

𝑝𝑟𝑒

𝑚𝑓,𝑓𝑖𝑙𝑚 𝑘 + 1

𝑝𝑟𝑒

𝑝𝑟𝑒

𝑚𝑟

𝑚𝑟𝑚𝑎𝑖𝑛 𝑘

Pressure (bar)

𝑚𝑎𝑖𝑛
𝑚𝑟,𝐼𝑛𝑒
𝑘

50
40
30
20
10
0

1

𝑘

𝑚𝑎𝑖𝑛
𝑚𝐼𝑉𝐶
𝑚𝑎𝑖𝑛
𝑚𝑓,𝑓𝑟𝑒𝑠ℎ

9
𝑝𝑟𝑒

6

2

𝑃𝑆𝑃𝐾

𝑝𝑟𝑒

𝑇𝑆𝑃𝐾

𝑚𝑎𝑖𝑛
𝑃𝑆𝑃𝐾
𝑚𝑎𝑖𝑛
𝑇𝑆𝑃𝐾

𝑝𝑟𝑒

7

3

𝑝𝑟𝑒

𝑃𝐸𝑂𝐸

8

𝑝𝑟𝑒

𝑇𝐸𝑂𝐸
𝑚𝑎𝑖𝑛
𝑃𝐸𝑂𝐶
𝑚𝑎𝑖𝑛
𝑇𝐸𝑂𝐶

𝑚𝑎𝑖𝑛
𝑇𝐸𝑉𝑂

𝑘+1

𝑚𝑟𝑚𝑎𝑖𝑛 𝑘 + 1

𝑚𝑎𝑖𝑛
𝑃𝐸𝑉𝑂

4

𝑚𝑟

5

𝑚𝑎𝑖𝑛
𝑚𝑟,𝐼𝑛𝑒
𝑘+1

SPK
EVC

-300

IVC

-200

INJ

-100

EVO

0

100
200
300
Crank position (deg ACTDC)

IVO

EVC

400

Figure 4.2: State-space TJI engine model architecture.

59

IVC

500

INJ

600

700

20
15
10
5
0
800

Valve lift (mm)

𝑚𝑓,𝑓𝑖𝑙𝑚 𝑘

According to the above two assumptions, the properties of the residual gas in the main-chamber
can be fully described by the following two states:
1. The mass of inert gas left from last cycle at EVC: mmain
r,Ine (k). Here, k represents the engine
cycle index.
2. The mass of residual gas left from last cycle at EVC: mmain
(k)
r
From EVC to IVC, the volume occupied by residual gas changes a little due to the pressure
main to Pmain .
change from PEVC
IVC
main
Vrmain (k) = VEVC

main (k)
PEVC

!1/n
(4.1)

main (k)
PIVC

main is the main-chamber volumes at EVC and IVC, respectively; n is the polytropic index.
where VEVC

The volume of the fresh air-fuel mixture flowing into the main-chamber then can be calculated
accordingly:
main
main
V fmain
resh (k) = VIVC −VEVC

main (k)
PEVC

!1/n

main (k)
PIVC

(4.2)

The mass of the fresh air-fuel mixture and the average in-cylinder temperature can be calculated
by ideal gas law:
main (k)V main (k)
PIVC
f resh
main
m f resh (k) =
RTInt (k)

(4.3)

main
main
main (k) = PIVC (k)VIVC
TIVC
mmain
IVC (k)R

(4.4)

where R is the gas constant and mmain
IVC is the total mass inside the main-chamber at IVC:
main (k) + mmain (k)
mmain
IVC (k) = mr
f resh

(4.5)

The mass of the fresh fuel flowing into the main-chamber can be calculated based on intake
manifold AFR, ϕInt :

main
mmain
f , f resh (k) = m f resh (k)

60

1
ϕInt (k) + 1

(4.6)

During the compression, combustion and expansion processes, the mass of the gas inside the
main-chamber is a constant. The temperature during combustion can be simplified to be a combination of polytropic process and heat release process [39]. Therefore, the temperature at end of
combustion (EOC) can be represented by:

main (k) = T main (k)
TEOC
IVC

main
VIVC

!n−1
+

main (k)
VEOC

QLHV mmain
f , f resh (k)
mmain
IVC (k)Cv

(4.7)

main can be calculated based on ideal gas law.
And PEOC

The cank angle position at EOC, θEOC , can be obtained by first substituting xb (θEOC ) = 0.9
into the first equation in Eqn. (3.24):
 1

θEOC − θign_b
ln(0.1) m+1
= Sol
= −
∆θd
a

(4.8)

Then by using the second equation in Eqn. (3.24) and Eqn. (3.26), θEOC can be calculated by:
θEOC = θign + ∆θd · Sol −

Z θ
EOC

β ṁtur (θ )dθ

θign

(4.9)

= θign + ∆θd · Sol − β mtur
where θ0 can be approximated by spark timing, and the total mass of the turbulent jet during
combustion, mtur , will be calculated in the next subsection.
After combustion, the total mass of the inert gas can be calculated by:
main
main
mmain
EOC,Ine (k) = mr,Ine (k) + m f , f resh (k)(1 + AFRs)

(4.10)

where AFRs is the stoichiometric AFR. After the combustion, the gas in the main-chamber polytropically expands to EVO:
main (k) = T main (k)
TEVO
EOC

main (k)
VEOC

!n−1

main
VEVO

(4.11)

main can be calculated by the ideal gas law. After EVO, the gas in the main-chamber
And PEVO

expands until its pressure is equal to the exhaust manifold pressure.
main (k) =
TEVC

main (k)
PEVC
main (k)
PEVO

61

! n−1
n

main (k)
TEVO

(4.12)

main can be approximated by P
where PEVC
Exh .

main can be approximated by P
PIVO
Exh . Then the mass of residual gas in the main-chamber can

be calculated by:
mmain
IVO (k + 1) =

main
PExh (k + 1)VIVO
 main  n−1
n
P
(k)
main (k)]
R[ IVO
TEVO
main

(4.13)

PEVO (k)

From IVO to EVC, part of the gas in the exhaust manifold flows back into the engine cylinder due
to the pressure difference between intake and exhaust manifolds. This part of the residual gas is
calculated based on the widely used Fox model[40].
mmain
back f low (k + 1) =
PInt −0.87 p
C1 OF
(
)
|PInt − PExh |mmain
IVC (k + 1)
Ω P

(4.14)

Exh

Since mmain
IVC (k + 1) cannot be determined without knowing the mass of residual gas, its value
is approximated by
main (k + 1)V main
PIVC
IVC
mmain
(k
+
1)
=
IVC
RTInt (k + 1)

(4.15)

Then the mass of residual gas left for next cycle can be calculated by:
main (k + 1)V main
PEVC
EVC
mmain
(k
+
1)
=
r
main (k)
RTEVC

(4.16)

The content in the residual gas is the same as that at EOC. Therefore, the inert gas left for the
next cycle can be calculated by:
mmain
EOC,Ine (k)
main
main
mr,Ine (k + 1) = mr (k + 1)
mmain
IVC (k)

4.3.2

(4.17)

Pre-chamber state-space model

In the pre-chamber, the AFR of the air-fuel mixture is always kept close to or richer than stoichiometry. Therefore, there is usually no air left after combustion in the pre-chamber. There might
be some vaporized fuel left in the residual gas inside the pre-chamber, but its amount is very small

62

compared with the fuel amount from vaporization and injection. The mass of the fuel film on the
pre-chamber wall is also very important because it determines how much fuel will vaporize into
the gas mixture before ignition. As a result, the states selected for pre-chamber are:
pre

1. The mass of inert gas left from last cycle at IVC: mr (k).
pre

2. The mass of fuel film on the pre-chamber wall at IVC: m f , f ilm (k)
Before IVC, the gas flow through the orifice connecting the two chambers mainly goes from
the pre-chamber to the main-chamber. Therefore, the gas mixture in the pre-chamber is mainly the
residual gas left from last cycle. At IVC, it is assumed that:
pre

main = P
PIVC = PIVC
Int

(4.18)

pre

where PIVC is the pre-chamber pressure at IVC.
pre

During the compression stroke, the residual gas in pre-chamber is compressed from PIVC , V pre
pre

pre

pre

main . The mass of the residual gas
to PSPK , Vr,SPK . Here, it is assumed that PSPK is the same as PSPK
pre

mr

keep unchanged.

pre

Vr,SPK (k) =

PInt (k)
main (k)
PSPK

!1/n
V pre

(4.19)

where the subscript SPK denotes the properties at spark timing and V pre is the pre-chamber volpre

ume. Then Tr−SPK can be obtained from ideal gas law.
The volume and mass of the gas coming from the main-chamber at SPK can be then calculated
by the following equations:
pre

pre

pre

Vmain,SPK (k) = V pre −Vr,SPK (k) −Vin j (k)
pre

(4.20)
pre

where Vin j is the total volume of the injected air and vaporized fuel. And mmain,SPK can be then
calculated.
The average temperature in pre-chamber at SPK is:
pre

P (k)V pre
pre
TSPK (k) = SPK
pre
mSPK (k)R
63

(4.21)

pre

pre

pre

pre

mSPK (k) = min j (k) + mr (k) + mmain,SPK (k)

(4.22)

pre

where min j is the total mass of the injected air and vaporized fuel in pre-chamber.
The combustion in the pre-chamber is very fast, so a constant volume combustion is used to
approximate the combustion process:
pre

pre
pre
TEOC (k) = TSPK (k) +

m f uel,vap (k) ¡ QLHV

pre

pre
PEOC (k) =

(4.23)

pre

mSPK (k)Cv
pre

mSPK (k)RTEOC (k)
V pre

(4.24)

pre

where m f uel,vap is the total vaporized fuel in the pre-chamber.
Then the gas in the pre-chamber expands until its pressure is similar to the main-chamber
pressure. This time is defined as end of expansion (EOE).
 n−1
n
main
PEOE (k)
pre
pre

TEOC (k)
TEOE (k) =  pre
PEOC(k)


(4.25)

main can be approximated by the pressure at
If SPK is controlled around the optimum value, PEOE
pre

top dead center when motoring the engine. Then the mass at EOE, mEOE can be calculated by
ideal gas law. The mass flow rate of the turbulent jet after EOE is usually very small, so mtur in
Eqn. (4.9) can be calculated by:
pre

pre

mtur (k) = mSPK (k) − mEOE (k)

(4.26)

After EOE, it is assumed that all the gas pushed into pre-chamber does not mix with the gas in
pre-chamber and then was pushed out as main-chamber pressure decrease. Therefore, the gas in
the pre-chamber is first polytropically compressed and then expands to IVC.
pre
TIVC (k + 1) =

main (k + 1)
PIVC
pre
PEOC (k)

! n−1
n

pre

TEOC (k)

(4.27)

The mass of the residual gas can then be calculated based on the ideal gas law.
P (k + 1)V pre
pre
mr (k + 1) = Int pre
RTIVC (k + 1)
64

(4.28)

The fuel film on the pre-chamber wall keeps vaporizing during the entire engine cycle. The
vaporization rate is shown in Eqn. (3.14). The mass of the fuel film after pre-chamber fuel injection
can be expressed by:
pre

pre
ṁ f , f ilm (t) = −

m f , f ilm (t)

pre

δf ρf

∆H

−
Ps0
R T (t)
M f hm
e m f
Rm T f (t)

(4.29)

pre

Here, T f (t) is approximated by at (TEOC + TIVC ), where at is a calibration parameter. The initial
pre

value of m f , f ilm (t), that is the fuel film mass right after the fuel injection, is:
pre

pre

pre

m f , f ilm (0) = m f , f ilm (k) + min j,w (k)

(4.30)

pre

where m f , f ilm (k) is the fuel film mass right before the fuel injection in engine cycle k, and
pre

min j,w (k) is the mass of the injected fuel shooting on the wall.
Then, the mass of the fuel film right before fuel injection during the next engine cycle can be
calculated by solving Eqn. (4.29).
h
i
pre
pre
pre
m f , f ilm (k + 1) = m f , f ilm (k) + min j,w (k) e−avap ∆t(Ω)

(4.31)

where ∆t(Ω) is the time duration of each engine cycle, Ω is engine speed, and avap is:
M f hm,ss Ps0 − R∆H
T
e m f
avap =
δ f ρ f Rm T f

(4.32)

where hm,ss replaced hm in Eqn. (4.29) as a calibration parameter in the state-space model.

4.3.3

CA50 calculation

The parameter, ∆θd , in Eqn. (4.9) is mainly determined by the main-chamber gas properties. According to the experiment results, ∆θd can be approximated by a linear function of PInt , Ω, λ and
main , where Îť is the normalized air-fuel ratio in the main-chamber and xmain is the main-chamber
xCSP
CSP

CSP mass fraction at IVC. Then, θCA50 can be finally obtained by
θCA50 = θSPK + ∆θ0 + a(PInt − PInt0 ) + b(Ω − Ω0 )
main − xmain ) + e(m − m
+c(λ − λ0 ) + d(xCSP
tur
tur0 )
CSP0

65

(4.33)

main , and m
where the parameters, PInt0 , Ω0 , λ0 , xCSP0
tur0 , represent the nominal operational values.

A linear least-square fitting process is used to determine the unknown parameters, ∆θ0 , a,
b, c, d and e, in Eqn. 4.33. Here, the following vectors are defined: P = [∆θd0 , a, b, c, d, e]T
main −
is the parameter vector; y j = [θCA50 − θign ]; and x j = [1, PInt − PInt0 , Ω − Ω0 , λ − λ0 , xCSP
main , m − m
T
xCSP0
tur
tur0 ] . and the least-square method calculates the parameter vector that minimize

the following expression using N data points.
N

∑ (y j − x j ∗ P)2

(4.34)

j=1

where j is data point index. The solution of the problem is
P = (X T ∗ X)−1 X T ∗Y

(4.35)

where the rows of X and Y are xTj and y j respectively. Then the parameters, ∆θd0 , a, b, c, d and e,
can be obtained.
Considering all the equations above, we finally got the nonlinear state-space equation in the
following form.
x(k + 1) = f [x(k), u(k)]
(4.36)
y(k) = h[x(k), u(k)]
where the function f and h are shown in the Appendix A, x(k) is the state vector; u(k) is the control
input vector; and y(k) is the output vector defined below.
pre

pre

T
x(k) =[mr (k), m f , f ilm (k), mmain
(k), mmain
r
r,CSP (k)]
pre

pre

T
u(k) =[θSPK (k), min j, f uel (k), min j,air (k), mmain
f uel (k), PInt , Ω(k)]

(4.37)

y(k) =[θCA50 (k)].
where Ω is engine speed.

4.4

Controller Development

It can be found that the burn duration in the main-chamber is determined by the gas properties in
both combustion chambers. In other words, the CA50 can be controlled by spark timing and other
66

Figure 4.3: Control system architecture.

variables in the pre-chamber. In this paper, three variables are used to control the main-chamber
pre

CA50. They are spark timing, θSPK ; pre-chamber fuel and air injection amounts, min j, f uel and
pre

min j,air . The other parameters, mmain
f uel (k), PInt (k) and Ω(k), are considered as disturbances. Fig. 4.3
shows the basic architecture of the control system. The feedforward controller is a model-based
controller. The feedback controller uses output constraint control method to regulate the CA50 to
the desired value. The main-chamber AFR is controlled by the main-chamber fuel injector using a
map-based feedforward control as a function of intake manifold pressure.

4.4.1

Feedforward control

The feedforward controller calculates the value of pre-chamber air and fuel injection amounts
and spark timing based on mmain
f uel (k), PInt (k) and Ω(k). The pre-chamber air injection quantity is
determined by a pre-calibrated map. The pre-chamber fuel injection quantity and spark timing is
obtained by solving


Ν pre = 1
SPK

(4.38)


θ
CA50 = θCA50,desired
pre

where ÎťSPK is the normalized air-fuel ratio in the pre-chamber at spark timing. The detailed
expression of these equations can be easily obtained based on the discussion in Section III, as

67

follows

 main main main
mIVC −m f uel −mr,CSP pre

pre


mmain,SPK + min j,air

main

mIVC



=


ϕ
s




 mmain

f uel pre
pre
pre


mmain,SPK + (1 − k)min j, f uel + ṁ f ,vap ∆θSPK

main

mIVC






 + (m pre + m pre )[1 − exp(−avap ∆θ
)]
f , f ilm

SPK0

in j,w

(4.39)










desired = θ


θCA50
SPK + ∆θ0 + a(PInt − PInt0 ) + b(Ω − Ω0 )





pre

main − xmain ) − e{m

+ c(λ − λ0 ) + d(xCSP

tur0 − mSpk [1
CSP0



pre

pre


PEOE V pre 1 pre pre m f uel,vap QLHV − 1

−(
) n (mSpk TSpk +
) n ]}
R
Cv
where
pre

pre

mSPK = mr

pre

pre

+ min j + mmain,SPK

(4.40)

pre

Since Eqn. (4.38) is a highly nonlinear equation of min j, f uel and θSPK , it is difficult to solve
the equations analytically. Numerical methods can be used, but the computational effort would be
too large for real time engine control. By investigating Eqn. (4.38), it can be found that there are
only two highly nonlinear terms that are functions of the two unknown variables. The other terms
pre

can be determined once the feedforward control inputs and min j,air are given. The two terms are

pre
pre


mmain,SPK = f1 (PInt , min j,air , θSPK )
pre
(4.41)
m f uel,vap QLHV 1

pre
pre
pre
pre
pre

−
(m T +
) n = f2 (PInt , min j,air , mr , min j, f uel , θSPK )
Spk Spk
Cv
since the expression of f1 and f2 can be easily derived based on results in Section III, thus will not
pre

be shown here. Then, min j, f uel and θSPK can be solved for any desired θCA50 by replacing the two
terms with their linear approximation

∂ f1


f1 =
(θ
− θSPK0 ) + f1 (θSPK0 )



∂ θSPK SPK



∂ f2
∂ f2
pre
f2 =
(θSPK − θSPK0 ) +
(min j, f uel
pre

∂ θSPK
∂ min j, f uel






pre
 − m pre
in j, f uel0 ) + f 2 (min j, f uel0 , θSPK0 )
68

(4.42)

pre

The values of mr

main come from the state estimator in the feedback controller.
and mCSP

The pre-chamber AFR needs to be maintained around the stoichiometry for good ignitability
and stable combustion in the pre-chamber. For the same main-chamber condition, the AFR in the
pre-chamber is determined by all three control inputs, but only two are independent. If two of them
are given, the remaining one can be calculated according to the first equation in Eqn. (4.38). As a
pre

result, min j,air is not used by feedback controller but calculated based on the other two controller
pre

outputs, θSPK and min j, f uel , in the ’Air Compensation’ block of Fig. 4.3.

4.4.2

Output covariance constraint control

The nonlinear state-space engine model was first linearized by Matlab around the operational condition with IMEP=6 bar, λ =1.85, and Ω=1500 rpm. The linearized TJI engine model can be
expressed as
ex (k) + Be
ew
eu (k) + G
e (k)
xe(k + 1) = Ae

(4.43)

ex (k) + De
e u (k) + H
ew
e (k) + ve(k)
ye(k) = Ce
pre
T
e
e
e in j, f uel (k)]T is the control input; and w
e (k) = [m
e main
where ue(k) = [θeSPK (k) , m
f uel (k) , PInt (k) , Ω (k)]

is considered as disturbance. The parameter with the overhead tilde denotes the deviation of the
corresponding parameter from its reference value, for example, x (k) = x0 (k) + xe(k), given x0 (k)
the reference value around which the model is linearized.
For the current control hardware, the value of CA50 needs one engine cycle to compute and
transmit to the engine controller. Therefore, the measured CA50 is delayed by one engine cycle,
resulting one more state,
xd (k + 1) = ye(k) ,

(4.44)

Accordingly, the actual plant model is
x p (k + 1) = A p x p (k) + B p u (k) + G p w p (k)
y p (k) = C p x p (k) + v p (k)

69

(4.45)

where



04×1







e
e
 A

 B 
Ap = 
 , Bp = 
,
e
Ce
0
D




e 4×1
 G 0

Gp = 
 , C p = 01×4 1
e
H
1

(4.46)

and the augmented state vector and the disturbance input vector are




e (k) 
 xe(k) 
 w
x p (k) = 
 , w p (k) = 
.
xd (k)
ve(k)

(4.47)

The controller is designed so that the CA50 tracks the desired value r (k). The tracking error is
defined as
e (k) = r (k) − y p (k) .

(4.48)

To eliminate the steady-state tracking error, an integral action is introduced by defining the
integration of the tracking error as
ei (k + 1) = ei (k) + Ts e (k) .

(4.49)
T


For simplicity, Ts is set to 1. The augmented state vector xi (k) =

x p (k) , ei (k)

leads the

following augmented state equation
xi (k + 1) = Ai xi (k) + Bi ui (k) + dr(k) + Gi wi (k)
(4.50)

yi (k) = Ci xi (k)
zi (k) = Mi xi (k) + vi (k)
where


 Ap
Ai = 
−C p


05×1









05×1


 Bp 

 , Bi = 
, d = 
1
0
1



5×1
1×5 1
 G 0

 0

Gi = 
 , Ci = 

0
1
Cp 0





and the vector zi (k) is the noisy measurements of yi (k). Therefore, Mi = Ci .

70

(4.51)

Suppose that the plant is controlled by a strictly proper output feedback stabilizing control law
given by
xc (k + 1) = Ac xc (k) + Fzi (k)

(4.52)

ui (k) = Kc xc (k) + Kr r(k)
Then, the closed-loop system becomes
x (k + 1) = Ax (k) + Dr r(k) + Dw (k)
(4.53)

yi (k) = Cx (k)

ui (k) = Cu x (k) + Kr r(k)




 xi (k) 
 wi (k) 
where x (k) = 
 , w (k) = 
. The closed-loop system matrices A, C, D, and Dr
xc (k)
vi (k)
can be easily obtained based upon the above equations.
Consider the closed-loop system above. Let Wi and Vi denote positive symmetric matrices with


dimension compatible to the process noise wi and vi . Defined W = block diag

Wi Vi

and X

the closed-loop controllability Gramian from the weighted disturbance input W −1/2 w. Since A is
stable, X is the positive semi-definite solution of the following Lyapunov equation
X = AXAT + DW DT .
(4.54)


 y1 
Rewrite the performance output vector yi into yi := 
, where y j = C j x for j = 1, 2.
y2
Then, the OCC problem is to find a output feedback controller for the plant that minimizes the
OCC cost
JOCC = trace RCu XCuT , R > 0,

(4.55)

subject to the output covariance constraints
Y j = C j XCTj < Ȳ j , j = 1, 2,

(4.56)

where Ȳ j > 0 ( j = 1, 2) are given.
The OCC problem is actually a linear quadratic (LQ) control problem with a special choice of
output-weighting matrix Q. The value of Q is obtained by the iteration algorithm proposed in [41].
71

According to [42] and [41], the feedback controller is in the following form.
xc (k + 1) = (Ai + Bi Kc − FM p )xc (k) + Fzi (k)

(4.57)

ui (k) = Kc xc (k) + Kr r (k)
Note that
F = Ai KMiT (V + Mi KMiT )−1 ,

(4.58)

where K is the solution to the algebraic Riccati equation
K = Ai [K − KMiT (Vi + Mi KMiT )−1 Mi K]ATi + DiWi DTi ;

(4.59)

control gain Kc is provided by
Kc = −(R + BTi SBi )−1 BTi SAi ,

(4.60)

where S is the solution to the algebraic Riccati equation
S = ATi [S − SBi (BTi SBi + R)−1 BTi S]Ai +CiT QCi ;

(4.61)

Kr = (BTi SBi + R)−1 BTi F1 ,

(4.62)

F1 = −[I − ATi + ATi SF3 F2 ]−1 ATi SF3 F2 d

(4.63)

F2 = BTi R−1 Bi , F3 = (I + F2 S)−1 .

(4.64)

and

where

and

The OCC problem has several different interpretations. In a stochastic point of view, first
assume that wi and vi are uncorrelated zero-mean white noises with intensity matrices Wi > 0 and
Vi > 0. Let E be the expectation operator and define E∞ [·] := limk→∞ E [·], it is easy to see that
h
i


OCC problem is to minimize E∞ ui (k) RuTi (k) subject to the OCCs Y j := E∞ y j (k) yTj (k) <
Ȳ j , j = 1, 2. These constraints may be interpreted as constraints on the variance of the performance
variables or lower bounds on the residence time of the performance variables [43].
72

In a deterministic point of view, first define `∞ and `2 norms as follows

 
2

y j 
 := supk≥0 yT (k) y j (k)
j
∞
kwk 22

T
:= ∑∞
k=0 w

,

(4.65)

(k) w (k)

and the `2 disturbance set





 −1/2 
 2
n
w
W := w ∈ R s.t. 
W
w
 ≤ 1 ,
2

(4.66)

where W > 0 is a real symmetric matrix. Then, for any w ∈ W , the following inequalities hold

 
2
 

y j 
 ≤ Y j , j = 1, 2,
∞

 
2 


u j 
 ≤ Cu XCuT , j = 1, 2.
jj
∞

(4.67)

where [¡] j j is the jth diagonal entry. Moreover, [44] and [45] show that the above bounds are the
least upper bounds that hold for any arbitrary signal w ∈ W .
In other word, if the control weighting matrix is defined as R := diag [r1 , r2 ...rnu ], where nu
is the dimension of u (in this paper, nu = 2), the OCC problem is to minimize the sum of the
worst-case `∞ norms on the control signals given by
(
nu

JOCC =

∑ rj

j=1

)


 
2
sup 
u j 
 ∞

(4.68)

w∈W

subject to the constraints on the worst-case `∞ norms on the performance variables of the form

 
2
sup 
y j 
 ∞ ≤ Ȳ j , j = 1, 2.

(4.69)

w∈W

Note that the corresponding cost function of the LQ controller is defined as
∞

J=

∑ [yi(k)T Qyi(k) + ui(k)T Rui(k)]

(4.70)

k=0

4.4.3

Baseline controllers

To verify the performance of designed controllers, several different baseline controllers were also
designed. The first kind of baseline controller is proportional and integral (PI) controller.
GPI (z) = KP +
73

KI
z−1

(4.71)

Table 4.1: PI controller parameters
P Gain

I Gain

Ziegler-Nichols

0.427

0.128

Modified Ziegler-Nichols

0.188

0.090

Tyreus-Luyben

0.293

0.033

The proportional and integral gains ("P" and "I") of the PI controller are determined by three
different methods: Ziegler-Nichols (ZN), Modified Ziegler-Nichols (MZN) and Tyreus-Luyben
(TL), see references [46] and [47]. Their values are shown in Table. 4.1. PI controllers use only
spark timing, θSpk , to control the CA50. An OCC baseline controller is also developed only using
θSpk as the controller output. The pre-chamber fuel and air injection amounts are determined by
the feedforward controller.

4.5
4.5.1

Experimental Validation
Model validation

Experiments were conducted using the single cylinder DM-TJI engine described in Section II.
A TMAP (temperature and manifold air pressure) sensor is mounted in the intake manifold to
measure the pressure and temperature. Two Kistler pressure sensors are used to measure the mainchamber and exhaust manifold pressures. Another Kistler pressure sensor integrated with the spark
plug is used to measure the pre-chamber pressure; see Fig. 3.1. The exhaust Îť is obtained by an
ECM lambda meter, and the CA50 is calculated by the combustion analysis system from A&D
technology in real time. The experimental data used to calibrate and validate the model are shown
in Fig. 4.4. The engine speed ranges from 1200 to 2000 rpm; and IMEP (indicated mean value
pressure) ranges from 4.5 to 7.2 bar. The pre-chamber fuel injection amount varies from 0.7 to
1.35 mg.
Direct measurement of the four states in the state-space model is not possible on the TJI engine.
In the authors’ previous research, [48] and [49], a crank-angle-resolved (CAR) TJI engine model
was developed. It is able to calculate all the relevant parameters for every crank angle degree.

74

IMEP (bar)

7
6
5
4
1000

1200

1400
1600
1800
Engine speed (rpm)

2000

2200

Figure 4.4: Experimental matrix.
mmain
(mg)
r

mmain
(mg)
r,CSP

60

State−space model

State−space model

35

50
40
30
30

40

50
CAR model
mpre
(mg)
r

20
15
15

20
25
CAR model
mpre
(mg)
f,film

30

35

0.1
0.15
CAR model

0.2

0.2
State−space model

State−space model

25

10
10

60

3
2.5
2
1.5
1

30

1

1.5

2
2.5
CAR model

3

0.15
0.1
0.05
0

0

0.05

Figure 4.5: State-space engine model validation.

Since the CAR engine model was calibrated and validated by the experimental data under different
steady-state operational conditions, the state-space model can be validated by the simulation results
of the CAR engine model. The simulation results are shown in Fig. 4.5.
Fig. 4.5 compares the values of the four states calculated from the state-space engine model
and those from the CAR engine model at different operational conditions. The two dash-dot lines

75

on each side of the data points show the 5% error bound. All of the four states calculated by the
state-space model match the CAR model very well.
The calculated CA50 values are compared with the experiment results in Fig. 4.6. This figure also shows how CA50 is influenced by different variables. The upper-left plot compares the
experimental and calculated CA50 under five different operational conditions. These conditions
main , and
differ from each other mainly by their main-chamber λ . The other parameters, PInt , Ω, xCSP

mtur in Eqn. (4.33) remain unchanged or within small variations. The upper-right and bottom-left
plots are obtained similarly by mainly changing PInt and mtur , respectively. However, varying one
parameter while maintaining other parameters unchanged is challenging during the experiment.
For example, all the parameters, except Ω, will be changed when changing PInt . This is why the
simulation points cannot form a straight line in the upper-right plot. However, the trend of the
curve is evident: increasing PInt leads to reduced burn duration. The last plot compares CA50
values for all the experimental conditions. The dash-dot line shows the 10% error bound.
Fig. 4.7 shows the CA50 values when the pre-chamber fuel injection amount changes from 1.9
mg to 1.1 mg. The experimental CA50 data have a lot of noise, so it is hard to see the transient
response of CA50 from a single set of experimental data. Therefore, the experimental data in
Fig. 4.7 are obtained by averaging the results of 14 identical experiments. It can be seen that
the CA50 does not reach its steady-state value immediately after reducing the pre-chamber fuel
injection. This is caused by the fuel vaporization from the fuel film on the pre-chamber wall.

4.5.2

Controller validation

In the experiment, the controller is implemented into the MotoTron controller shown in Fig. 4.8.
The pressure and other signals from the engine are captured by the Phoenix AM module from A&D
Technology. The Phoenix RT module calculates the CA50 in real-time and sends it to MotoTron
controller through the CAN bus.
Fig. 4.9 compares the experimental results when the controller is turned on and off. When the
OCC controller is turned on, CA50 is adjusted close to the optimal value. The standard deviation

76

40

35

Spk−CA50 (deg)

Spk−CA50 (deg)

40

30
25
Simulation
Experiment

20
15

1.4

1.6
Lambda

1.8

35

30

25
0.7

2

40

0.8
0.9
PInt (bar)
Spk−CA50 (deg)

1

40
Simulation
Experiment

35

35

Simulation

Spk−CA50 (deg)

Simulation
Experiment

30

30
25
20

25
10

11

12
mtur (mg)

13

14

20

25
30
35
Experiment

40

CA50 (deg ATDC)

Figure 4.6: Comparison between experimental and calculated CA50 under different operational
conditions.

12
10
8
6

Experiment
Simulation

4
0

5

10

15

20

25
30
Engine cycle

35

40

45

50

Figure 4.7: Comparison between experimental and calculated CA50 when pre-chamber fuel is
changed.

77

Figure 4.8: Test bench structure.

→ Control on

15

→ Control off

Std=1.12

10

IMEP (bar)

CA50
(deg ATDC)

20

Std=1.58

5
7
6

Std=0.058
→ Control on

5
0

50

100

Std=0.068
→ Control off

150
200
Engine cycle

250

300

350

Figure 4.9: Experimental results with controller on and off.

(S.D.) of CA50 is reduced from 1.58 to 1.12 deg after the controller is turned on. The S.D. of
IMEP is also reduced accordingly. This indicates that the combustion phasing controller is able to
improve the combustion stability by regulating the CA50 close to its optimal value.
The OCC controller is compared with the PI baseline controllers in Fig. 4.10. The settling
times for the MZN and TL PI controllers are very long compared with the OCC controller. Since
the ZN controller has higher PI gains than these of MZN and TL controllers, the settling time is
improved. However, The high gain results in large steady state CA50 variations. This is because
the PI controller is not model-based and does not consider the engine cycle-to-cycle dynamics,
disturbances and measurement noises. However, OCC controller utilizes all the plant information,
and thus, has better performance.

78

CA50
(deg ATDC)

20

Measurement
Reference

PI TL
S.D.=1.91

10

CA50
(deg ATDC)

0
20
PI MZN
S.D.=2.05

10

CA50
(deg ATDC)

0
20
PI ZN
S.D.=2.31

10

CA50
(deg ATDC)

0
20
OCC
S.D.=1.34

10
0

0

50

100

150

200
250
Engine cycle

300

350

400

450

Figure 4.10: Experimental results of OCC and PI controllers.

In Fig. 4.11, the OCC controller is compared with the OCC baseline controller, which does not
pre

use pre-chamber fuel injection amount, min j, f uel , as the control input. The CA50 takes a longer
time to track the desired value. Although it takes a similar number of engine cycles to track the
desired CA50 for both cases, the combustion is less stable when only using spark timing as the
control input, which can be observed from the standard deviation of CA50 in Fig. 4.11. The reason
pre

is quite obvious. Since min j, f uel is not used for closed loop control, the controller relies only on
adjusting spark timing to track the desired CA50. This is why the variation of spark timing of the
baseline OCC controller is much larger than the proposed OCC controller; see the bottom plot in
Fig. 4.11. To achieve similar settling time, the baseline OCC controller needs a larger controller
gain than the proposed OCC controller, resulting in significant CA50 variations. It should also
be noted that the OCC baseline controller has better performance than other PI controllers. This
further demonstrates the benefit of using the model-based OCC control scheme.

79

CA50
(deg ATDC)

20

Measurement
Reference

OCC baseline
S.D.=1.80

10

CA50
(deg ATDC)

0
20
OCC
S.D.=1.34

10

Spark time
(deg ATDC)

mpre
(mg)
fuel

0
OCC

1.2
1
0.8
35

OCC
OCC baseline

30
25
20
0

50

100

150

200
250
Engine cycle

300

350

400

450

Figure 4.11: Experimental results of OCC controllers with and without using pre-chamber fuel.

Finally, Fig. 4.12 shows the experiment results under different operational conditions (engine
speed: 1200 and 2000 rpm; IMEP: 6.5 and 4.55 bar). Overall, the OCC controller is able to track
the desired CA50 within 10 engine cycles. Fig. 4.13 compared the standard deviations of the CA50
with different controllers discussed in this paper. The proposed and baseline OCC controllers have
overall lower S.D. values than the PID controllers. The proposed OCC controller has much lower
value than the baseline OCC controller because of utilization of pre-chamber fuel as a control
input.

4.6

Conclusions

In this chapter, a cycle-by-cycle state-space TJI engine model is first developed and validated
using both experimental data and high fidelity model. The validation results show that the CA50
(crank location when 50% fuel is burned) modeling error is less than 10%. An OCC (output co-

80

CA50
(deg ATDC)

20

Measurement
Reference

OCC 1200rpm
10

S.D.=1.35

CA50
(deg ATDC)

0
20
OCC 2000rpm
S.D.=1.48

10
0

0

50

100

150

200
250
Engine cycle

300

350

400

450

Figure 4.12: Experimental results of OCC controllers under different conditions.

Figure 4.13: Standard deviation of CA50 with different controllers.

variance constraint) combustion phasing controller is developed based on the statespace engine
model. The controller uses both spark timing and pre-chamber fuel quantity as the control inputs. Experimental results show that the proposed OCC controller has better performance than the
baseline controllers with a minimal improvement of over 22% for CA50 standard deviation. This
demonstrates that utilizing multiple pre-chamber control variables and OCC control scheme is able
to improve the combustion phasing control performance. The future work is to study if further performance improvement can be achieved by utilizing linear-parameter-varying (LPV) control based
on an LPV model that captures the engine dynamics accurately.

81

CHAPTER 5
CONCLUSIONS AND FUTURE WORK

5.1

Conclusions

In this dissertation, a control-oriented TJI combustion model is firstly developed for a rapid
compression machine (RCM) configured for TJI combustion, and then, the RCM model was extended to a TJI combustion model for the Dual-Mode (DM) TJI engine at Michigan State University. To achieve stable combustion with high thermal efficiency and low engine-out emissions, an
optimal combustion phasing controller is developed for the DM-TJI engine. The conclusions are
summarized as follows.
• A newly proposed parameter-varying Wiebe combustion model is used to link the combustion processes in both pre- and main-combustion chambers. The developed model can be
calibrated using a simple and systematic calibration procedure based on the experimental
data. The model validation process shows a good agreement between the model and experimental pressure traces, which indicates that the developed model is capable of accurately
capturing the TJI combustion dynamics. The validation results also indicate that the model
is able to predict the combustion process that is not used to calibrate the model parameters.
This shows that the developed model has the potential to be used for studying TJI combustion
engines and developing the associated control strategies.
• The TJI engine model was calibrated using experimental data at engine speeds between 1200
and 2000 rpm with indicated mean effective pressure between 4.2 to 7.2 bar, and then, the
calibrated model is validated using a set of data different from that used for calibration. The
simulation and experimental results show a fairly good agreement. The simulation results
also shows the effectiveness of the fuel vaporization model and the combustion model based
on the proposed parameter-varying Wiebe function that models the coupling combustion

82

effect between the pre- and main-chambers. The relative modeling error of the simulated
pressure curves is less than 8%. For most of the other pressure related variables, like indicated mean effective pressure and main-chamber burn duration, the relative errors are within
5%.
• The cycle-by-cycle state-space TJI engine model is developed and validated using both experimental data and high fidelity model. The validation results show that the CA50 (crank
location when 50% fuel is burned) modeling error is less than 10%. An OCC (output covariance constraint) combustion phasing controller is developed based on the state-space engine
model. The controller uses both spark timing and pre-chamber fuel quantity as the control
inputs. Experimental results show that the proposed OCC controller has better performance
than the baseline controllers with a minimal improvement of over 22% for CA50 standard
deviation. This demonstrates that utilizing multiple pre-chamber control variables and OCC
control scheme is able to improve the combustion phasing control performance.

5.2

Recommendations for Future Work

Currently, the parameter-varying Wiebe function only uses the mass flow rate of turbulent jets
to regulate the main-chamber burn rate. The modeling error increases as operational conditions
become far away from the nominal condition. Other properties of the turbulent jet should be
studied and used for main-chamber burn rate calculation. Possible properties include: turbulence,
temperature, kinetic/internal energy, etc.
In this research, the engine model used for control design is linearized around single operational
condition. Although the simulation results under different operational conditions showed similar
performance, further study need to be conducted to see if further performance improvement can be
achieved by linear-parameter-varying (LPV) control based on LPV model that captures the engine
dynamics accurately.
For the current TJI engine, pre-chamber combustion is not optimized through feedback control.
In the future, the pre-chamber pressure can be studied and used as the feedback signal. Ionization
83

signal is another good candidate. It is beneficial to study whether ionization signal could give more
information about the pre-chamber combustion and thus used for feedback combustion optimation.
In the future, engine subsystems, such as variable valve timing actuator and supercharger, will
be added and the control system will evolve correspondingly. Especially, main-chamber AFR/EGR
should be controled so that maximum dilution rate should be achieved while maintaining a stable
combustion.

84

APPENDICES

85

APPENDIX A
NONLINEAR STATE-SPACE MODEL

1

1

1 main
−1
−1
main
main main
mmain
(k + 1) = PEVC
(k + 1)VEVC
PEVC (k + 1) n R− n VEOC
(k) n mmain
r
IVC (k)

1

1

1
−1
main
main main
n −1 − n main
mmain
R VEOC (k) n
r,Ine (k + 1) = PEVC (k + 1)VEVC PEVC (k + 1)

"

main (k)V main
PIVC
IVC
R



"

main (k)V main
PIVC
IVC
R

main
VIVC
main (k)
VEOC

n−1
+



main
VIVC
main
VEOC (k)

n−1
+

QLHV mmain
f uel (k)

#− 1n

Cv

QLHV mmain
f uel (k)

#− 1n
(A.1)

Cv

h
i
main
mmain
r,Ine (k) + m f uel (k) (1 + AFRs)

(A.2)



− 1 
pre


main (k)
1
n
1
mIVC
(k)
VSPK
R
QLHV R pre
pre
pre
pre
pre
pre
pre
pre
pre
pre
main
min
(k)
+
m
(k)
+
V
−
V
−
mrpre (k + 1) = PIVC
(k)V pre +
m f uel,vap (k)
m
(k)
T
(k)
+
m
(k)
T
(k)
(k + 1) n V pre n PSPK
pre
r
j
in j,air
in j,air
in j, f uel
in j, f uel
main (k)
main
Cv
PSPK
(k)
VSPK
VIVC
(A.3)
"

#
h
i
M f hm,ss Ps0
∆H
pre
pre
m pre
(k
+
1)
=
m
(k
+
1)
+
m
(k)
exp∆t
−
exp
f , f ilm
f , f ilm
in j,wall
δ f ρ f Rm T f
Rm T f

where:
main
mmain
(k) +
IVC (k) = mr

"
main
VEOC
(k) = f θ SPK + 0.5 ∗ ∆θ 0

main
main
mmain
IVC (k) − m f uel (k) − mr,Ine (k)

Îť=

mmain
f uel (k)

(
mtur =

pre
1 − PEOE

1
n

(k) V

pre 1n

=



PInt
PInt0

a 

rpm
rpm0

b 

Îť
Îť0

c 



 main
 1n
main (k)
PIVC
main
main PEVC (k)
VIVC

−VEVC
main (k)
main (k)
RTInt
PIVC

mIne
mIne0

main
main (k) − mmain (k)
mmain
r,Ine
f resh (k) − m f uel (k) + mr

mmain
f uel (k)

(A.4)

#

d

"

− β mtur = f θ SPK + 0.5 ∗ ∆θ 0

−1
= mmain
f uel (k)



PInt
PInt0

(A.5)

a 

rpm
rpm0

b 

Îť
Îť0

c 

mIne
mIne0

#

d
− β mtur





 main
1
 Pmain (k)

PEVC (k) n
main
main
main
main
main
IVC
VIVC −VEVC
 − m f uel (k) + mr (k) − mr,Ine (k)
main (k)
main (k)
 RTInt

PIVC

(A.6)

(A.7)


− 1 ) 

pre


n
mIVC
(k)
V main (k)
QLHV R pre
R
pre
pre
pre
pre
pre
pre
pre
PSPK (k)V +
m f uel,vap (k)
min
V pre − SPKmain V pre − pre
min
(k) Tinpre
(k) + min
(A.8)
j (k) + mr (k) + main
j,air (k) Tin j,air (k)
j,
f
uel
j,
f
uel
Cv
PSPK (k)
VSPK (k)
VIVC

main n
VIVC
pre
main
main
PSPK
(k) = PSPK
(k) = PIVC
(k)
main
VSPK (k)

86

(A.9)

1

n
1
QLHV R pre
pre
m f uel,vap (k)
V pre1− n PSPK
(k)V pre +
Cv


pre

−1
mIVC
(k)
V main (k)
R
pre
pre
pre
pre
pre
pre
min
V pre − SPKmain V pre − pre
min
(k)
T
(k)
+
m
(k)
T
(k)
j (k) + mr (k) + main
in
j,air
in
j,air
j, f uel
in j, f uel
PSPK (k)
VSPK (k)
VIVC
1− 1n

pre
main
TIVC
(k + 1) =R−1 PIVC
(k + 1)

pre

pre
TEOC
(k) =

pre
PSPK
(k)V pre QLHV m f uel,vap (k)
+
R
Cv

!

pre
pre
min
j (k) + mr (k) +


pre
−1

main (k)
mIVC
(k)
VSPK
R
pre
pre
pre
pre
pre
pre
(k)
+
m
(k)
T
(k)
V
−
V
−
(k)
T
m
pre
in j,air
in j,air
in j, f uel
in j, f uel
main (k)
main
PSPK
(k)
VSPK
VIVC

87

(A.10)

(A.11)

APPENDIX B
ENGINE MODEL PARAMETERS

Table B.1: Parameters for Combustion Phase Model
Symbol

Eqn.

Value

Symbol

Eqn.

Value

C1

(4.14)

0.85

b

(4.33)

0.00105 s

OF

(4.14)

1.20 deg/m

c

(4.33)

17 deg

δf

(3.14)

0.03 mm

d

(4.33)

-660 deg

a

(4.33)

-60.9 deg/bar

e

(4.33)

-3.06 deg/mg

88

BIBLIOGRAPHY

89

BIBLIOGRAPHY

[1] Lino Guzzella and Christopher Onder. Introduction to modeling and control of internal combustion engine systems. Springer Science & Business Media, 2009.
[2] Benjamin Ellies, Charles Schenk, and Paul Dekraker. Benchmarking and hardware-in-theloop operation of a 2014 mazda skyactiv 2.0 l 13: 1 compression ratio engine. Technical
report, SAE Technical Paper, 2016.
[3] Reese R. Progress (and challenges) along the path to 2025. In Engine Research Center 2015
Symposium. University of Wisconsin, Madison, 2015.
[4] Timothy Johnson. Vehicular emissions in review. SAE International Journal of Engines, 9
(2016-01-0919):1258–1275, 2016.
[5] Neil Fraser, Hugh Blaxill, Grant Lumsden, and Michael Bassett. Challenges for increased
efficiency through gasoline engine downsizing. SAE International Journal of Engines, 2
(2009-01-1053):991–1008, 2009.
[6] Shupeng Zhang, Guoming Zhu, and Zongxuan Sun. A control-oriented charge mixing and
two-zone hcci combustion model. IEEE transactions on vehicular technology, 63(3):1079–
1090, 2014.
[7] Xiaojian Yang and Guoming Zhu. SI and HCCI combustion mode transition control of an
HCCI capable SI engine. Control Systems Technology, IEEE Transactions on, 21(5):1558–
1569, Sept 2013. ISSN 1063-6536. doi: 10.1109/TCST.2012.2201719.
[8] Rudolf H Stanglmaier and Charles E Roberts. Homogeneous Charge Compression Ignition
(HCCI): Benefits, Compromises, and Future Engine Applications. Technical report, SAE
Technical Paper, 1999.
[9] LA Gussak, Michael C Turkish, and Donald C Siegla. High chemical activity of incomplete
combustion products and a method of prechamber torch ignition for avalanche activation of
combustion in internal combustion engines. Technical report, SAE Technical Paper, 1975.
[10] Elisa Toulson, Andrew Huisjen, Xuefei Chen, Cody Squibb, Guoming Zhu, Harold Schock,
and William P Attard. Visualization of Propane and Natural Gas Spark Ignition and Turbulent
Jet Ignition Combustion. Technical report, SAE Technical Paper, 2012.
[11] William P Attard, Neil Fraser, Patrick Parsons, and Elisa Toulson. A turbulent jet ignition
pre-chamber combustion system for large fuel economy improvements in a modern vehicle
powertrain. SAE International Journal of Engines, 3(2010-01-1457):20–37, 2010.
[12] Michael C Turkish. 3-Valve Stratified Charge Engines: Evolvement, Analysis and Progression. Technical report, SAE Technical Paper, 1974.

90

[13] Tasuku Date, Shizuo Yagi, Akira Ishizuya, and Isao Fujii. Research and Development of the
Honda CVCC engine. Technical report, SAE Technical Paper, 1974.
[14] Elisa Toulson, Harold J Schock, and William P Attard. A review of pre-chamber initiated jet
ignition combustion systems. Technical report, SAE Technical Paper, 2010.
[15] Ruitao Song, Gerald Gentz, Guoming Zhu, Elisa Toulson, and Harold Schock. A controloriented jet ignition combustion model for an si engine. In ASME 2015 Dynamic Systems and
Control Conference, pages V001T11A001–V001T11A001. American Society of Mechanical
Engineers, 2015.
[16] Elisa Toulson, Harry C Watson, and William P Attard. Modeling alternative prechamber fuels
in jet assisted ignition of gasoline and lpg. Technical report, SAE Technical Paper, 2009.
[17] A Ghorbani, G Steinhilber, D Markus, and U Maas. Numerical investigation of ignition in a
transient turbulent jet by means of a pdf method. Combustion Science and Technology, 186
(10-11):1582–1596, 2014.
[18] M. Canova, R. Garcin, S. Midlam-Mohler, Y. Guezennec, and G. Rizzoni. A controloriented model of combustion process in a hcci diesel engine. In American Control
Conference, 2005. Proceedings of the 2005, pages 4446–4451 vol. 7, June 2005. doi:
10.1109/ACC.2005.1470696.
[19] Michael E Mortenson. Mathematics for computer graphics applications. Industrial Press
Inc., 1999.
[20] John B Heywood. Internal combustion engine fundamentals, volume 930. Mcgraw-hill New
York, 1988.
[21] JI Ghojel. Review of the development and applications of the wiebe function: a tribute to the
contribution of ivan wiebe to engine research. International Journal of Engine Research, 11
(4):297–312, 2010.
[22] Hua Zhao and Nicos Ladommatos. Engine combustion instrumentation and diagnostics. Warrendale, PA: Society of Automotive Engineers, 2001. 842, 2001.
[23] Erik HellstrĂśm, Anna Stefanopoulou, and Li Jiang. A linear least-squares algorithm for
double-wiebe functions applied to spark-assisted compression ignition. Journal of Engineering for Gas Turbines and Power, 136(9):091514, 2014.
[24] Asghar Ghorbani, Gerd Steinhilber, Detlef Markus, and Ulrich Maas. Ignition by transient
hot turbulent jets: An investigation of ignition mechanisms by means of a pdf/redim method.
Proceedings of the Combustion Institute, 35(2):2191–2198, 2015.
[25] Stephen R Turns et al. An introduction to combustion, volume 287. McGraw-hill New York,
1996.
[26] Yanjun Huang, Amir Khajepour, Tianjun Zhu, and Haitao Ding. A supervisory energy-saving
controller for a novel anti-idling system of service vehicles. IEEE/ASME Transactions on
Mechatronics, 22(2):1037–1046, 2017.
91

[27] Ravi Teja Vedula, Ruitao Song, Thomas Stuecken, Guoming G Zhu, and Harold
Schock. Thermal efficiency of a dual-mode turbulent jet ignition engine under lean
and near-stoichiometric operation.
International Journal of Engine Research, page
1468087417699979, 2017.
[28] Nicolas Carpentiers.
Turbulent Jet Ignition: In the antechamber of F1 power.
http://en.f1i.com/magazine/76167-tji-antechamber-f1-power.html, 2016.
[29] Hui Zhang and Junmin Wang. Adaptive sliding-mode observer design for a selective catalytic
reduction system of ground-vehicle diesel engines. IEEE/ASME Transactions on Mechatronics, 21(4):2027–2038, 2016.
[30] Enrico Corti and Claudio Forte. Spark advance real-time optimization based on combustion
analysis. Journal of Engineering for Gas Turbines and Power, 133(9):092804, 2011.
[31] Baitao Xiao. Adaptive model based combustion phasing control for multi fuel spark ignition
engines. Proquest Dissertations and Theses, 2013.
[32] Peter G Scotson and Peter E Wellstead. Self-tuning optimization of spark ignition automotive
engines. IEEE Control Systems Magazine, 10(3):94–101, 1990.
[33] Ibrahim Haskara, Guoming G Zhu, and Jim Winkelman. Multivariable egr/spark timing
control for ic engines via extremum seeking. In American Control Conference, 2006, pages
6–pp. IEEE, 2006.
[34] Enrico Corti, Claudio Forte, Giorgio Mancini, and Davide Moro. Automatic combustion
phase calibration with extremum seeking approach. Journal of Engineering for Gas Turbines
and Power, 136(9):091402, 2014.
[35] Qilun Zhu, Robert Prucka, Shu Wang, Michael Prucka, and Hussein Dourra. Model-based
optimal combustion phasing control strategy for spark ignition engines. SAE International
Journal of Engines, 9(2016-01-0818):1170–1179, 2016.
[36] Maru Yoon, Kangyoon Lee, and Myoungho Sunwoo. A method for combustion phasing
control using cylinder pressure measurement in a crdi diesel engine. Mechatronics, 17(9):
469–479, 2007.
[37] Jinwu Gao, Yuhu Wu, and Tielong Shen. Combustion phase control of si gasoline engines
using hypothesis test. IFAC-PapersOnLine, 48(15):153–158, 2015.
[38] Guoming G Zhu, Ibrahim Haskara, and Jim Winkelman. Closed-loop ignition timing control for si engines using ionization current feedback. IEEE transactions on control systems
technology, 15(3):416–427, 2007.
[39] Shupeng Zhang, Guoming Zhu, and Zongxuan Sun. A control-oriented charge mixing and
two-zone hcci combustion model. IEEE transactions on vehicular technology, 63(3):1079–
1090, 2014.
[40] Jonathan W Fox, Wai K Cheng, and John B Heywood. A model for predicting residual gas
fraction in spark-ignition engines. Technical report, SAE Technical Paper, 1993.
92

[41] Guoming Zhu, MA Rotea, and R Skelton. A convergent algorithm for the output covariance
constraint control problem. SIAM Journal on Control and Optimization, 35(1):341–361,
1997.
[42] Xuefei Chen, Yueyun Wang, Ibrahim Haskara, and Guoming Zhu. Optimal air-to-fuel ratio
tracking control with adaptive biofuel content estimation for lnt regeneration. IEEE Transactions on Control Systems Technology, 22(2):428–439, 2014.
[43] SM Meerkov and T Runolfsson. Output residence time control. IEEE Transactions on Automatic Control, 34(11):1171–1176, 1989.
[44] G Zhu, M Corless, and R Skelton. Robustness properties of covariance controllers. In Proceedings of Allerton Conference, Monticello, IL, 1989.
[45] G Zhu and R Skelton. Mixed l2 and l∞ problems by weight selection in quadratic optimal
control. International Journal of Control, 53(5):1161–1176, 1991.
[46] Michael L Luyben and William L Luyben. Essentials of process control. McGraw-Hill„
1997.
[47] Chang C Hang, Karl Johan Åström, and Weng Khuen Ho. Refinements of the ziegler nichols
tuning formula. In IEE Proceedings D Control Theory and Applications, 1991.
[48] Ruitao Song, Gerald Gentz, Guoming Zhu, Elisa Toulson, and Harald Schock. A controloriented model of turbulent jet ignition combustion in a rapid compression machine. Proceedings of the Institution of Mechanical Engineers, Part D: Journal of Automobile Engineering,
page 0954407016670303, 2016.
[49] Ruitao Song and Guoming G Zhu. Control-oriented crank-resolved and state-space models
for a gasoline turbulent jet ignition engine. In American Control Conference (ACC), 2017.
IEEE, 2017.

93