Application of the Schur Complement in Sum of Squares Optimisation
Elias August
1 a
, Sigurdur Hafstein
2 b
, Jacopo Piccini
1 c
, Stefania Andersen
2 d
,
and Anna Bavarsad
1 e,∗
1
Reykjavik University, Department of Engineering, Menntavegur 1, 102 Reykjavik, Iceland
2
University of Iceland, Faculty of Physical Sciences, Dunhagi 5, 107 Reykjavik, Iceland
Keywords:
Schur Complement, Sum of Squared Polynomials, Stochastic Differential Equation, Lyapunov Function,
Gain Matrix, Numerical Method.
Abstract:
In this paper, we use the Schur Complement in combination with the sum of squares decomposition, first, to
determine whether a nonlinear stochastic dynamical systems has a stable equilibrium and, second, to find a
stabilising gain matrix for nonlinear dynamical systems. In both cases, we consider systems whose dynamics
can be described using polynomial vector fields. Using many different examples, we highlight the effectivity
of using our approaches. In some cases, we manage to obtain results that surpass previous ones. We believe
that the presented approaches have many potential applications, for example, in the fields of aerospace and
quantum control.
1 INTRODUCTION
The Schur complement, a fundamental concept in lin-
ear algebra and matrix theory, has found wide-ranging
applications across various fields of science and en-
gineering. Named after Issai Schur, who introduced
it in the early 20th century, the Schur complement
provides a powerful tool for analysing and simplify-
ing complex matrix inequalities (Carlson et al., 1986).
Pablo Parrilo connected the question whether a poly-
nomial consists of a sum of squares (SOS) to mod-
ern optimisation via linear matrix inequalities and
semidefinite programming (Parrilo, 2003). For dy-
namical systems consisting of polynomials of any de-
gree, this strengthens the requirement of positivity to
the condition that the polynomial function is a SOS.
Crucially, to solve linear matrix inequalities, the re-
sults by Schur, are often applied (Boyd and Vanden-
berghe, 2004). The results by Schur are now a century
old and those by Pablo Parillo just over two decades.
a
https://orcid.org/0000-0001-9018-5624
b
https://orcid.org/0000-0003-0073-2765
c
https://orcid.org/0000-0002-4180-8140
d
https://orcid.org/0000-0001-6747-775X
e
https://orcid.org/0000-0002-6530-2689
∗
This work was supported in part by the Icelandic
Research Fund under Grant 228725-051 and has received
funding from European Union’s Horizon 2020 research and
innovation programme under grant agreement no. 965417.
Nevertheless, we continue to find novel applications
of those to many problems in science and engineer-
ing. In this paper, we combine a novel use of the
Schur complement and SOS decomposition methods
to obtain global asymptotic stability (GAS) certifi-
cates for stochastic dynamical systems whose dynam-
ics can be modelled through polynomial functions,
and to design controllers guaranteeing GAS for poly-
nomial nonlinear dynamical systems.
Natural as well as engineered systems have often
complex dynamics that are affected by noise and dis-
turbances. Additionally, they are difficult, if not im-
possible, to model perfectly, which leads to parameter
uncertainty and approximations in the mathematical
modelling. For biological systems, at times, noise car-
ries information about the underlying network (Mun-
sky et al., 2009; August, 2012). For this reason, one
often requires models based on stochastic differen-
tial equations (SDE). Moreover, in many engineering
fields, system modelling requires the use of nonlinear
ordinary differential equations (Lavretsky and Wise,
2013; Slotine and Li, 1991), which often makes the
design of a satisfactory feedback control law chal-
lenging (Iqbal et al., 2017).
1.1 Stochastic Differential Equations
One often requires models that are described by SDE,
for path planning, modelling aerodynamic forces act-
424
August, E., Hafstein, S., Piccini, J., Andersen, S. and Bavarsad, A.
Application of the Schur Complement in Sum of Squares Optimisation.
DOI: 10.5220/0013016200003822
Paper published under CC license (CC BY-NC-ND 4.0)
In Proceedings of the 21st International Conference on Informatics in Control, Automation and Robotics (ICINCO 2024) - Volume 1, pages 424-431
ISBN: 978-989-758-717-7; ISSN: 2184-2809
Proceedings Copyright © 2024 by SCITEPRESS – Science and Technology Publications, Lda.

