Adaptive Pressure Control for Use in Variable-Thrust Rocket
Development
Anil Alan and Yildiray Yildiz
Bilkent University, Ankara, Turkey

arXiv:1705.01849v1 [] 4 May 2017

Umit Poyraz
Roketsan INC., Ankara, Turkey

Abstract
The precise control of gas generator pressure is important to obtain variable thrust in ducted
rockets. A delay resistant closed loop reference model adaptive control is proposed in this paper
to address this problem. The proposed controller combines delay compensation and adaptation
with improved transient response. The controller is successfully implemented using an industrial
grade cold air test setup, which is a milestone towards obtaining a fully developed throttleable
rocket gas generator controller. Simulation and experimental comparisons with alternative adaptive approaches and a fixed controller demonstrate improved performance and effective handling of
time delays and uncertainties. A step by step design methodology, covering robustfying schemes,
selection of adaptation rates and initial controller parameters, is also provided to facilitate implementations.

1

I.

INTRODUCTION

Precise control of the gas pressure in the gas generator is the key element of having
variable thrust for a ducted rocket, and a great example of a general control problem where
time delays, uncertainties and nonlinearities need to be carefully addressed. A ducted rocket
whose thrust is controlled is called a throttleable ducted rocket (TDR, see Fig. 1). Controlling the thrust gives the ability to alter the speed during a mission and allow the missile to
make maneuvers more efficiently. This way, the missile enlarges its no-escape zone, which is
the maximum range that it can outrun its target [1].
As seen in Figure 1, TDR propulsion systems contain two main parts, the gas generator
(GG) and the ram combuster (RC). Fuel-rich (or oxidizer deficient) solid propellant in GG is
ignited and partially burned to obtain fuel in the gaseous form which is then sent to the RC
to meet with air taken in and compressed by the air-intakes and have further combustion to
produce thrust. Variable thrust is achieved by modifying the fuel mass flow rate from the
GG to the RC, which can be done using several different methods including changing the
burning area of the GG propellant in a controlled manner, changing the throat area between
the GG and the RC using a mechanical valve (see Fig. 1), introducing secondary injection
to the GG chamber to control burning rate and utilizing a vortex valve which introduces a
swirl to the flow in order to control the effective throat area [2, 3]. Among these methods
of regulating the fuel mass flow rate, changing the throat area using a mechanical valve is
the most commonly used approach [3–10].
Valve

Gas Generator

Ram Combuster

Fuel enriched
solid propellant
.

.

.

Air flow

.

.

.

.

.

.

.

.
.

.

.

.

Metallic particles

.

. .

Fuel flow

Air Intake
Throat

FIG. 1: Throttleable ducted rocket components

A typical thrust/speed control loop has a hierarchical structure where an outer loop is
2

driven by the error between the desired and the measured thrust/speed and determines the
required GG pressure to minimize this error. Required GG pressure is then provided to
the inner feedback loop as a reference, which controls the mass flow rate to achieve the
desired pressure (see Fig. 2) [1, 11–13]. The GG pressure control loop is the focus of
this paper. Studies in the literature on the GG pressure control problem mainly involves
linearization of the nonlinear system dynamics around certain operating points and designing
linear controllers for the resultant linear dynamics. Among others, proportional integral (PI)
type controllers appear to be the most common controller used for this purpose [12, 14, 15].
Another challenge for the GG control problem, apart from the nonlinear dynamics, is that as
the fuel burns, the free volume inside GG increases with time, which makes the system time
varying. Moreover, the GG dynamics contain several uncertainties emanating from metallic
particles inside the solid fuel, such as deposition at nozzle throat and ablation of mechanical
elements. One of the more sophisticated control approaches to address these issues is gain
scheduling [11, 13]. In this approach, a full knowledge of the controlled plant in terms of a
high fidelity mathematical model is employed to prepare look-up tables which are then used
to assign appropriate controller gains at different operation modes and conditions. Another
method, which is employed for the flight performance evaluation study of the Meteor missile,
is called as the ‘performance funnel’ [16, 17]. In this approach, a proportional controller is
utilized with a time varying gain, which is adjusted online to keep the error of the closed
loop system within a predefined performance funnel.
In this study, we build upon the above mentioned successful approaches by eliminating
the need for a precise system model for classical approaches and eliminating the conservatism
introduced by robust methods. This is achieved by developing and utilizing a delay resistant
adaptive controller which includes a unique combination of the elements of the Adaptive
Posicast Controller (APC) studied by the author Yildiz and his coworkers in [18] and closed
loop reference model (CRM) adaptive controller proposed in [19–21].
APC is an adaptive controller developed for time delay systems which extends the ideas
from the Smith Predictor [22], finite spectrum assignment [23] and adaptation [24–26]. The
ability of the APC to accommodate large delays has been successfully validated through
several simulation and experimental studies presented in [27–32]. Other notable studies on
the adaptive control of time delay systems can be seen in [33], where unknown input delays
and [34], where both state and input delays are addressed. Also, extension of predictor
3

REQUIRED PRESSURE

PRESSURE
CONTROLLER

REQUIRED
THROAT AREA

THROAT AREA

ACTUATOR &
VALVE CONTROL
MECHANISM

FUEL FLOW TO
RAM COMBUSTER
GAS GENERATOR
DYNAMICS

MEASURED PRESSURE

DESIRED
SPEED/THRUST

REQUIRED PRESSURE

SPEED/THRUST
CONTROLLER

RAM COMBUSTER &
MISSILE DYNAMICS

MEASURED SPEED/THRUST

GAS GENERATOR LOOP

FUEL FLOW TO RAM COMBUSTER

FIG. 2: Thrust/speed control loop hierarchical structure

feedback to nonlinear and delay adaptive systems with actuator dynamics modeled by partial
differential equations can be found in [35]. Closed loop reference model (CRM) adaptive
control proposed in [19–21, 36, 37] introduces an error feedback modification of the reference
model which is shown to improve the transient response. Similar approaches where the
reference model is modified to obtain better transients can be found in [38–41]. Therefore, in
this paper, APC and CRM adaptive control constitute the main ingredients of the proposed
controller design which can handle large delays and provide improved transient response.
Furthermore, the projection algorithm [42] is used to prevent parameter drift and a procedure
to select the adaptation and CRM feedback gains, inspired from [27] and [19], are provided to
reduce the tuning time and effort for experimental implementations. The effectiveness of this
approach is validated through extensive simulations together with experiments conducted
using an industry grade cold air test setup, where comparisons with conventional MRAC
and CRM adaptive control are provided.
In TDR research, cold air test setup (CATS) is widely appointed as the crucial step in
validation of subsystems and methods. CATS is used, for example, to validate the numerical
simulation results of flow characteristics, to test the structures that are used to change the
throat area and to characterize the materials that are planned to be used in the construction
[43–48]. Controller is also an important subsystem that is required to be qualified and CATS
4

is utilized to conduct comparative analysis of alternative control systems [49], which is then
used to acquire a proper control methodology based on the gained insight. In this study,
using the facilities provided by Roketsan Inc., we setup and utilize an industrial grade CATS
to validate the proposed controller and compare it by various competing methods.
Organization of the paper is as follows: Section II describes the mathematical modeling of
the system together with model enhancement steps using experimental data. In Section III,
the designs of MRAC, CRM and the proposed DR-CRM adaptive controllers and a PI controller are described with implementation enhancements. Simulation results are presented
in Section IV for two different scenarios. Experimental setup is described and experiment
results for the same scenarios are given in Section V. A summary is provided in Section VI.

II.

SYSTEM MODEL

Overall cold air test setup (CATS) consists of a control volume (pressure chamber), an
actuator, a valve mechanism, drive-train elements, a gas supply and a pressure regulator
(see Fig. 3). A continuous flow of gas is provided by a nitrogen source to the plant from the
inlet and flow is adjusted by a pressure regulator. The output of the model, which is the
pressure inside the control volume, is controlled through changing the exit throat area of
the flow, which is the model input. The throat area is increased/decreased using the linear
motion of a pintle at the exit throat. Drive-train elements are used to convert the rotational
motion of the actuator, which is a brushless direct current motor, to translational motion of
the pintle with required amount of reduction.

A.

Plant

Assuming ideal gas conditions, the difference between the mass flow rates going into the
control volume (pressure chamber), ṁin [kg/sec], and going out of the control volume, ṁout
[kg/sec], is given as
ṁin − ṁout =

Ṗ V
RT

(1)

where R [J/(kg K)] is the specific gas constant and P [Pa], T [K] and V [m3 ] are the pressure,
the temperature and the volume of the gas inside the control volume. The process is assumed
to be isothermal with no change in the control volume, hence V̇ = Ṫ = 0.
5

