Constrained Model Predictive Control using Feedforward Steady-State
Target Optimization tracking approach for an Isolated Power System
Paulo R. Loma Marconi
1
Abstract— A constrained Linear Quadratic Model Predictive
Model Control (LQ-MPC) with offset-free tracking and dis-
turbance rejection using the Feedforward Steady-State Target
Optimization (SSTO) approach is presented in this paper.
The MPC control law uses Finite-Receding-Horizon method
in dual-mode. Moreover, the stability and recursive feasibility
are guaranteed by off-line dual-mode and on-line optimization.
The results are presented by simulation.
I. STATE-SPACE MODEL OF THE SYSTEM
The isolated power system [1] to be studied is a linear
Single-Input Single-Output (SISO) model with the addition
of the disturbance in the output, Fig. 1.
Fig. 1: Isolated power system
The linear continuous-time state-space model is as follows,
∆ ˙ω
∆ ˙p
m
∆ ˙p
v
=
−D/M 1/M 0
0 −1/T
t
−1/T
t
−1/(RT
g
) 0 −1/T
g
∆ω
∆p
∆p
v
+
0
0
1/T
g
∆p
ref
+
−1/M
0
0
∆p
d
∆f =
50 0 0
∆ω
∆p
∆p
v
+ ∆p
d
(1)
where ∆f is the output values, ∆p
ref
is the input, and ∆p
d
is the disturbance, the rest of the parameters can be consulted
at [1].
Thus, the Discrete Linear Time-Invariant (LTI) state-space
system is
x(k + 1) = A x(k) + B u(k) + B
d
d(k)
y(k) = C x(k) + D
d
d(k), k = 0, 1, 2, . . .
(2)
where x ∈ R
n
, u ∈ R
m
, y ∈ R
p
, and d ∈ R
q
are the states,
input, output, and disturbance vectors respectively. A, B, C,
are the state, input, and output matrices, and B
d
, D
d
are the
disturbance matrices in the input, and in the output.
1
The author is with the Department of Automatic Con-
trol Systems Engineering, The University of Sheffield, UK
prlomamarconi1@sheffiel.ac.uk
II. FORMULATION OF THE LQ-MPC TRACKING
PROBLEM
A. Feedforward tracking with SSTO
The Feedfoward tracking approach by SSTO requires a
target optimizer that calculates the steady-state values for
the states x
ss
and the input u
ss
given a reference r in the
presence of a disturbance d, Fig. 2.
Fig. 2: MPC and SSTO model
B. Cost function
Because it is desired that x → x
ss
and u → u
ss
, new
variables can be defined as ξ ≡ x − x
ss
and v ≡ u − u
ss
,
where ξ → 0 and v → 0.
Therefore, the constrained LQ-MPC problem is to mini-
mize the cost function
IP
N
(v(k)) : min
v (k)
N−1
X
j=0
(ξ
|
(k + j|k) Q ξ(k + j|k)
+ v
|
(k + j|k) R v(k + j|k))
+ ξ
|
(k + N|k) P ξ(k + N|k)
(3)
subject to, for j = 0, 1, 2, . . . , N − 1
ξ(k + 1 + j|k) = A ξ(k + j|k) + B v(k + j|k)
ξ(k|k) = ξ(k)
P
x
ξ(k + j|k) ≤ q
x
− P
x
x
ss
(4)
P
u
v(k + j|k) ≤ q
u
− P
u
u
ss
(5)
P
x
N
ξ(k + N|k) ≤ ˜q
x
N
(6)
where Q and R are weight matrices, P is the terminal cost
matrix, (4), (5), (6) are the linear inequality (polyhedra)
constraints, and N is the finite-horizon value. Note that
B
d
, D
d
and d are not present because the SSTO will manage
the disturbance rejection.
C. Defining the inequality constraints
Given |∆p
ref
| = |u| ≤ u
min,max
, and |∆f| = |y| ≤
y
min,max
, for y = C x,
C
−C
x(k + j|k) ≤
y
max
−y
min
I
−I
u(k + j|k) ≤
u
max
−u
min
can be written in terms of ξ and v as follows,
C
−C
| {z }
P
x
ξ(k + j|k) ≤
y
max
−y
min
| {z }
q
x
−
C
−C
| {z }
P
x
x
ss
I
−I
|{z}
P
u
v(k + j|k) ≤
u
max
−u
min
| {z }
q
u
−
I
−I
|{z}
P
u
u
ss
(7)
and using a proper deadbeat mode-2 terminal inequality
constraint to ensure a fast convergence to zero,
I
n×n
⊗
P
x
P
u
K
∞
(A + BK
∞
)
0
.
.
.
(A + BK
∞
)
N−1
| {z }
P
x
N
ξ(k + N|k)
≤ 1
n
⊗
q
x
− P
x
x
ss
q
u
− P
u
u
ss
| {z }
˜q
x
N
(8)
where K
∞
is the deadbeat mode-2 gain calculated using
Linear Quadratic Regulator (LQR) method.
D. Target equilibrium optimization under constraints (SSTO)
The offset-free equilibrium points (x
ss
, u
ss
) that satisfies
the constraints exits if only if x
ss
= Ax
ss
+ Bu
ss
+ B
d
d
and y
ss
= Cx
ss
+ D
d
d = r.
Thus, the SSTO problem is
min
{x
ss
,u
ss
}
kC x
ss
+ D
d
d − rk + ρku
ss
k (9)
subject to,
I − A −B
C 0
| {z }
A
eq
=T
x
ss
u
ss
|{z}
x
eq
=
B
d
r − D
d
d
| {z }
b
eq
P
x
0
0 P
u
| {z }
A
neq
x
ss
u
ss
|{z}
x
neq
≤
q
x
q
u
|{z}
b
neq
(10)
III. SOLVING THE LQ-MPC OPTIMIZATION
PROBLEM
The solution of the LQ-MPC problem (3) in the compress
form can written as,
min
v (k)
1
2
v
|
(k) H v(k) + ξ
|
(k) L
|
v(k) + ξ
|
(k) M ξ(k)
(11)
subject to
P
c
v(k) ≤ q
c
+ S
c
ξ(k) (12)
where H, L, M are the cost matrices given Q and
R. P
c
, q
c
, S
c
are the stacked inequality constraints for
P
x
, P
u
, P
x
N
, q
u
, q
x
, ˜q
x
N
given the predictions matrices x =
F x(k) + Gu(k), [1]. Note that the prediction matrices are
constructed using the states and input without the tracking
model.
The solution of (11) is the optimal control input sequence,
v
∗
= {v
∗
(k|k), v
∗
(k + 1|k), . . . , v
∗
(k + N − 1|k)}
which can be solved numerically.
Because MPC uses the receding-horizon principle, only
the first optimized control input is needed
v(k) = v
∗
(k|k) = K
N
(ξ(k))
but the real control input applied to the system is
u(k) = u
ss
+ v
∗
(k|k) (13)
A. Stability
In order to ensure stability it is necessary to compute a
terminal cost matrix P that satisfies the Lyapunov equation
(A + B K)
|
P (A + B K) − P + (Q + K
|
R K) = 0 (14)
which solution can be computed numerically. A good ap-
proach is to use a deadbeat mode-2 control gain K that
converges the states to the origin x = 0 as long as (A+BK)
is stable, which means that the eigenvalues λ are inside the
unit circle |λ| < 1.
IV. ALGORITHM FOR THE IMPLEMENTATION
Now that the solution of the LQ-MPC is established and
the stability can be guaranteed, there are basic conditions to
meet such as,
• (A, B) is stabilizable which can be proved by the
controllability, in Matlab rank(ctrb(A,B)).
• Q 0, can be C
|
C
• R 0, changed later for tuning.
The following steps reflect the Matlab implementation
algorithm
Off-line mode:
1. Discretization of the model (1) using
zero holder (zoh), in Matlab will be
sys=c2d(idss(Ac,Bc,C,D,Bd,x0,0),Ts),
where T
s
is the sampling time.
2. Construct the target equilibrium constraints (10)
3. Define the optimization equation as z=[xss;uss]
and fun=@(z) C(1)
*
z(1)+Dd
*
d-r+rho
*
z(4)
4. Solve the optimization problem for SSTO using
f=fmincon(fun,z0,Aneq,bneq,Aeq,beq) in
Matlab, where z0=[r,0,0,0], fun is (9), and
f=[xss,uss] is the optimized vector solution.
5. Construct the LQ-MPC inequality constraints (11), and
compute the deadbeat mode-2 terminal constraint using
K=-dlqr(A,B,Q,R).
6. Construct the prediction matrices F, G using
[F,G]=predict_mats(A,B,N) and the
constraint matrices (12) using [Pc,qc,Sc]=...
constraint_mats(F,G,Pu,qu,Px,qx,PxN,qxN);
7. Compute the deadbeat mode-2 control Kl using
Kl=-acker(A,B,[0,0,0])
8. Calculate the Lyapunov solution of (14) with
P=dlyap((A+B
*
K)’,Q+K’
*
R
*
K), which gives a
terminal cost matrix P that guaranties stability.
9. Finally before the on-line mode is initiated, set the
initial conditions as x=[0;0;0;]
On-line mode: iteration for k steps
1. Recall steps 3-4 of the off-line mode for the SSTO
optimization but at each iteration z0=[x;us(1,k)],
where us(1,k) is the first optimized control input.
2. Calculate ξ(k) and v(k) using xi(:,k)=x-xss
v(:,k)=us(:,k)-uss.
3. Compute the constraints (10) with the new x
ss
and u
ss
4. Compute (12) like it was done in the off-line mode
step 6.
5. Store the states xs(:,k)=x for the iteration process.
6. Solve the LQ-MPC optimization problem
(11) without using the matrix M as follows,
[v,fval,flag]=quadprog(H,L
*
xi(:,k),...
Pc,qc+Sc
*
xi(:,k)) where v is the optimal control
input v(k)
7. Check feasibility using flag, if it is < 1 the solution
is infeasible.
8. Use (13) as the real control input applied to the sys-
tem, us(:,k)=v(1)+uss using only the first
optimized value (receding principle).
9. Close the loop (2), x=A
*
x+B
*
us(1,k)+Bd
*
d
y=C
*
xs+Dd
*
d
10. Calculate the cost value of (11) using the matrix M .
11. Return to step 1.
V. DESIGN REQUIREMENTS
For the simulation the following specifications are de-
manded,
• Reference r = 0.3
• Frequency settles to |∆f | = |y| ≤ 0.01 within 2 [s] of
the disturbance initiaion.
• |∆p
ref
| = |u| ≤ 0.5, and |∆f| = |y| ≤ 0.5.
• Zero steady-state error in the output.
VI. SIMULATION AND RESULTS
The tuning and design parameters, as well as the simula-
tion value results are shown in TableI
Fig. 3a shows that the input constraint is satisfied, and in
Fig. 3b it can be seen that the zero offset-free tracking is also
satisfied rejecting disturbances at different times. The low
value of R is to minimize the overshoot as much as possible
which is 0.6%, the best Q is C
|
C where the settling time
t
s
is around 0.7. Also, the selection of N = 3 was lowest
possible to increase performance, Fig. 3c.
Parameter Value
Q C
|
C
R 0.08
N 3
d
1
0.5 at t = 30[s]
d
2
−0.2 at t = 60[s]
r 0.3
x
0
[0; 0; 0]
|
t
s
0.7[s]
overshoot 0.6%
TABLE I: Tuning and design parameters
(a) u(k)
(b) y(k) vs r(k)
(c) y(k) vs r(k)
Fig. 3: Input and output
Fig. 4 shows that the constraints are satisfied for the three
states. And the cost value shown in Fig. 5 is relatively low
due the low value of R.
(a) States vs step (k)
(b) x1(k) vs x2(k)
(c) x1(k) vs x3(k)
(d) x2 (k) vs x3(k)
Fig. 4: States and phase plots
Fig. 5: Cost value J(k) function
One of the drawbacks of the SSTO approach is that it does
not have robustness because it needs a perfect model of the
plant. Also, when d > 0.5, the solution is infeasible, in order
to handle bigger disturbances it is necessary to change the
cost function, [2] shows an alternative solution.
REFERENCES
[1] P. Trodden, Lecture Notes ACS616 2018/19
[2] S. Dughman, J.A. Rossiter, A survey of guaranteeing feasibility and
stability in MPC during target changes 9th International Symposium
on Advanced Control of Chemical Processes, June 2015.