ing on aerial vehicles, or the design of ascent phase
control for reusable launch vehicles that involves at-
titude manoeuvring through a wide range of flight
conditions, including stochastic disturbances or non-
deterministic disturbances, engine failure, and aero-
surface locks (Berning Jr, 2020; Rodnishchev and So-
mov, 2018; Xu and Xin, 2011). Thus, modelling us-
ing SDE is important in aerospace science and engi-
neering for the design of control strategies – for exam-
ple, the interpolation between points of operation in
gain scheduling can become otherwise prohibitively
costly (to compute and store) – and path/mission plan-
ning, where it can provide safety certificates, to name
but two.
Another example is finding feedback control laws
that stabilise a pure state of a quantum system by
means of continuous quantum non-demolition mea-
surements. The preparation of such a state is crucial
for quantum technologies. Quantum non-demolition
measurements can be seen as “classical” ones (they
do not prevent the wave function from collapsing).
However, they do not destroy the particle under ob-
servation: In the double-slit experiment, information
about a photon’s trajectory is gained without captur-
ing/destroying the photon. Significantly, analysis and
control of a stochastic dynamical system is a more
challenging problem compared with its determinis-
tic counterpart, particularly, for nonlinear stochas-
tic dynamical systems (Ludyk, 2018; Cardona et al.,
2018; Wiseman, 1994). Crucially, control strategies
for quantum system pure state stabilisation are still in
their infancy.
1.2 Nonlinear Control
Many effective techniques exist to design a controller
for a nonlinear dynamical system. However, they typ-
ically ensure asymptotic stability only near the op-
erating point of interest (local asymptotic stability).
One common method is to linearise the set of dif-
ferential equations either through Jacobian lineari-
sation or feedback linearisation. Alternatively, one
can try to solve the State-Dependent Riccati Equa-
tion, which is often very difficult, to derive an op-
timal control law (Cloutier, 1997). For systems de-
scribed by polynomial vector fields, a powerful ap-
proach to controller design was developed in (August
and Papachristodoulou, 2022). This approach, based
on SOS decomposition, ensures asymptotic stability
of the closed loop system.
1.3 Organisation of the Paper
The organisation of the paper is the following. In Sec-
tion 2, we present the different methods used in this
paper. More precisely, in Section 2.1, we present the
Schur Complement and in Section 2.2 the SOS de-
composition. In Section 2.3, we show how the two
can be used to determine the stability of stochastic
dynamical systems and similarly, in Section 2.4 how
they can be used to design a stabilising controller for
a nonlinear dynamical system, whose dynamics are
described by means of polynomial functions. The
main contributions of this paper are the convex prob-
lem relaxations presented in Sections 2.3.1 and 2.4.1
and their application in Section 3. For instance,
we demonstrate the usefulness of the presented ap-
proaches in determining a stabilising controller for a
quantum system as well as an aerospace system. Fi-
nally, we conclude the paper in Section 4.
2 METHODS
In this section, we present the different methods, in-
cluding their novel application, used to obtain the re-
sults presented in Section 3.
2.1 Schur Complement
Thereom 1. Consider,
F =

A B
B
T
C

.
If matrices A and C are symmetric then the following
statements are equivalent:
i F ≻ 0,
ii A ≻ 0 and C −B
T
A
−1
B ≻ 0,
iii C ≻ 0 and A − BC
−1
B
T
≻ 0.
Similarly, the following statements are equivalent:
iv F ≺ 0,
v A ≺ 0 and C −B
T
A
−1
B ≺ 0,
vi C ≺ 0 and A − BC
−1
B
T
≺ 0.
This theorem can be traced back to Schur’s orginal
paper (Schur, 1917) and similar results for semidefi-
nite F exist. To see that (i) ⇔ (ii), note that if M is
nonsingular then M
T
FM ≻ 0 ⇔ F ≻ 0 and that

I −A
−1
B
0 I

T

A B
B
T
C

I −A
−1
B
0 I

=

A 0
0 C −B
T
A
−1
B

.
All other equivalencies can be shown similarly. Ma-
trices C − B
T
A
−1
B and A − BC
−1
B
T
are called the
Schur Complements.
Application of the Schur Complement in Sum of Squares Optimisation
425

2.2 Sum of Squares Decomposition
Testing for non-negativity of real-valued polynomial
function F(x) of degree 2d is NP-hard (Murty and
Kabadi, 1987), x ∈ R
n
. However, a sufficient con-
dition for F(x) to be nonnegative is that it can be de-
composed into SOS (Parrilo, 2003):
F(x) =
∑
i
f
2
i
(x) ≥ 0,
where f
i
are polynomial functions. Now, F(x) is SOS
if and only if there exists a matrix R such that
F(x) =
∑
i
f
2
i
(x) = χ
T
Rχ, R = R
T
⪰ 0, (1)
χ
T
=
h
1 x
(1)
... x
(n)
x
(1)
x
(2)
... x
d
(n)
i
.
The entries of vector χ consist of all monomial com-
binations of the elements of vector x up to degree d
(including x
0
(i)
= 1) and, thus, its length is ℓ =

n+d
d