Air flow outlet
Actuator & valve
mechanism
Throat with variable
open throat area
Gas source

Pintle
Spindle &
driving elements

Pressure regulator

Brushless DC motor

Air flow inlet

CONTROL VOLUME
(PRESSURE CHAMBER)

Pressure
sensor

FIG. 3: Schematic of cold air test setup

Mass flow coming out of the control volume is a function of the throat area (At [mm2 ])
and the pressure inside the control volume (P ). Assuming chocked flow conditions, the
relationship between the throat area and the resulting mass flow rate out of the control
volume is given by [11]
ṁout

γ
2 ( γ−1
)
= P At ((
)
)
γ+1

r

γ
P At
= ∗ ,
∗
RT
c

(2)

where γ is the specific heat ratio of air, T ∗ is the temperature at the throat and c∗ [m/s]
is the characteristic velocity of the gas inside the control volume. Using (1) and (2), it is
obtained that
Ṗ = c1 ṁin − c2 P At = M − c2 P At ,
where c1 =

B.

RT
,
V

c2 =

c1
c∗

(3)

and M = c1 ṁin .

Actuator

A brushless DC motor is used in position controller mode as the actuator, with its driver
card and an encoder to measure the position/speed of the rotor. The closed loop actuator

6

dynamics is approximated as a first order linear system, whose transfer function is given as
Wact (s) =

1
θmes (t)
=
θcom (t)
τact s + 1

(4)

where θmes [quadrature] is the measured actual rotational position of the rotor while θcom
[quadrature] is the commanded rotational position (4000 quadratures (qc) correspond to 1
rotation) and τact is the actuator closed loop time constant.

C.

Valve geometry and drive-train elements

The drive-train elements consist of a gear box and a spindle to convert the rotational
motion of the motor into the translational motion of a pintle at the throat area (see Fig. 4).
The linear position of the pintle determines the throat opening.

Air flow direction

Pintle

a

y

Pintle moving direction

x

CONTROL
VOLUME

y0

r0

Air flow direction

FIG. 4: Valve geometry

The pintle has one degree of freedom in x direction and open throat area changes as the
pintle moves along the x axis due to its conical surface. The cross sectional area of the
cylindrical part at the back of the pintle is smaller than the fixed throat area, which makes
sure that the open throat area At is always larger than zero and protects the system from
rapid pressure build up. Below, we provide the valve and the drive train models.
7

1.

Valve Geometry

There exist a nontrival relationship between the movement of the pintle and the minimum
throat area, where the choked flow conditions occur, due to their complex geometries [43, 45].
The size and the location of the minimum throat area is hard to estimate analytically due
to the fact that location of the choked flow line, where the throat area is minimum, shifts
towards the upstream as the pintle moves into the throat [45]. In this study, the size of the
open throat area is approximated as the projection of the real area on the vertical surface
that is perpendicular to the pintle center line. Movement of the pintle along the x axis
reduces the projected throat area by
y = y0 − tan(α)x,

(5)

where y0 is the radius of the pintle at cylindrical part and α is the half of the cone angle at
the tip of the pintle (see Fig. 4). The projected open throat area is then calculated as
At = (r02 − y 2 )π.

2.

(6)

Drive-Train Elements

A gear box with a reduction ratio of 1 : R1 is used to increase the torque output of the
actuator. The spindle has an R2 [mm] thread pitch, i.e. one turn of rotation corresponds to
R2 [mm] translational motion. Therefore, the relationship between the actuator rotational
position, θ [quadrature], and the linear position of the pintle, x [mm], can be calculated as
x=

θR2
R1 × 4000

(7)

where 4000 quadratures (qc) correspond to 1 rotation.
Using (5-7), it is obtained that
At = (a1 + a2 θ + a3 θ2 )π
where a1 = r02 − y02 , a2 =

D.

2y0 tan(α)R2
R1 ×4000

(8)

2 2
and a3 = −( tan(α)R
).
R1 ×4000

Model Enhancements Using Experimental Data

To improve the fidelity of the system model, open loop experimental tests are performed
and the obtained experimental data is used to adjust model parameters. Firstly, actua8

4

#105
Reference
Experimental result
Model output

3.5

Angular Position [qc]

3
2.5
2
1.5
1
0.5
0
-0.5
0

10

20

30

40

50

60

70

Time [sec]

FIG. 5: Open loop test results and updated model results of the actuator

tor model is updated: Brushless DC motor is commanded to track inputs in the position
controller mode and based on the response of the actuator the time constant τact in (4)
is updated. Experiments also revealed that a considerable amount of time delay exists in
the actuator control loop, which is due to the communication and computation lags. After
adjusting the time constant and incorporating a time delay, the enhanced actuator model
output is compared with the experimental results and the outcomes are presented in Fig. 5,
which shows that the updated model has a good agreement with the test data. To improve
the system model further, parameters in (3) are considered next: R, T , V and c∗ are available for the test conditions with good accuracy, and therefore the values of these parameters
are easily obtained. However mass flow rate (ṁin ) is not always feasible to measure, especially for relatively small flow rate values. Therefore, the mass flow rate going into the
CATS plant is calculated via (3) using steady state pressure values at different operating
points and corresponding throat areas. Several values for ṁin at different operating points
are plotted in Fig. 6 together with a polynomial fit. At low plant pressure, mass flow rates
are nearly constant. However, mass flow rate decreases at higher pressure, because high
back pressure overcomes the mechanical force in the pressure regulator and reduces the flow
rate. Using the polynomial that is fitted to the data in Fig. 6, (3) is updated as
Ṗ = c1 (c3 P 3 + c4 P 2 + c5 P + c6 ) − c2 P At .

9

(9)

1.05
Curve fit
Experimental data

mdotin (normalized)

1
0.95
0.9
0.85
0.8
0.75
0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

1

1.1

Pressure (normalized)

FIG. 6: ṁin calculated from (3) in tests and curve fitted to the data (9)

1.4
Initial model
Updated model
Experimental results

Pressure (normalized)

1.2
1
0.8
0.6
0.4
0.2
0
25

30

35

40

45

50

Time (sec)

FIG. 7: Open loop responses of experimental setup and simulation with initial plant
model (3) and updated plant model (9)

Open loop simulation results with the overall updated system model along with experimental results, which are obtained for a range of operating points, are given in Fig. 7. It is
noted that the model enhancements can be improved further by making comparisons at several other operating points followed by further tuning of the parameters but it is determined
that this level of fidelity is enough for simulation evaluation purposes.
10

E.

Modeling For Controller Design

As explained earlier, the nonlinear model of the CATS developed in the earlier sections
is used to evaluate controller alternatives in the simulation environment. To facilitate the
controller design, a simpler model is developed following the steps listed below:
• Nonlinear plant model and valve equation in (3) and (8), respectively, are linearized.
• Actuator dynamics are ignored due to small time constants compared to the plant.
Linearizing (3) around an equilibrium point (P, At ) = (P0 , At0 ), it is obtained that
Ṗ = ∆Ṗ = (−c2 At0 )∆P − (c2 P0 )∆At .

(10)

where ∆P = P − P0 and ∆At = At − At0 . Defining ap ≡ −c2 At0 and bp ≡ −c2 P0 , (10) can
be rewritten as
∆Ṗ = ap ∆P + bp ∆At .

(11)

The value of a3 in (8) is much smaller than a1 and a2 for meaningful physical parameters,
and therefore (8) is approximated as
At ≈ (a1 + a2 θ)π.

(12)

It is noted that the valve equation (12) is used to convert the required throat area determined
by the pressure controller to the required actuator rotational position, which is provided to
the actuator/valve control mechanism as a reference (see Fig. 2).
Therefore, together with the actuator time lag, τ , the system model used to develop the
controller is obtained that
bp e−sτ
∆P (t)
=
Wp (s) =
∆At,com (t)
s − ap
III.

(13)

CONTROLLER DESIGN

The structure of the closed loop control system of a conventional throttleable ducted
rocket is provided in Fig. 2. The outer loop determines the required gas generator (GG)
pressure to obtain a desired thrust/speed profile, and the required pressure becomes the
reference for the inner loop pressure controller. The pressure controller is the main focus of
this paper and will be referred to as “the controller” in the following sections. 4 different
11

controllers are evaluated both in the simulation environment and in the experimental set-up:
Model reference adaptive controller (MRAC), closed loop reference model (CRM) adaptive
controller, the delay resistant closed loop reference model (DR-CRM) adaptive controller,
which is introduced in this paper, and a PI controller. These controllers are explained below
in the given order. It is noted that although MRAC is a very well known controller, its
development is also provided here to explain the sequential controller improvement steps
starting from MRAC and ending with the DR-CRM adaptive control.

A.

Model reference adaptive controller (MRAC)

Consider following plant dynamics
yp (t) = Wp (s)u(t),

Wp (s) =

kp Zp (s)
,
Rp (s)

(14)

where yp ∈ < and u ∈ < are the measured output and the control input of the system,
respectively. Zp (s) and Rp (s) are monic polynomials with orders of m and n and kp ∈ < is
the constant gain of the plant. Following assumptions are made for the plant [50]:
• System order n is known along with the relative degree n∗ = n − m.
• Sign of kp is known.
• Polynomial Zp (s) is Hurwitz.
The reference model, which gives the desired response of the closed loop system, is given
as
ym (t) = Wm (s)r(t),

Wm (s) =

km Zm (s)
,
Rm (s)

(15)

where ym ∈ < and r ∈ < are the output of the reference model and bounded reference
signal, respectively. Wm (s) is chosen as strictly positive real with relative degree equal to
the relative degree of the plant.
State space description of the plant and the signal generators for the output feedback
problem with controllable (Λ, bλ ) pair are given as
ẋp (t) = Ap xp (t) + bp u(t),
ω̇1 (t) = Λω1 (t) + bλ u(t)
ω̇2 (t) = Λω2 (t) + bλ yp (t)
12

yp (t) = hTp xp (t)
(16)

where xp ∈ <n , ω1 ∈ <n−1 , ω2 ∈ <n−1 , Ap ∈ <n×n , bp ∈ <n , hp ∈ <n , Λ ∈ <(n−1)×(n−1) is
Hurwitz and bλ ∈ <n−1 .
It can be shown there exist constant parameters θ0 ∈ <, θ1 ∈ <n−1 , θ2 ∈ <n−1 and θr ∈ <
such that the controller given as
u(t) = θ0 yp (t) + θ1T ω1 (t) + θ2T ω2 (t) + θr r(t)

(17)

satisfies the desired reference model response characteristics [50]. When the plant parameters
are unknown, for n∗ = 1, the following adaptation law
Θ̇(t) = −Γsign(kp )e1 (t)Ω(t)

(18)

where


θ0 (t)


θ (t)
1


Θ(t) = 
,
θ2 (t)


θr (t)



yp (t)





ω (t)
1


Ω(t) = 
,
ω2 (t)


r(t)

(19)

e1 (t) = yp (t) − ym (t) is the tracking error, Γ ∈ <2n×2n is a diagonal matrix with positive
elements, stabilizes the closed loop system and ensures that e1 → 0 as t → ∞ [50].

B.

Closed-loop reference model (CRM) adaptive controller

The reference model in the classical model reference adaptive control (15) is unaffected
by the tracking error. In CRM adaptive controller, however, the tracking error e1 (t) =
yp (t) − ym (t) is fed back to the reference model for the purpose of improving the transient
dynamics [19]. Consider the following state space representation of the reference model
dynamics given in (15)
ẋm (t) = Am xm (t) + bm r(t),

ym (t) = hTm xm (t),

(20)

where Am ∈ <n×n , bm ∈ <n and hm ∈ <n . In classical model reference adaptive control,
Zm (s)
Am , bm and hm are chosen such that the transfer function hTm (sI−Am )bm = Wm (s) = km R
m (s)

becomes strictly positive real. In CRM adaptive controller, the reference model is modified
as
ẋm (t) = Am xm (t) + bm r(t) + L(yp (t) − ym (t)),
13

ym (t) = hTm xm (t),

(21)

where L ∈ <n is a design parameter vector. The relationship between the reference model
output ym , the reference r and the tracking error e1 then becomes
ym (t) = Wm (s)r(t) + WL (s)e1 (t)

(22)

where
hTm (sI − Am )L = kL

ZL (s)
= WL (s).
Rm (s)

(23)

Polynomial ZL (s) is order of n − 1. Boundedness of all the signals in the closed loop system
along with the convergence of the tracking error as in classical MRAC is valid for the CRM
approach, using the same controller structure (17) and the adaptive law (18), as long as L
is chosen such that the transfer function
We =

Zm (s)
Rm (s) − kL ZL (s)

(24)

is strictly positive real (see [36]).

C.

Delay resistant closed loop reference model (DR-CRM) adaptive controller

Consider the following plant with an input time delay
yp (t) = kp

Zp (s)
u(t − τ ) = Wp (s)u(t − τ )
Rp (s)

(25)

where yp ∈ < is the measured output, u ∈ < is the control signal, τ is the known time delay,
Zp (s) and Rp (s) are monic coprime polynomials with orders of m and n, respectively, and
kp ∈ < is the constant gain of the plant. Following assumptions are made for the plant:
• System order n is known along with the relative degree n∗ = n − m.
• Sign of kp is known.
• Polynomial Zp (s) is Hurwitz.
The reference model dynamics are given with the closed loop reference model structure
as
ẋm (t) = Am xm (t) + bm r(t − τ ) + L(yp (t) − ym (t))
ym (t) =

hTm xm (t)
14

(26)

where xm ∈ <n , ym ∈ <, r ∈ <, Am ∈ <n×n and bm , L, hm ∈ <n . Input-output relationship
of the closed loop reference model is given as
ym (t) = Wm (s)r(t − τ ) + WL (s)e1 (t)

(27)

where e1 = yp − ym is the tracking error. The transfer functions describing the closed loop
reference model are

Zm (s)
= Wm (s)
Rm (s)
ZL (s)
hTm (sI − Am )L = kL
= WL (s)
Rm (s)

hTm (sI − Am )bm = km

(28)

where Rm (s) is a monic polynomial with order n while Zm (s) and ZL (s) are two monic
polynomials with order n − 1, km ∈ < and kL ∈ < are the gains of the transfer functions.
It is noted that under model matching conditions the tracking error becomes zero, which
reduces the reference model (27) to
ym (t) = Wm (s)r(t − τ ).

(29)

State space description of the plant (25) and the signal generators for the output feedback
problem with the controllable (F, g) pair are given as
yp (t) = hTp x(t)

ẋp (t) = Ap xp (t) + bp u(t − τ ),

(30)

ω̇1 (t) = F ω1 (t) + gu(t − τ )
ω̇2 (t) = F ω2 (t) + gyp (t)

where xp ∈ <n is the state vector, ω1 , ω2 ∈ <n , Ap ∈ <n×n , bp , hp ∈ <n , F ∈ <n×n is
Hurwitz and g ∈ <n . It can be shown that there exists constant controller parameters
β1∗ ∈ <n , β2∗ ∈ <n and k ∗ ∈ < such that the controller
u(t) = β1∗T ω̄1 (t) + β2∗T ω̄2 (t) + k ∗ r(t)

(31)

where ω̄1 (t) , ω1 (t + τ ), ω̄2 (t) , ω2 (t + τ ) together with the rewritten plant dynamics
x̄˙ p (t) = Ap x̄p (t) + bp u(t),
ω̄˙ 1 (t) = F ω̄1 (t) + gu(t)

ȳp (t) = hTp x̄p (t)
(32)

ω̄˙ 2 (t) = F ω̄2 (t) + g ȳp (t)
where x̄p (t) , xp (t + τ ) and ȳp (t) , yp (t + τ ), satisfy the model matching conditions [18].
15

It is shown in [50] that the plant output, yp (t), can be expressed as a linear combination
of ω1 (t), ω2 (t) as
yp (t) = cT ω1 (t) + dT ω2 (t)

(33)

where c, d ∈ <n . When (33) is substituted into (30), it is obtained that




ω̇1 (t)
ω1 (t)

 = A
 + bu(t − τ ),
ω̇2 (t)
ω2 (t)


(34)



 
F
0
g
 and b =  . Non-casual
where A ∈ <2n×2n and b ∈ <2n are given as A = 
gcT F + gdT
0
terms in (31) can then be calculated as




Z 0
ω̄1 (t)
ω1 (t)
Aτ

=e 
+
eAη bu(t + η)dη.
(35)
−τ
ω̄2 (t)
ω2 (t)
When (35) is substituted into (31), the control signal becomes
Z 0
∗T
∗T
φ∗ (η)u(t + η)dη + k ∗ r(t)
u(t) = α1 ω1 (t) + α2 ω2 (t) +

(36)

−τ
∗
∈ <n , φ∗ (η) ∈ < are the corresponding controller parameters which eliminate the
where α1,2