.
While R is usually non-unique, (1) poses certain con-
straints on it of the form
tr(A
j
R) = c
j
, j = 1, 2,. ..,m,
where tr(M) denotes the trace of square matrix M and
matrices A
j
and constants c
j
are problem dependent.
As illustration consider
F(x) = 2x
4
1
+ 2x
3
1
x
2
− x
2
1
x
2
2
+ 5x
4
2
=


x
2
1
x
2
2
x
1
x
2


T


q
11
q
12
q
13
q
12
q
22
q
23
q
13
q
23
q
33




x
2
1
x
2
2
x
1
x
2


= q
11
x
4
1
+ q
22
x
4
2
+ (q
33
+ 2q
12
)x
2
1
x
2
2
+ 2q
13
x
3
1
x
2
+ 2q
23
x
1
x
3
2
.
For χ
1
= x
2
1
, χ
2
= x
2
2
, and χ
3
= x
1
x
2
, we obtain the
following linear equalities:
q
11
= 2, q
22
= 5, q
33
+2q
12
= −1, 2q
13
= 2, 2q
23
= 0.
Thus, j = 1, .. .,5, and, for example, c
1
= 2, c
3
= −1,
A
1
=


1 0 0
0 0 0
0 0 0


, and A
3
=


0 1 0
1 0 0
0 0 1


.
In order to find R, we solve the optimisation prob-
lem associated with the following semidefinite pro-
gramme, where matrix A
0
can be chosen to select a
particular solution R:
min tr(A
0
R)
s.t. tr(A
j
R) = c
j
, j = 1, .. ., m
R = R
T
⪰ 0. (2)
In this paper, to solve (2), we use SOSTOOLS (Pa-
pachristodoulou et al., 2021) and the SeDuMi
solver (Sturm, 1999).
2.3 Stochastic Differential Equations
Consider the following SDE,
dX = b(X)dt + σ(X)dW, (3)
where X is the random state vector, W is an m-
dimensional Wiener Process, b : R
n
→ R
n
the drift
function and σ : R
n
→ R
n×m
the state-dependent dif-
fusion matrix (Ventsel and Freidlin, 1970; Gihman
and Skorohod, 1972).
1
Moreover, we assume that
both maps, b and σ, satisfy the usual conditions for
the existence and uniqueness of solutions (El-Samad
and Khammash, 2004). For the reminder of this pa-
per, we use the short-form notation σ for σ(X). Let
the origin be an equilibrium point of system (3). The
associated infinitesimal generator is given by
AV (X) =
n
∑
i=1
"
b
i
(X)
∂V (X)
∂X
i
+
1
2
n
∑
j=1
(σσ
T
)
i, j
∂
2
V (X)
∂X
i
X
j
#
for any given function V ∈ C(R
n
) ∩ C
2
(R
n
\ {0}).
If V (0) = 0, AV (0) = 0, and for all X ∈ R
n
\ {0},
V (X) > 0 and AV (X) < 0 then the expected value of
V (X) decreases with time, where X evolves according
to (3) (Øksendal, 2003).
Now, let
dX = A(X)Xdt + G(X)XdW, (4)
where dW is 1-dimensional and G(x) and A(x) ∈
R[x]
n×n
are state-dependent matrices, and positive
definite matrix Q be such that
V (X) =

X
T
QX

p
2
, p > 0.
Then, in the following, we provide conditions for the
associated AV (X) to be negative definite (Hafstein
et al., 2018), where for clarity we use x instead of X.
First note that
b
i
(x) =
n
∑
j=1
A
i, j
(x)x
j
,
V (x) =
n
∑
i=1
x
i
n
∑
j=1
Q
i, j
x
j
!
p
2
,
∂V (x)
∂x
i
=
p
2

x
T
Qx

p−2
2
2
n
∑
j=1
Q
i, j
x
j
. (5)
1
When noise in (3) is intrinsic, is it also termed vanish-
ing stochastic perturbations, since then σ(0) = 0.
ICINCO 2024 - 21st International Conference on Informatics in Control, Automation and Robotics
426

Then, it follows from (5) that
n
∑
i=1
b
i
(x)
∂V (x)
∂x
i
= p

x
T
Qx

p−2
2
n
∑
i=1
n
∑
k=1
n
∑
j=1
x
k
Q
k,i
A
i, j
(x)x
j
= p

x
T
Qx

p−2
2
x
T
QA(x)x
= p

x
T
Qx

p−4
2
x
T
Qxx
T
QA(x)x
= −
1
2
p

x
T
Qx

p−4
2
H
0
(x).
Now, for q = 2 − p,
∂
2
V (x)
∂x
i
x
j
= p
∂
∂x
j

x
T
Qx

p−2
2
n
∑
k=1
Q
i,k
x
k
= −qp

x
T
Qx