non-causality in controller (31) using (35).
In the case of unknown plant parameters, the control input, u(t), can be split into two
sub-signals as
u(t) = u1 (t) + u2 (t)
where
u1 (t) =

α1∗T ω1 (t)

+

α2∗T ω2 (t)

(37)

0

Z
+

φ∗ (η)u(t + η)dη + k ∗ r(t)

(38)

φ̃(t, η)u(t + η)dη + k̃(t)r(t)

(39)

−τ

and
u2 (t) =

α̃1T (t)ω1 (t)

+

α̃2T (t)ω2 (t)

Z

0

+
−τ

where α̃i (t) = αi (t) − αi∗ for i = 1, 2; φ̃(t, η) = φ(t, η) − φ∗ (η) and k̃(t) = k(t) − k ∗ .
Substituting (37), using (38) and (39), into (30), the closed loop dynamics is obtained as
"
#
Z 0
Ẋp (t) = Amn Xp (t) + bmn θ̃(t − τ )w(t − τ ) +
φ̃(t − τ, η)u(t − τ + η)dη + k ∗ r(t − τ )
−τ

yp (t) = hTmn Xp (t)
(40)
16


 
∗
∗
A
b
β
b
β
p 1
p 2

 bp 
 p
h
i

 

where Amn =  0 F + gβ1∗ gβ2∗  , bmn =  g  , hTmn = hTp 0 0 ,

 

T
ghp
0
F
0
h
iT
h
iT
h
i
Xp (t) = xTp (t) ω1T (t) ω2T (t) , w(t) = ω1T (t) ω2T (t) r(t) and θ̃(t) = α̃1 (t) α̃2 (t) k̃(t) .
˜ = 0, the closed loop dynamics
It is noted that when the parameter errors are zero, i.e. (·)


represented by (40) becomes equivalent to that of the reference model dynamics (29), which
shows that the system formed by (Amn , bmn , hmn ) is a non-minimal representation of the
reference model:
kp
Zm (s)
=
Wm (s).
Rm (s)
km
Using (40) and (41), the plant output is obtained as
#
"
Z 0
kp
φ̃(t − τ, η)u(t − τ + η)dη + k ∗ r(t − τ ) .
yp (t) =
Wm (s) θ̃(t − τ )w(t − τ ) +
km
−τ
hTmn (sI − Amn )bmn ≡ kp

Subtracting (27) from (42), it can be obtained as
#
"
Z 0
kp
φ̃(t − τ, η)u(t − τ + η)dη − WL e1 (t).
e1 (t) =
Wm (s) θ̃(t − τ )w(t − τ ) +
km
−τ
Solving (43), the tracking error can be found as
#
"
Z 0
Zm (s)
e1 (t) = kp
φ̃(t − τ, η)u(t − τ + η)dη
θ̃(t − τ )w(t − τ ) +
Rm (s) − kL ZL (s)
−τ
where

Zm (s)
Rm (s)−kL ZL (s)

(41)

(42)

(43)

(44)

= We (s) has sufficient degrees of freedom, in terms of the design param-

eter vector L, to be determined as a SPR transfer function (see Section III B).
The closed loop reference model in (26) can be rewritten as,
ẋmn (t) = Amn xmn (t) + bmn k ∗ r(t − τ ) + GL(yp (t) − ym (t))

(45)

hTmn xmn (t)

for xmn

ym (t) =
h
iT
∗T
∗T
∗T
= xp (t) ω1 (t) ω2 (t) where x∗p (t), ω1∗ (t) and ω2∗ (t) are the signals in the refer-

ence model, corresponding to the signals xp (t), ω1 (t) and ω2 (t) in the closed loop dynamics,
respectively. G ∈ <3n×n is the constant matrix to transform xm to the controllable subspace
in xmn (see [36]). Error dynamics, e(t) = Xp (t) − xmn (t), in non-minimal form, is found by
subtracting (45) from (40) as
"

Z

ė(t) = Ae e(t) + bmn θ̃(t − τ )w(t − τ ) +

φ̃(t − τ, η)u(t − τ + η)dη
−τ

e1 (t) = hTmn e(t)
17

#

0

(46)

where
Ae = Amn − GLhTmn .

(47)