p−4
2
n
∑
k=1
Q
i,k
x
k
n
∑
ℓ=1
Q
j,ℓ
x
ℓ
+ px
T
Qx

x
T
Qx

p−4
2
Q
i, j
.
Next, consider, with y := G(x)x,

G(x)xx
T
G(x)
T

i, j
=

yy
T

i, j
= y
i
y
j
.
Then,
1
2
n
∑
i=1
n
∑
j=1
(σσ
T
)
i, j
∂
2
V (x)
∂x
i
x
j
= −
1
2
p

x
T
Qx

p−4
2
(H
1
(x) + H
2
(x)),
where
H
1
(x) = q
n
∑
k=1
n
∑
i=1
x
k
Q
k,i
y
i
n
∑
ℓ=1
n
∑
j=1
x
ℓ
Q
ℓ,i
y
j
= q

x
T
QG(x)x

2
and
H
2
(x) = −x
T
Qx
n
∑
i=1
n
∑
j=1
y
i
Q
i, j
y
j
= −x
T
Qxy
T
Qy.
Thus, AV (x) < 0 if H(x) > 0 for all x ∈ R
n
, x ̸= 0,
where
H(x) = H
0
(x) + H
1
(x) + H
2
(x)
= −x
T

2QA(x) + G(x)
T
QG(x)

xx
T
Qx
+q

x
T
QG(x)x

2
. (6)
2.3.1 Problem Relaxation
Note that, given c ∈ R, for unknown Q, unless p = 2,
H(x) ≥ c

x
T
x

2
is a bilinear matrix inequality, which
is a non-convex problem. To make the problem con-
vex, we first require that Q ⪯ ¯cI
n
, ¯c ∈ R, and then
make use of the following inequality

x
T
QG(x)x

2
x
T
Qx
≥

x
T
QG(x)x

2
¯c
,
which implies that H(x) ≥
−x
T

2QA(x) + G(x)
T
QG(x)

xx
T
x ¯c+q

x
T
QG(x)x

2
.
Next, we replace

x
T
QG(x)x

2
by z
T
Xz, where z con-
tains only those monomials that are necessary for
equality
z
T
Xz =

x
T
QG(x)x

2
, (7)
to hold (if such a matrix X exists), and require that
z
T
Xz ≥

x
T
QG(x)x

2
, which can be written as con-
vex constraint using the Schur Complement (see (8)).
Thus, if state-dependent matrix G(x) has only polyno-
mial entries then to find a feasible solution of (6) such
that H(x) ≥ c

x
T
x

2
, we solve the following SOS op-
timisation problem, where we seek to enforce (7) by
minimising the trace of matrix X:
given A(x),G(x) ∈ R[x]
n×n
, ¯c ∈ R, c ∈ R, q
min. tr(X)
s. t.
¯
H(x) − c

x
T
x

2
is SOS ∀x ∈ R
n
¯
H(x) = −x
T
(2QA(x) + G(x)
T
QG(x))xx
T
x ¯c
+qz
T
Xz
Q ⪰ I
n
, Q ⪯ ¯cI
n
, X ⪰ 0
y
T
M(x)y is SOS ∀x ∈ R
n
& ∀y ∈ R
2
M(x) =

z
T
Xz x
T
QG(x)x
x
T
QG(x)x 1

. (8)
If the solution of (8) is such that (7) holds for matrix
X then H(x) ≥
¯
H(x) and the expected value of (4)
decreases with time.
Remark 1. Note that if z contains not only those
monomials that are necessary for equality z
T
Xz =

x
T
QG(x)x

2
to hold then there might exist a ma-
trix
˜
X, for which tr(
˜
X) < tr(X) holds, while z
T
˜
Xz >
z
T
Xz =

x
T
QG(x)x

2
also holds, as the following ex-
ample shows.
For instance, consider
Q =

1.300976011928209 0.176426076324787
0.176426076324787 1.103417407153368

and
G =

0 0
0 −2

.
Then, z
T
Xz = (x
T
QGx)
2
, where z
T
=

x
2
1
x
1
x
2

,
X =

0.124504641629441 0.778686414770154
0.778686414770154 4.870119897636248

Application of the Schur Complement in Sum of Squares Optimisation
427

and tr(X) = 4.9946. However, if we allow “unneces-
sary” monomials in z to appear – that is, those that
would appear if matrix G were non-singular – then
z
T
=

x
2
1
x
1
x
2
x
2
2

, and, for example, for
˜
X =



0.005258248861173 0.023272688728787 0.154806237483972
0.023272688728787 0.110381268124464 0.732647367566061
0.154806237483972 0.732647367566061 4.878012769236294



,
where tr(
˜
X) = 4.99365, the following holds,
z
T
˜
Xz > (x
T
QGx)
2
.
2.4 Nonlinear Control
In this paper, we represent nonlinear dynamical sys-
tems, whose dynamics are governed by polynomial
function, using the following notation:
˙x = A(x)x +Bu,
x(t) ∈ R
n
, u(t) ∈ R
m
,
A(x) ∈ R[x]
n×n
, B ∈ R
n×m
,
where x is the system state, u the input, A(x) the state-
dependent polynomial system matrix representing the
system dynamics, B is the input matrix, and we de-
note the set of all polynomials by R[x]. Particularly,
A(x) ∈ R[x]
n×n
means that the following holds for its
(i, j)-th entry, A
(i, j)
∈ R[x], where i, j ∈ {1,2, .. ., n}.
Although matrix B can also be state-dependent, for
simplicity, we consider constant input matrices only.
Now, consider
˙x = A(x)x + BK(x,k)x,
K(x, k) = −B
T
Q f (x,k),
f (x, k) =
k
∑
i=0
(x
T
x)
2k
. (9)
If the origin is an asymptotically stable equilibrium
point of (9) then we seek a positive definite matrix Q
such that we can prove this by means of the Lyapunov
function given by V (x) = x
T
Qx, which implies that
the following inequality must hold for all x , x ̸= 0,
x
T

QA(x)x − QBB
T
Q f (x,k)

x < 0. (10)
2.4.1 Problem Relaxation
Note that (10) is a bilinear matrix inequality and, thus,
a non-convex problem. To make this problem convex,
we replace matrix QBB
T
Q by matrix X and require
that X ⪰ QBB
T
Q, which can be written as convex con-
straint using the Schur Complement. Then, we try to
obtain positive definite matrix Q by solving the fol-
lowing SOS optimisation problem, where we seek to
enforce X = QBB
T
Q by minimising the trace of ma-
trix X:
given A(x) ∈ R[x]
n×n
, B ∈ R
m×n
, c ∈ R, k
minimise tr(X)
subject to x
T
(X f (x, k) − QA(x))x is SOS ∀x ∈ R
n

X QB
B
T
Q I
m

⪰ 0, Q ⪰ cI
n
(11)
If (11) has a feasible solution given by Q and X =
QBB
T
Q holds then matrix Q stabilises system (9).
3 RESULTS
In the following, we demonstrate the use of the ap-
proaches presented in this paper by applying them
to different exemplary systems taken from the liter-
ature (Hafstein et al., 2018; Cardona et al., 2018; Au-
gust and Papachristodoulou, 2022).
3.1 Linear SDE
3.1.1 Example 1: 2D System
Consider
A =

0 2.5
−2.5 −0.9

, G =

0 0
0 −2

.
First, note that, for this system, the method presented
in (Hafstein et al., 2018) was unable to find a matrix
Q, for which H(x) > 0 holds. Now, for p = 0.1, c =
0.3, and ¯c = 1.41, we solve (8) and determine that
H(x) ≥ c

x
T
x

2
for
Q =

1.3010 0.1764
0.1764 1.1034

by means of the following SOSTOOLS programme:
A =[0 2.5 ; -2.5 -.9];
G =[0 0; 0 -2];
p =0 .1;
q =2 - p ;
c =0 .3;
cb = 1 .41;
solver_opt . solve r =' sedu m i' ;
pva r x1 x2 v1 v2
x =[ x1 x2 ]';
v =[ v1 v2 ]';
p = sosprogra m ([ x;v ]) ;
[p , Q ]= s o s p o l y m at r i x v a r ( p , m o n o m i als (
x ,0 ) ,[2 2] ,' symm e t r i c ' ) ;
[p , X ]= s o s p o l y m at r i x v a r ( p , m o n o m i als (
x ,0 ) ,[2 2] ,' symm e t r i c ' ) ;
ICINCO 2024 - 21st International Conference on Informatics in Control, Automation and Robotics
428