It can be shown [18], utilizing the error dynamics (46), that the controller (37)-(39), along
with the adaptation laws given as
˙ = θ̇(t) = −sign(k )Γ e (t)ω(t − τ )
θ̃(t)
p
θ 1
˙ η) = φ̇(t) = −sign(k )Γ e (t)u(t − τ + η)
φ̃(t,
p
φ 1

(48)
−τ ≤η ≤0

stabilizes the closed loop system as long as L is chosen to ensure that We (s) is SPR. Furthermore, tracking error e1 (t) converges to zero. Γθ ∈ <(2n+1)×(2n+1) is a diagonal matrix
with positive elements and Γφ ∈ <+ .

D.

Proportional-integral controller (PI)

A proportional-integral controller is designed for the plant dynamics given in (14). The
controller transfer function is given by
GP I (s) = Kp (1 +

1
)
Ti s

(49)

where scalar controller parameters Kp and Ti are selected to make the closed loop dynamics
provide a similar response with the reference model (15) used for MRAC. Frequency response
plots of the reference model and the compensated system are provided in Fig. 8.

E.

Implementation Enhancements

The implementation of the adaptive controllers requires some modifications to address
the issues that were not taken into consideration during the initial design but arise in real
experimental tests. Below we explain these experimental requirements and how they are
addressed.

1.

Disturbance rejection

The goal of the model reference adaptive control designs presented in Sections III A, III B,
III C is to force the plant output follow the reference model output while keeping all the
18

10
Compensated closed loop system
Reference model

Magnitude (dB)

5
0
-5
-10
-15
-20
10-2

10-1

100

101

102

101

102

Frequency (rad/sec)

0

Phase (deg)

-20

-40

-60

-80
10-2

10-1

100

Frequency (rad/sec)

FIG. 8: Comparison of the Bode plots of the compensated closed loop system with
the designed PI controller and the reference model

system signals bounded. The disturbances, on the other hand, are not explicitly taken into
account. In throttleable ducted rocket propulsion systems, the solid propellant in the gas
generator contains metallic particles (see Fig. 1), which can cause deposition or ablation at
the throat between the gas generator and the ram combuster, which acts like an additive
disturbance on the effective throat area. Therefore, it is important to incorporate disturbance rejection capabilities into the gas generator pressure controller. Since the plant model
(13) developed for the controller design is first order, the disturbance rejection modification
explained below utilizes a scalar plant model, which simplifies the process.
Consider the plant
yp (t) = ap yp (t) + bp (u(t) + d0 )
19

(50)

where yp is the system output, d0 is the unknown constant disturbance, u is the plant input
and ap and bp are the plant parameters. An additional adaptive controller parameter, θ3 (t),
is introduced to the controller structure for MRAC and CRM adaptive controller, which
results in the following control signal
u(t) = θ0 (t)yp (t) + θr (t)r(t) + θ3 (t)

(51)

where theta0 , thetar and θ3 are the adaptive control parameters to be determined and r is
a bounded reference signal. The adaptation laws for the controller parameters are given as
˙
Θ̄(t)
= −sign(bp )Γ̄e1 (t)Ω̄(t)
(52)
h
i
h
i
T
T
where Θ̄ (t) = θ0 (t) θr (t) θ3 (t) , Ω̄ = yp (t) r(t) 1 , Γ̄ > 0 is the diagonal adaptation
rate matrix and e1 is the tracking error given by e1 = yp − ym . The reference model for the
MRAC is appointed as
ẏm (t) = am ym (t) + bm r(t),

(53)

and the reference model for CRM adaptive controller is selected as,
ẏm (t) = am ym (t) + bm r(t) + l(yp (t) − ym (t)),

(54)

where am , bm , l ∈ <, am = −bm for unity DC gain and am < 0 and l > 0. It is noted that
for a first order plant, choosing l > 0 satisfies the stability conditions.
Similarly, DR-CRM adaptive control signal is modified as
Z 0
u(t) = αy (t)yp (t) +
λ(t, η)u(t + η)dη + k(t)r(t) + θ3 (t).

(55)

−τ

The adaptive laws for the controller parameters αy , λ, k, θ3 ∈ < are given as
θ̇(t) = −sign(bp )Γθ e1 (t)ω(t − τ )
(56)
∂λ
(t, η) = −sign(bp )γλ e1 (t)u(t + η − τ ) − τ ≤ η ≤ 0
∂t
h
i
h
i
where θT (t) = αy (t) k(t) θ3 (t) , ω T (t) = yp (t) r(t) 1 , Γθ > 0 and γλ ∈ <+ are adaptation rate constants.
The reference model for the DR-CRM adaptive controller is determined as
ẏm (t) = am ym (t) + bm r(t − τ ) + l(yp (t) − ym (t)),
where am , bm , l ∈ <, am = −bm for unity DC gain and am < 0 and l > 0.
20

(57)

2.

Robustness

The design of the DR-CRM adaptive controller presented in Section III C portrayed an
idealized case, where the delay free part of the plant dynamics are assumed to be linear and
time-invariant with unknown but constant parameters. Furthermore, the measurements are
assumed to be perfect. However, it is known that in reality, no plant is truly linear or finite
dimensional. Parameters may vary with time and operating conditions, and measurements
are always contaminated with noise. The plant model used for the controller design is always
an approximation of reality. Therefore, we need a robustfying modification against possible
parameter drifts in adaptive controller parameters, due to these non-ideal situations. One
common remedy utilized to prevent parameter drift is the projection algorithm [51], which
is explained below.
Consider a controller parameter vector Θ(t) ∈ <k . The upper bound on the controller
parameter vector norm ||Θ||, to activate the projection modification is defined as Θmax ∈ <+ .
A continuous and differentiable convex function f : <k → < is introduced as
f (Θ) =

||Θ||2 − Θ2max
Θ2max

(58)

where  > 0 is a constant which defines the projection tolerance region where ||Θ|| ≤
√
Θmax ( 1 + ) forms a hard bound on the parameter norm. The projection operator, Proj :
<k × <k → <k , is described as

T

y − ∇f (Θ)(∇f (Θ))
yf (Θ), if ||Θ|| > Θmax ∧ y T ∇f (Θ) > 0
||∇f (Θ)||2
Proj(Θ, y) ,

y,
otherwise

(59)

where y ∈ <k , ∇f (Θ) ∈ <k is the gradient vector of f , evaluated at Θ.
The adaptive control law is then modified as
Θ̇ = Proj(Θ, −sign(kp )Γe1 Ω)

(60)

where Γ ∈ <k×k is a diagonal matrix with positive elements, e1 ∈ < is the tracking error,
Ω ∈ <k is the regressor vector containing system signals. It can be shown that with the
adaptive law (60) utilizing projection, the parameter drift is prevented [51].
Although we observed in the experiments that the introduction of the projection algorithm prevents parameter drift, there are no well defined procedures to determine the upper
21

bound Θmax on the controller parameter vector Θ. One method is to calculate this bound
using the worst case uncertainty case. Another method is to conduct several experiments
without projection and observe the variation of controller parameters, which can help define a reasonable upper bound. A third approach is setting initial values for the controller
parameters that would satisfy the matching conditions for the nominal plant dynamics and
then determining the upper bound for the parameters as a certain percentage higher than
these initial values. In our experiments, we utilized the second method.

3.

Digital implementation of the integral term in DR-CRM adaptive controller

For computer implementation purposes, the finite integral term in (55) is approximated
as,
Z

0

λ(t, η)u(t + η)dη =
−τ

m
X

λi (t)u(t − idt) = λ̄T (t)ū(t)

(61)

i=1

where dt = 50 ms is the sampling interval of the implemented software, τ = 300 ms is
τ
= 6. λ̄ ∈ <m is the vector containing parameters λ̄T (t) =
the time delay and m = dt
h
i
h
i
m
T
and
ū
∈
<
is
the
delayed
input
vector
ū
(t)
=
λ1 (t) . . . λm (t)
u(t − dt) . . . u(t − mdt) .

Adaptation law (56) is updated as
˙ = −sign(b )Γ̄ e (t)ω̄(t)
θ̄(t)
p
θ̄ 1

(62)

where


αy (t)





yp (t)







 λ (t) 
 u(t − dt) 
 1 


 . 


..
 .. 


.




θ̄(t) = 
 ω̄(t) = 

λm (t)
u(t − mdt)








 k(t) 


r(t)




θ3 (t)

(63)

1

and Γ̄θ̄ > 0 is a diagonal adaptation rate matrix.

4.

Initialization of the controller parameters

One way to initialize the adaptive controller is to use initial values that would satisfy
model matching for the nominal plant dynamics. However, in MRAC and CRM adaptive
22

controller designs, initial parameters θ0 (0) and θr (0) had to be lowered to obtain the best
performance since the presence of the time delay prevents exact model matching due to the
lack of exact delay compensation terms. On the other hand, the initial parameters for DRCRM adaptive controller are set to satisfy the model matching conditions for the nominal
plant dynamics. Initial parameter θ3 (0) is selected zero for all adaptive controllers.

5.

Selecting the design parameters for MRAC

Adaptation rate for a particular parameter θi in MRAC is chosen according to the empirical formula [52]
Γ̄ii =

|θi∗ |
3τm (r̄)2

(64)

where θi∗ is the ideal value of the controller parameter, τm is the smallest time constant of
the reference model and r̄ is the maximum possible amplitude of the reference signal which
is equal to the amplitude of the operating range. Since the ideal controller parameters, θi∗ ,
are assumed to be unknown, nominal values calculated using the matching conditions are
used instead.
Adaptation rates obtained from (64) are calculated for the worst case scenario which is
usually valid at the beginning of the operation when the tracking error and the system states
are of the same order of magnitude with the reference signal. Moreover, the calculation (64)
requires the estimation of the ideal control parameters. Due to these approximations, we
include a fine-tuning matrix W for MRAC to fine-tune the adaptation rates as Γ̄W = Γ̄W
where



p1 0 0


W = 0 p2 0


0 0 p3

(65)

and adjustable constants p1,2,3 are used for the fine-tuning process. It was experienced during
simulations that selection of p1 = p2 = 1 and p3 > 1 provide faster and more robust system
response.

23

6.

Selecting the design parameters for CRM and DR-CRM adaptive controllers

CRM gain, `, helps suppress the oscillations in the case of high adaptation rates. However,
a numerically large ` can cause so-called “peak phenomena” in non-zero initial tracking
error [53]. Therefore, a procedure to determine the optimum value of the CRM gain and
adaptation rates is needed to reduce the time and effort spent for the controller tuning.
Below steps, inspired from [19], are followed separately for CRM and DR-CRM adaptive
controllers:
1. Find the adaptation rates using (64) and define the adaptation rate vector, γ̄.
2. Find the norm of the adaptation rate vector, ||γ̄||.
3. Define the CRM gain, `, equal to ||γ̄||.
4. Increase ` and γ̄ together by keeping ` = ||γ̄|| until a desired tracking performance is
obtained in numerical simulations.
The same adaptation rate, γλ , which is determined through simulations, is chosen for
each λi since they have same order of magnitude.

F.

Step by step controller design procedure

A clear, step by step procedure is provided below to facilitate the adaptive controller
design.
Step 1: Determine the reference model dynamics by choosing appropriate (Am , bm , hm ) for
the performance specifications of the closed loop system. (Am , bm , hm ) should form a
strictly positive real (SPR) transfer function for MRAC.
Step 2: Signals ω1,2 (t) in (16) and (30) are generated by choosing (F, g) and (Λ, bΛ ) so that they
form controllable pairs. Since these signals operate like state observers, eigenvalues of
the matrices F and Λ should be faster than the reference model dynamics. For a first
order plant, these signals are not required.
Step 3: Set the initial conditions of the controller parameters using Section III E 4.

24

Step 4: Find the adaptation rates and CRM gain, `, using Section III E 5-III E 6. CRM gain
should be selected such that We (s) in (24) is an SPR transfer function for both CRM
and DR-CRM adaptive controllers.
Step 5: Tune the parameters γλ and p3 in (63) and (65) using numerical simulations and then
fine-tune again during experiments.
Step 6: Integrate the projection algorithm provided in (60) to the adaptation laws.
All of the adaptive controllers considered in this paper requires minimal amounts of
computational resources and memory. DR-CRM, having the highest number of terms in
the control signal, for example, needs only 256 bytes of memory for data storage. It has
116 operations per cycle, which corresponds to around 2320 floating point operations per
second (flops). See the Appendix for detailed memory requirement and computational load
calculations.

IV.

SIMULATIONS

Simulation results employing the full nonlinear model, which is developed in Section II,
are presented in this section. All the parameters defining the controller (initial conditions of
the controller parameters, adaptation rates, CRM gain, (F, g) and (Λ, bΛ ) pairs) are obtained
using Section III F. Reference model is chosen to satisfy the performance specifications listed
in Table I. The controller gains of the PI controller is designed so that the closed loop system
provide a similar response with the reference model, as explained in Section III D.
TABLE I: Specifications of chosen reference model
Steady state error Rise time Settling time (5 %) Maximum overshoot
0%

0.6 sec

1.5 sec

10 %

Firstly, performance of MRAC is compared with the PI controller. Then, MRAC, CRM
adaptive controller and DR-CRM adaptive controller are comparatively evaluated, where
demanding tracking tasks are utilized to reveal the performance differences between these
controllers. Numerical simulations are carried out using Matlab
of 50 ms.
25

R

with a sampling interval

A.

MRAC vs. PI controller

Reference
PI
MRAC

Pressure (normalized)

1.2

1

0.8

0.6

0.4

0.2
25

30

35

40

45

50

55

60

65

Time (sec)

FIG. 9: Tracking curves for the PI controller and MRAC at 3 different operating
points.

In Fig. 9, simulation results are given demonstrating the performances of MRAC and
the PI controller. It is noted that although the PI controller and MRAC show very similar
performances around the normalized nominal pressure, between 0.6 and 0.8, the closed loop
system with the PI controller presents an oscillatory response for higher pressure operating
conditions and a slow response for lower pressures. On the other hand, MRAC can adapt
itself to changing operating conditions and provides a more consistent performance across
operating points.

B.

Comparative evaluation of adaptive controllers

The promise of CRM and delay compensation modifications for the conventional MRAC
controller is that they provide higher performance without causing excessive oscillations. To
demonstrate the difference that these modifications can make, more demanding reference
trajectories, compared to those of the case presented in Fig. 9, are used for the reference
tracking tests. Also, time constant of the reference model, τm , is halved to obtain a fast
response. Adaptation rates of MRAC are updated using (64) and a new p3 , employed in
26

(65), is obtained empirically. The results of these simulations are given in Fig. 10-12.
1.6
MRAC
CRM
Reference
DR-CRM

Pressure (normalized)

1.4
1.2
1
0.8
0.6
0.4
0.2
20

30

40

50

60

70

Time (sec)

FIG. 10: Reference tracking of MRAC, CRM adaptive control and DR-CRM
adaptive control in simulations.

Throat area command (normalized)

1.2
CRM
MRAC
DR-CRM

1
0.8
0.6
0.4
0.2
0
20

30

40

50

60

70

Time (sec)

FIG. 11: Evolution of control inputs of MRAC, CRM adaptive control and
DR-CRM adaptive control in simulations.

MRAC, CRM adaptive controller and DR-CRM adaptive controller provide similar performances for small variations in the pressure demand but MRAC response becomes os27

12
CRM
DR-CRM

Projection tolerance region
defined by 0

11

# max

10

MRAC

||#||

9
8
7
6
5
4
3
20

30

40

50

60

70

Time (sec)

FIG. 12: Controller parameters of MRAC, CRM adaptive control and DR-CRM
adaptive control in simulations with the projection boundary

cillatory once the demanded variation is increased by 3 times. Although CRM adaptive
controller provides a considerably more damped response compared to MRAC, the best response is obtained for the case where DR-CRM adaptive controller is employed. Figure 11
shows that DR-CRM controller provides the smoothest control input. In Fig. 12, adaptive
control parameters are shown where MRAC parameters hit the projection boundary and
stays within the projection tolerance limits.

V.

EXPERIMENTS

The experimental results are obtained using a cold air test setup (CATS), which is designed and manufactured by Roketsan Missiles Inc. (See Fig. 3 for a schematic of the overall
system and Fig. 13 for the real experimental system). The pressure chamber in CATS has
two ports, which are called the inlet port and the exit port. The inlet port is connected
to a nitrogen gas source of 230 bars through a pressure regulator. The pressure regulator
ensures safe test conditions by adjusting the inflow pressure. There is a solenoid between the
pressure chamber and the pressure regulator to stop the flow in case of emergencies. Outlet
port of the pressure chamber has a shape of a nozzle whose effective throat area is continuously altered during the operation by the actuator and valve mechanism. An EC-max 30,
28

Pressure
chamber

Exit port
Inlet port
Solenoid

Pressure
transducer

Pressure
regulator
Nitrogen
gas source

Master
computer
Power
source

Compact Rio
Actuator
driver

Ethernet
connection

CAN
connection
NI 9478, NI 9881 and NI 9203
FIG. 13: Cold Air Test Setup.

60 Watt, 24 Volt brushless direct current motor (Maxon Motor Company

R

) with a EPOS2

70/10 driver is used as the actuator. Output shaft of the motor is connected to a gear box
and a spindle, respectively. Other side of the spindle is connected to a conical pintle which is
located such that its linear position determines the effective throat area at the outlet of the
pressure chamber (see Fig. 4). There is a pressure transducer inside the pressure chamber,
which provides real time pressure data to the controller. A slave Compact Rio computer
(National Instrument

R

), running LabVIEW software, is used to collect data from the

pressure transducer, run the pressure controller cycles to calculate the necessary effective
throat area and send this data to the actuator driver which is responsible for controlling the
actuator position. The communication between the slave Compact Rio computer and the
actuator driver card is through NI 9881 card using a CAN bus protocol. Pressure transducer
29

Solenoid

Actuator

Laptop (master)
Running Labview

NI 9881

Ethernet
(FIFO)

NI 9478
Compact Rio (slave)
Running Labview

CAN

EPOS
(driver)

Motor
Encoder

Data acquisition card
(NI 9203)

Pressure transducer

FIG. 14: Schematic of the CATS hardware and data communication.

is connected to the data acquisition card (NI 9203) that is also connected to the Compact
Rio. All the algorithms for data acquisition, pressure controller calculations, data sending
and corresponding communication phases are prepared in the master computer prior to the
experiments using LabVIEW and the code is embedded to the slave Compact Rio computer
through ethernet connection. Communication between the master and slave computers is
realized with FIFO: First-In-First-Output methodology. Master computer monitors the experiment in real time and is able to intervene in the experimental process in case of a safety
hazard. The solenoid located in between the pressure chamber and the pressure regulator is
controlled by the master computer through the digital card NI 9478. A detailed schematic
of the hardware and data communication are given in Fig. 14.
The same scenarios used for the simulations are employed for the experimental tests.
First, better performance of MRAC over the PI controller is demonstrated by performing
experiments at three different operation points with the same controller gains and parameters
used in the simulations. Then, a comparative evaluation showing the advantage of DR-CRM
adaptive controller over other adaptive controllers are presented. Finally, an experiment is
conducted for a larger period of time to show the effectiveness of the projection algorithm.
All the numbers in the figures are normalized.

30

1.2

MRAC
PI
Reference

Pressure (normalized)

1.1
1
0.9
0.8
0.7
0.6
0.5
0.4
10

15

20

25

30

35

40

Time (sec)

FIG. 15: Test results of PI controller and MRAC for three different operating
conditions. The adaptation rates used for these experiments are the same as the
ones used in simulations.

A.

MRAC vs. PI controller

Test results are given in Fig. 15. PI controller shows acceptable performance at the
nominal operating point, the linearized model of which was used for the controller design.
However, as the operating point deviates from the nominal design conditions, advantage of
the adaptive controller, which provides consistent transient performance at different operating conditions, is observed.

B.

Comparative evaluation of adaptive controllers

To demonstrate the performance differences between the adaptive controllers, a more
challenging reference signal, compared to that of the previous subsection, is used and time
constant of the reference model, τm , is halved. Same controller design parameters are used
as in the simulations.
Experiment results are given in Fig. 16-18, exhibiting similar trends with the simulations with CRM and DR-CRM adaptive controller having slightly more damped responses.
DR-CRM adaptive controller is able handle the demanding operation conditions in the ex-

31

DR-CRM
Reference
CRM
MRAC

Pressure (normalized

1.4
1.2
1
0.8
0.6
0.4
0.2
20

30

40

50

60

70

Time (sec)

FIG. 16: Reference tracking of MRAC, CRM adaptive control and DR-CRM
adaptive control in experiments.

Throat area commanded (normalized)

1.4
DR-CRM
CRM
MRAC

1.2
1
0.8
0.6
0.4
0.2
0
-0.2
20

30

40

50

60

70

Time (sec)

FIG. 17: Evolution of control inputs of MRAC, CRM adaptive control and
DR-CRM adaptive control in experiments.

periments and provides a reasonable performance. CRM adaptive controller, damping most
of the oscillations as intended, results in undesired high amplitude overshoots. MRAC’s
response is similar to that of the other two adaptive controllers for small pressure deviation
demands but becomes oscillatory for larger deviations in the reference signal. DR-CRM
32

12

# max

Projection tolerance region
defined by 0

MRAC
CRM
DR-CRM

10

||#||

8

6

4

2
20

30

40

50

60

70

Time (sec)

FIG. 18: Controller parameters of MRAC, CRM adaptive control and DR-CRM
adaptive control in experiments and projection boundary

adaptive controller has the smoothest controller input, exhibited in Fig. 17. Fig. 18 shows
that the controller parameters of MRAC hit the projection boundary, but then prevented
to grow further.
Effect of the projection algorithm can be observed more clearly in the experiments that
are conducted for longer times. Fig. 19 presents the evolution of the norm of the MRAC
controller parameters in a longer test whose tracking curve is depicted in Fig. 20. Controller parameters tend to increase due to non-ideal situations such as unmodeled dynamics,
disturbances and noise, but the projection algorithm keeps them within a predefined bound.

VI.

SUMMARY AND FUTURE WORK

The pressure control of gas generators in ducted rockets are addressed in this paper. To
solve this control problem, which includes time delays, uncertainties and nonlinear dynamics,
a delay resistant closed loop reference model (DR-CRM) adaptive controller is proposed. The
controller merges the benefits of two different approaches: The first approach is the adaptive
posicast controller (APC) which compensates the time delay by making use of positively
forecasted output of the plant. The second approach is the closed-loop reference model
(CRM) modification, which damps the high frequency oscillations due to high adaptive
33

4.5
4
MRAC

||#||

3.5

Projection tolerance region
defined by 0

# max

3
2.5
2
1.5
20

40

60

80

100

120

140

160

180

200

Time [sec]

FIG. 19: Controller parameters in the long-term test

1.2
System pressure
Reference

Pressure (normalized)

1

0.8

0.6

0.4

0.2
0

20

40

60

80

100

120

140

160

180

200

Time (sec)

FIG. 20: Pressure tracking in the long-term test

learning rates by altering the reference model structure with tracking error feedback.
DR-CRM adaptive controller is tested using a cold air test setup (CATS) which is utilized
as a test bed for throttleable ducted rocket development. A detailed nonlinear mathematical model of the CATS is presented. Model reference adaptive control (MRAC) and CRM
adaptive control along with a proportional integral (PI) controller are also designed for a
comparative evaluation with the proposed controller. Implementation enhancements are pre34

sented to ensure the robustness of the controllers against non-ideal experimental conditions.
A step by step controller design procedure is also provided.
Performances of all controllers are evaluated in simulations and experiments, and similar
conclusions were drawn for both cases. MRAC clearly demonstrates its advantage over constant gain PI controller. More demanding conditions reveal the advantage of DR-CRM adaptive controller over CRM adaptive controller and MRAC. Effect of the projection algorithm
to ensure the boundedness of the controller parameters are also shown in the experiments.
The results of this research showed that combining different control approaches is effective
in addressing complex control problems such as gas generator pressure control. As a future
work, the authors plan to continue this research with experimental studies on real gas
generators and continue to improve/modify the proposed approach according to possibly new
implementation requirements. These requirements may include time-varying or uncertain
time delays, the need for discrete domain design due to non-ideal sampling or different kind
of disturbances that require more sophisticated disturbance rejection methods.

35

Appendix: Memory requirement and computational load of the implementation of
the DR-CRM adaptive controller

Input time delay is observed to be 300 ms on the average, whereas sampling interval of
the controller cycle is chosen as 50 ms, which adds 6 controller parameters, λi , and 6 states,
u(t − mdt), to the controller structure (see (63)). There are 9 states overall, 9 controller
parameters, 9 multiplication results, which adds up to 1 controller signal. We need 9 terms
to the define adaptation laws, 4 terms to define reference model (am , bm , ` and ym ), 1 tracking
error term and 9 adaptation rate terms. In addition, we need 13 terms for the projection
algorithm (Θmax , , ||Θ||, f and ∇f ). Overall, we have 64 float variables that needs 256 bytes
of memory space.
9 multiplication and 8 summation operations are needed to define the controller signal.
There exist 18 multiplications in the adaptive law calculations and 9 summation operations
are required to update the controller parameters. In addition, 4 summation and 3 multiplication operations are needed to be employed to form the reference model output and
the tracking error. Furthermore, in the projection algorithm, there are 2 comparison and 1
logical operations, 18 summation and 43 multiplication operations along with a square root
operator. In total, 116 floating point operations are conducted per controller cycle that runs
with a sampling rate of 50 ms, which results in 2320 floating point operations per second
(FLOPS).

[1] Besser, H.-L. and Kurth, G., “METEOR - European Air Dominance Missile Powered by High
Energy Throttleable Ducted Rocket,” in “Proc. RTO-MP-AVT-208,” , 2012, pp. 1–17.
[2] Goldman, C. and Gany, A., “Thrust modulation of ram-rockets by a vortex valve,” in “AIAA,
ASME, SAE, and ASEE, Joint Propulsion Conference and Exhibit,” Lake Buena Vista, FL,
1996.
[3] Miller, W., McClendon, S., and Burkes, W., “Design Approaches for Variable Flow Ducted
Rockets,” in “Proc. AIAA/SAE/ASME 17th Joint Propulsion Conference,” Colorado Springs,
Colorado, 1981.
[4] Davis, C. A. and Gerards, A. B., “VARIABLE THRUST SOLID PROPULSION CONTROL
USING LABVIEW,” in “Proc. AIAA/ASME/SAE/ASEE 39th Joint Propulsion Conference

36

and Exhibit,” Huntsville, Alabama, AIAA 2003-5241, 2003.
[5] Chang, J., Li, B., Bao, W., Niu, W., and Yu, D., “Thrust control system design of ducted
rockets,” Acta Astronautica, Vol. 69, No. 1, 2011, pp. 86–95.
[6] Burroughs, S., “Status of Army Pintle Technology for Controllable Thrust Propulsion,” in
“37th AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit,” Salt Lake City,
Utah, 2001.
[7] Ostrander, M. and Thomas, M., “Air Turbo-Rocket solid propellant development and testing,”
in “Proc. AIAA/ASME/SAE/ASEE 33th Joint Propulsion Conference and Exhibit,” , 1997.
[8] Bao, W., Li, B., Chang, J., Niu, W., and Yu, D., “Switching control of thrust regulation and
inlet buzz protection for ducted rocket,” Acta Astronautica, Vol. 67, No. 7, 2010, pp. 764–773.
[9] Bao, W., Qi, Y., and Chang, J., “Multi-objective regulating and protecting control for ducted
rocket using a bumpless transfer scheme,” Proceedings of the IMechE, Part G:Journal of
Aerospace Engineering, pp. 311–325.
[10] Qi, Y., Bao, W., Zhao, J., and Chang, J., “Coordinated control for regulation/protection
mode-switching of ducted rockets,” Acta Astronautica, Vol. 98, 2014, pp. 138–146.
[11] Thomaier, D., “Speed control of a missile with throttleable ducted rocket propulsion,” in
“Proc. Advances in Air-Launched Weapon Guidance and Control 15 p (SEE N88-19553 1215),” , 1987.
[12] Sreeriatha, A. G. and Bhardwaj, N., “Mach Number Control-ler for a Flight Vehicle with
Ramjet Propulsion,” in “Proc. AIAA/ASME/SAE/ASEE 35th Joint Propulsion Conference
and Exhibit,” Los Angeles, CA, AIAA 99-294, 1999.
[13] Pinto, P. and Kurth, G., “Robust Propulsion Control in all Flight Stages of a Throtteable
Ducted Rocket,” in “Proc. AIAA/ASME/ASEE Joint Propulsion Conference & Exhibit,” San
Diego, California, AIAA-2011-5611, 2011, pp. 1–12.
[14] Bao, W., Niu, W., J.T., C., Cui, T., and Yu, D., “Control system design and experiment
of needle-type gas regulating system for ducted rocket,” Proceedings of the Institution of
Mechanical Engineers, Part G, Journal of Aerospace Engineering, Vol. 224, No. 5, 2010, pp.
563–573.
[15] Niu, W. Y., Bao, W., Chang, J., Cui, T., , and Yu, D. R., “Control system design and
experiment of needle-type gas regulating system for ducted rocket,” Journal of Aerospace
Engineering, Vol. 224, No. 5, 2010, pp. 563–573.

37

[16] Bauer, C., Hopfe, N., Caldas-Pinto, P., Davenne, F., and Kurth, G., “Advanced Flight Performance Evaluation Methods of Supersonic Air- Breathing Propulsion System by a Highly
Integrated Model Based Approach,” in “Proc. RTO-MP-AVT-208,” , 2012, pp. 1–14.
[17] Ilchmann, A., Ryan, E. P., and Sangwin, C. J., “Tracking with prescribed transient behaviour,” ESAIM: Control, Optimisation and Calculus of Variations, Vol. 7, 2002, pp. 471–
493.
[18] Yildiz, Y., Annaswamy, A., Kolmanovsky, I., and Yanakiev, D., “Adaptive Posicast Controller
for Time-delay Systems with Relative Degree n∗ ≤ 2,” Automatica, Vol. 46, No. 2, 2010, pp.
279–289.
[19] T.Gibson, A.Annaswamy, and E.Lavretsky, “Adaptive systems with closed-loop reference
models: Stability, robustness, and transient performance,” arXiv:1201.4897. Ithaca, NY: Cornell University Library. Retrieved from http://arxiv.org/abs/1201.4897.
[20] T.Gibson, A.Annaswamy, and E.Lavretsky, “Improved Transient Response in Adaptive Control Using Projection Algorithms and Closed Loop Reference Models,” in “AIAA Guidance,
Navigation, and Control Conference,” Minneapolis, Minnesota, AIAA-2012-4775, 2012.
[21] T.Gibson, A.Annaswamy, and E.Lavretsky, “On Adaptive Control with Closed-Loop Reference Models: Transients, Oscillations, and Peaking,” IEEE Access, Vol. 1, 2013, pp. 709–717.
[22] Smith, O. J., “A controller to overcome dead time,” ISA Journal, Vol. 6.
[23] Manitius, A. Z. and Olbrot, A. W., “Finite spectrum assignement problem for systems with
delays,” IEEE Transactions on Automatic Control, Vol. 24, No. 4.
[24] Ichikawa, K., “Frequency-domain pole assignement and exact model-matching for delay systems,” International Journal of Control, Vol. 41, 1985, pp. 1015–1024.
[25] Ortega, R. and Lozano, R., “Globally stable adaptive controller for systems with delay,”
International Journal of Control, Vol. 47, No. 1, 1988, pp. 17–23.
[26] Niculescu, S.-I. and Annaswamy, A. M., “An adaptive Smith-controller for time-delay systems
with relative degree n∗ ≤ 2,” Systems and Control Letters, Vol. 49, 2003, pp. 347–358.
[27] Yildiz, Y., Annaswamy, A., Yanakiev, D., and Kolmanovsky, I., “Spark Ignition Engine Idle
Speed Control: An Adaptive Control Approach,” IEEE Transactions On Control Systems
Technology, Vol. 19, No. 5, 2011, pp. 990–1002.
[28] Yildiz, Y., Annaswamy, A., Yanakiev, D., and Kolmanovsky, I., “Spark Ignition Engine Fuelto-Air Ratio Control: An Adaptive Control Approach,” Control Engineering Practice, Vol. 18,

38

No. 12, 2010, pp. 1369–1378.
[29] Yildiz, Y., Annaswamy, A., Yanakiev, D., and Kolmanovsky, I., “Adaptive Idle Speed Control
for Internal Combustion Engines,” in “Proc. Amer. Control Conf.”, New York City, 2007, pp.
3700–3705.
[30] Yildiz, Y., Annaswamy, A., Yanakiev, D., and Kolmanovsky, I., “Automotive Powertrain
Control Problems Involving Time Delay: An Adaptive Control Approach,” in “Proc. ASME
Dynamic Systems and Control Conference,” Ann Arbor, Michigan, 2008.
[31] Yildiz, Y., Annaswamy, A., Yanakiev, D., and Kolmanovsky, I., “Adaptive Air Fuel Ratio
Control for Internal Combustion Engines,” in “Proc. Amer. Control Conf.”, Seattle, Washington, 2008, pp. 2058–2063.
[32] Z.T.Dydek, A.Annaswamy, J.J.E.Slotine, and E.Lavretsky, “Composite adaptive posicast control for a class of LTI plants with known delay,” Automatica, Vol. 49, No. 6, 2013, pp. 1914–
1924.
[33] Bresch-Pietri, D. and Krstic, M., “Adaptive trajectory tracking despite unknown input delay
and plant parameters,” Automatica, Vol. 45, No. 9, 2009, pp. 2074–2081.
[34] Bekiaris-Liberis, N. and Krstic, M., “Delay-adaptive feedback for linear feedforward systems,”
System and Control Letters, Vol. 59, No. 5, 2010, pp. 277–283.
[35] Krstic, M., Delay Compensation for Nonlinear, Adaptive, and PDE Systems, Birkhauser,
Boston, 2009.
[36] T.Gibson, A.Annaswamy, and E.Lavretsky, “Closed-Loop Reference Models for OutputFeedback Adaptive Systems,” in “European Control Conference (ECC),” IEEE, Zurich, 2013,
pp. 365–370.
[37] T.Gibson, A.Annaswamy, and E.Lavretsky, “Adaptive Systems with Closedloop Reference
Models, Part I: Transient Performance,” in “American Control Conference,” IEEE, Piscataway, NJ, 2013, pp. 3376–3383.
[38] Lavretsky, E., “Adaptive output feedback design using asymptotic properties of lqg/ltr controllers,” IEEE Trans. Automat. Contr., Vol. 57, No. 6.
[39] Stepanyan, V. and Krishnakumar, K., “Mrac revisited: guaranteed perforamance with reference model modification,” in “American Control Conference,” , 2010.
[40] Stepanyan, V. and Krishnakumar, K., “Mmrac for nonlinear systems with bounded disturbances,” in “Conference on Decision and Control,” , 2011.

39

[41] Yucelen, T., Torre, G. D. L., and Johnson, E. N., “Improving transient performance of adaptive
control architectures using frequency-limited system error dynamics,” International Journal
of Control, Vol. DOI:10.1080/00207179.2014.922702.
[42] Pomet, J.-B. and Praly, L., “Adaptive Nonlinear Regulation: Estimation from the Lyapunov
Equation,” IEEE TRANSACTIONS ON AUTOMATIC CONTROL, Vol. 37, No. 6, 1992, pp.
729–740.
[43] Lee, J. H., Park, B. H., and Yoon, W., “Parametric investigation of the pintle-perturbed
conical nozzle flows,” Aerospace Science and Technology, Vol. 26, No. 5, 2013, pp. 268–279.
[44] Deng, R., Setoguchi, T., and Kim, H. D., “Computational Study on the Thrust Performance
of a Supersonic Pintle Nozzle,” in “International Symposium of Turbulence and Shear Flow
Phenomena,” Melbourne, Australia, 2015.
[45] Heo, J., Jeong, K., and Sung, H.-G., “Numerical Study of the Dynamic Characteristics of
Pintle Nozzles for Variable Thrust,” JOURNAL OF PROPULSION AND POWER, Vol. 31,
No. 1, 2015, pp. 230–237.
[46] Ko, H., Lee, J.-H., Chang, H.-B., and Yoon, W.-S., “Cold Tests and the Dynamic Characteristics of the Pintle Type Solid Rocket Motor,” in “49th AIAA/ASME/SAE/ASEE Joint
Propulsion Conference,” San Jose, CA, 2013.
[47] Lee, J. H., Kim, J. K., Jang, H. B., and Oh, J. Y., “Experimental and Theoretical Investigations of Thrust Variation with Pintle Positions Using Cold Gas,” in “44th
AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit,” Hartford, CT, 2008.
[48] Verma, S. B. and Haidn, O., “Cold Gas Testing of Thrust-Optimized Parabolic Nozzle in
a High-Altitude Test Facility,” Journal of Propulsion and Power, Vol. 27, No. 6, 2011, pp.
1238–1246.
[49] Peterson, Z. W., Eilers, S. D., and Whitmorey, S. A., “Closed-loop Thrust and Pressure Profile Throttling of a Nitrous-Oxide HTPB Hybrid Rocket Motor,” in “48th
AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit,” Atlanta, Georgia, 2012.
[50] Narendra, K. S. and Annaswamy, A. M., Stable adaptive systems, Dover Publications, New
York, 2005.
[51] Lavretsky, E. and Wise, K. A., Robust Adaptive Control, Springer, London, 2013.
[52] Yildiz, Y., Annaswamy, A. M., Yanakiev, D., and Kolmanovsky, I., “Spark-Ignition-Engine
Idle Speed Control: An Adaptive Control Approach,” IEEE TRANSACTIONS ON CON-

40

TROL SYSTEMS TECHNOLOGY, Vol. 19, No. 5, 2011, pp. 990–1002.
[53] Gibson, T. E., Closed-Loop Reference Model Adaptive Control: with Application to Very Flexible Aircraft, Phd thesis, Massachusetts Institute of Technology, 2014.

41