[~ ,z ,~ ]= findso s ((x ' * r and (2) *G*x) ˆ2)
p = sosi n e q (p ,x ' *(Q - e ye (2) )*x);
p = sosi n e q (p ,x ' *( cb * eye (2 ) -Q)* x ) ;
H = -x ' *(2* Q* A + G ' *Q*G ) * x *x ' *x * cb...
+ q *z ' *X*z ;
p = sosi n e q (p ,H - c*( x '*x ) ˆ2) ;
p = sosi n e q (p ,x ' *X*x ) ;
p = sosi n e q (p ,v ' *[z '*X * z x'* Q * G *x;...
x '*Q*G * x 1]* v ) ;
p = sosseto b j (p , trac e ( X ) ) ;
[p ,~ ]= sossolve ( p , s o l v e r _ opt );
Q = dou b le ( so s g e t s ol (p , Q ) )
sosgetso l (p , z ' *X * z ) -( x ' *Q * G * x)ˆ2
3.1.2 Example 2: 3D System
Consider
A =


0.2759 0.0831 −0.9603
0.8794 −0.8281 0.2348
−0.3879 −1.8186 −0.1508


,
and
G =


3.3849 −0.3554 0.1794
−0.3554 2.2541 0.3234
0.1794 0.3234 2.8609


.
For p = 0.1, c = 0.1, and ¯c = 1, using a SOSTOOLS
programme similar to the one used in Example 3.1.1
to solve (8), we determine that H(x) ≥ c

x
T
x

2
for
Q = I.
3.2 Density Matrix of a Two-State
Quantum System
A two-state system is a quantum system that consists
of two independent – that is, physically distinguish-
able – quantum states and all their quantum superpo-
sitions. For example, the Stern-Gerlach experiment
is such a system, where the distribution of magnetic
moments is not a continuous one but is limited to two
values. Moreover, there are matrices, so-called ob-
servables, associated with such a system. If A is an
observable then the measurement of a physical entity
corresponding to observable A must be an eigenvalue
λ
i
of A (Ludyk, 2018). Now, for the eigenvalues to be
always real, observable A must be Hermitian and, it
follows from the above, that the state variable can be
described by a linear combination of its eigenvectors.
The density matrix of a two-state quantum system
describes an ensemble of states. Let it be given by
density operator
ρ =
1
2
(I + σ
x
x + σ
y
y + σ
z
z) =
1
2

1 + z x − iy
x + iy 1 − z

,
where σ
x
, σ
y
, and σ
x
are the Pauli matrices given by
σ
x
=

0 1
1 0

, σ
y
=

0 −i
i 0

, σ
z
=

1 0
0 −1

,
and real scalars x, y, and z satisfy
x
2
+ y
2
+ z
2
≤ 1.
We would like to drive ρ towards a pure state, that is,
towards z = 1 or z = −1 by means of non-demolition
measurement and using feedback; for example, by
first entangling the system and then using the (pos-
sibly destructive) measurement of the entangled sys-
tem for feedback. For instance, such feedback can be
used for quantum-state preparation, which is crucial
for quantum technologies.
The quantum closed system evolves as in the fol-
lowing:
˙
ρ = −i[H,ρ] = −i(Hρ − ρH),
where H is a Hermitian matrix and [H,ρ] = Hρ− ρH.
The quantum open (measured) system with feedback
evolves as in the following (Cardona et al., 2018;
Wiseman, 1994):
dρ = −i f [σ
y
,ρ]dt + D(σ
z
− ikσ
y
,ρ)dt
+M(σ
z
,ρ)dW − iκ[σ
y
,ρ]dW, (12)
where κ and f are tuneable control parameters,
tr(ρ) = 1, ρ =

ρ
1
ρ
2
ρ
∗
2
ρ
3

,
D(L, ρ) = LρL
†
−
1
2
L
†
Lρ −
1
2
ρL
†
L,
and
M(L, ρ) = Lρ + ρL
†
− tr((L + L
†
)ρ)ρ.
Now, consider the following change of coordinates,
x = tr(ρσ
x
), z = tr(ρσ
z
) ⇒ x = ρ
2
+ ρ
∗
2
, z = ρ
1
− ρ
3
.
Then,
dx = 2(2κ −xκ
2
− x + f z)dt + 2(κ − x)zdW
and
dz = −2(zκ
2
+ f x)dt + 2(1 −z
2
− κx)dW.
Next, we use also the following coordinate trans-
formation: ¯z = z − 1. For control parameters κ = −¯z
and f = (¯z + 2)¯z, as suggested in (Cardona et al.,
2018), we obtain,
dx = −2

2¯z + x¯z
2
+ x − (¯z + 2)(¯z + 1)¯z

dt
−2(¯z + x)(¯z + 1)dW
= −2

x¯z
2
+ x − ¯z
3
− 3¯z
2

dt
−2

¯z + 1 ¯z + 1

χdW
= −2

¯z
2
+ 1 −¯z
2
− 3¯z

χdt
−2

¯z + 1 ¯z + 1

χdW (13)
Application of the Schur Complement in Sum of Squares Optimisation
429

and
dz = −2

¯z
3
+ ¯z
2
+ ¯z
2
x + 2¯zx

dt
−2(¯z
2
+ 2¯z − ¯zx)dW
= −2

¯z
2
+ 2¯z ¯z
2
+ ¯z

χdt
−2

−¯z ¯z + 2

χdW, (14)
where χ =

x
¯z

. Thus,
A(χ) = −2

¯z
2
+ 1 −¯z
2
− 3¯z
¯z
2
+ 2¯z ¯z
2
+ ¯z

and
G(χ) = −2

¯z + 1 ¯z + 1
−¯z ¯z + 2

.
To check for exponential stability of (13)–(14), for
p = 0.1, c = 0.065, and ¯c = 1, we use a SOSTOOLS
programme similar to the one used in Example 3.1.1
to solve (8) and determine that H ≥ c

x
T
x

2
for Q =
I. Note that, in this section, matrices A(χ) and G(χ)
are state dependent, as opposed to constant matrices A
and G in Section 3.1. Finally, exponential stability of
(13)–(14) has been already shown in (Cardona et al.,
2018) and we do not claim to improve on these results,
we rather aim at exemplifying our approach.
3.3 Nonlinear Control
The system dynamics of the examples in this section
are all given by ˙x = A(x)x + BKx.
3.3.1 Tunnel Diode Circuit
The dynamics of a tunnel diode circuit are defined by
A(x) =

−0.5g(x
1
) 0.5
−0.2 −0.3

, B =

0
0.2

,
g(x
1
) = 17.76 − 103.79x
1
+ 229.62x
2
1
− 226.31x
3
1
+ 83.72x
4
1
.
Here, x
1
is the voltage across the capacitor and x
2
the
current through the inductor. For c = 1 and k = 0,
we solve (11) and obtain stabilising gain matrix K =
−

0 0.2281

through the following SOSTOOLS
programme (that also confirms that X − QBB
T
Q ≈ 0):
solver_opt . solve r =' sedu m i' ;
pva r x1 x2 v
x =[ x1 ; x2 ];
w =[ x ; v];
p = sosprogra m ( w );
g =17 . 7 6 -10 3 . 79* x1 + 2 2 9 . 62* x1 ˆ2 ...
-22 6.31* x1 ˆ3+ 8 3 . 7 2 * x1 ˆ4;
A =[ -0.5 * g 0.5 ; -0.2 -0.3];
B =[0;0 . 2 ] ;
[p , Q ]= s o s p o l y m a tr i x v a r ( p , m o n o m i als (
x ,0) ,[2 2] ,' symm e t r i c ' ) ;
[p , X ]= s o s p o l y m a tr i x v a r ( p , m o n o m i als (
x ,0) ,[2 2] ,' symm e t r i c ' ) ;
p = sosi n e q (p ,x ' *(Q - e ye (2) )*x);
p = sosi n e q (p ,x ' *(X -Q * A ) * x);
p = sosi n e q (p ,w ' *[X Q * B;B '*Q 1]* w) ;
p = sosseto b j (p , trac e ( X ) ) ;
[p ,~ ]= sossolve ( p , s o l v e r _ opt );
X = sosgets o l (p , X);
Q = sosgets o l (p , Q);
X - Q*B*B'* Q
K = -B ' *Q
3.3.2 Air-Breathing Hypersonic Flight Vehicle
The longitudinal dynamics of a simplified version of
an air-breathing hypersonic flight vehicle are defined
by A
11
= 0,
A
12
= −(0.645x
2
+ 0.01921)(0.1247x
1
+ 0.7370)
2
,
A
21
= −0.0002706, A
22
= −0.0574 − 0.009716x
1
,
B =

0.014
0

.
Here, x
1
is the velocity, x
2
the angle of attack, and we
use the fact that for small angles sin(x
2
) ≈ x
2
. Further-
more, u is the thrust input. For c = 1 and k = 2, using
a SOSTOOLS programme similar to the one used in Ex-
ample 3.3.1 to solve (11), we obtain stabilising gain
matrix K = −

0.0140 0


1 + x
T
x +

x
T
x

2

.
3.3.3 Lorenz System
Consider the Lorenz system defined by
A(x) =


10 −10 0
28 −1 −x
(1)
x
(2)
0 8


, B =


1 0
1 0
0 1


.
For c = 11 and k = 0, using a SOSTOOLS programme
similar to the one used in Example 3.3.1 to solve (11),
we obtain stabilising gain matrix K = −11B
T
.
4 CONCLUSIONS
In this paper, we presented two related approaches to
either determine whether a nonlinear stochastic dy-
namical system has a stable equilibrium or to find a
stabilising gain matrix for a nonlinear dynamical sys-
tem, where we considered systems whose dynamics
ICINCO 2024 - 21st International Conference on Informatics in Control, Automation and Robotics
430

can be described using polynomial vector fields. Sig-
nificantly, by using the Schur Complement in com-
bination with the sum of squares decomposition, we
provided convex alternatives to bilinear matrix in-
equalities. Using different examples, we highlighted
the effectivity of using our approaches, which also
managed to obtain results that surpassed previous
ones. We believe that the presented approaches have
many potential applications, for example, in the fields
of aerospace and quantum control.
REFERENCES
August, E. (2012). Using noise for model-testing. Journal
of Computational Biology, 19:968–977.
August, E. and Papachristodoulou, A. (2022). Feedback
control design using sum of squares optimisation. Eu-
ropean Journal of Control, page 100683.
Berning Jr, A. (2020). Control and Optimization for
Aerospace Systems with Stochastic Disturbances, Un-
certainties, and Constraints. PhD thesis, University
of Michigan.
Boyd, S. and Vandenberghe, L. (2004). Convex Optimiza-
tion. Cambridge University Press, Cambridge, UK.
Cardona, G., Sarlette, A., and Rouchon, P. (2018). Expo-
nential stochastic stabilization of a two-level quantum
system via strict lyapunov control. In 2018 IEEE Con-
ference on Decision and Control, pages 6591–6596.
Carlson, D., Markham, T. L., and Uhlig, F. (1986). Em-
ilie haynsworth, 1916-1985. Linear Algebra and its
Applications, 75:269–276.
Cloutier, J. (1997). State-dependent riccati equation tech-
niques: an overview. In Proceedings of the 1997
American Control Conference (Cat. No.97CH36041),
pages 932–936.
El-Samad, H. and Khammash, M. (2004). Stochastic stabil-
ity and its application to the analysis of gene regula-
tory networks. Proceedings of the American Control
Conference, pages 3001–3006.
Gihman, I. I. and Skorohod, A. V. (1972). Stochastic Dif-
ferential Equations. Springer-Verlag, Berlin.
Hafstein, S., Gudmundsson, S., Giesl, P., and Scalas,
E. (2018). Lyapunov function computation for au-
tonomous linear stochastic differential equations us-
ing sum-of-squares programming. Discrete and Con-
tinuous Dynamical Systems - B, 23(2):939–956.
Iqbal, J., Ullah, M., Khan, S. G., Khelifa, B., and
´
Cukovi
´
c,
S. (2017). Nonlinear control systems - a brief
overview of historical and recent advances. Nonlin-
ear Engineering, 6(4):301–312.
Lavretsky, E. and Wise, K. (2013). Robust and Adap-
tive Control: With Aerospace Applications. Springer-
Verlag London, London, United Kingdom.
Ludyk, G. (2018). Quantum Mechanics in Matrix Form.
Springer, Cham, Switzerland.
Munsky, B., Trinh, B., and Khammash, M. (2009). Lis-
tening to the noise: random fluctuations reveal gene
network parameters. Molecular Systems Biology, 5.
Murty, K. G. and Kabadi, S. N. (1987). Some NP-complete
problems in quadratic and nonlinear programming.
Math. Program., 39:117–129.
Øksendal, B. (2003). Stochastic Differential Equation:
An Introduction with Applications. Springer-Verlag,
Berlin, 6 edition.
Papachristodoulou, A., Anderson, J., Valmorbida, G., Pra-
jna, S., Seiler, P., Parrilo, P. A., Peet, M. M., and Jagt,
D. (2021). SOSTOOLS: Sum of squares optimization
toolbox for MATLAB. http://arxiv.org/abs/1310.4716.
https://github.com/oxfordcontrol/SOSTOOLS.
Parrilo, P. A. (2003). Semidefinite programming relaxations
for semialgebraic problems. Math. Program., Ser. B,
96:293–320.
Rodnishchev, N. and Somov, Y. (2018). Control opti-
mization in aerospace engineering at stochastic per-
turbations and stream of faults. In 2018 5th IEEE
International Workshop on Metrology for AeroSpace
(MetroAeroSpace), pages 171–175.
Schur, I. (1917).
¨
Uber Potenzreihen, die im Inneren des Ein-
heitskreises beschr
¨
ankt sind. J. Reine Angew. Math.,
147:205–232.
Slotine, J. J. E. and Li, W. (1991). Applied Nonlinear Con-
trol. Prentice Hall, Englewood Cliffs, NJ, USA.
Sturm, J. F. (1999). Using SeDuMi 1.02, a MATLAB tool-
box for optimization over symmetric cones. Optimiza-
tion Methods and Software, 11–12:625–653. Avail-
able at http://sedumi.ie.lehigh.edu.
Ventsel, A. D. and Freidlin, M. I. (1970). On
Small Random Perturbations of Dynamical Systems.
Russ. Math. Surv., 25:1.
Wiseman, H. M. (1994). Quantum theory of continuous
feedback. Phys. Rev. A, 49:2133–2150.
Xu, Y. and Xin, M. (2011). Nonlinear stochastic con-
trol for space launch vehicles. IEEE Transactions on
Aerospace and Electronic Systems, 47(1):98–108.
Application of the Schur Complement in Sum of Squares Optimisation
431