Ellipsoidal Calculus

Basic Notions

We start with basic definitions.

Definition. Ellipsoid {\mathcal E}(q,Q) in {\bf R}^n with center q and shape matrix Q is the set

(1){\mathcal E}(q,Q) = \{ x \in {\bf R}^n ~|~ \langle (x-q), Q^{-1}(x-q)\rangle\leqslant1 \},

wherein Q is positive definite (Q=Q^T and \langle x, Qx\rangle>0 for all nonzero x\in{\bf R}^n). Here \langle\cdot,\cdot\rangle denotes inner product.

Definition. The support function of a set {\mathcal X}\subseteq{\bf R}^n is

\rho(l~|~{\mathcal X}) = \sup_{x\in{\mathcal X}} \langle l,x\rangle.

In particular, the support function of the ellipsoid (1) is

(2)\rho(l~|~{\mathcal E}(q,Q)) = \langle l, q\rangle + \langle l, Ql\rangle^{1/2}.

Although in (1) Q is assumed to be positive definite, in practice we may deal with situations when Q is singular, that is, with degenerate ellipsoids flat in those directions for which the corresponding eigenvalues are zero. Therefore, it is useful to give an alternative definition of an ellipsoid using the expression (2).

Definition. Ellipsoid {\mathcal E}(q,Q) in {\bf R}^n with center q and shape matrix Q is the set

(3){\mathcal E}(q,Q) = \{ x \in {\bf R}^n ~|~
\langle l,x\rangle\leqslant\langle l,q\rangle + \langle l,Ql\rangle^{1/2}
\mbox{ for all } l\in{\bf R}^n \},

wherein matrix Q is positive semidefinite (Q=Q^T and \langle x, Qx\rangle\geqslant0 for all x\in{\bf R}^n). The volume of ellipsoid {\mathcal E}(q,Q) is

(4){\bf Vol}(E(q,Q)) = {\bf Vol}_{\langle x,x\rangle\leqslant1}\sqrt{\det Q},

where {\bf Vol}_{\langle x,x\rangle\leqslant1} is the volume of the unit ball in {\bf R}^n:

(5){\bf Vol}_{\langle x,x\rangle\leqslant1} = \left\{\begin{array}{ll}
\frac{\pi^{n/2}}{(n/2)!}, &
\mbox{ for even } n,\\
\frac{2^n\pi^{(n-1)/2}\left((n-1)/2\right)!}{n!}, &
\mbox{ for odd } n. \end{array}\right.

The distance from {\mathcal E}(q,Q) to the fixed point a is

(6){\bf dist}({\mathcal E}(q,Q),a) = \max_{\langle l,l\rangle=1}\left(\langle l,a\rangle -
\rho(l ~|~ {\mathcal E}(q,Q)) \right) =
\max_{\langle l,l\rangle=1}\left(\langle l,a\rangle - \langle l,q\rangle -
\langle l,Ql\rangle^{1/2}\right).

If {\bf dist}({\mathcal E}(q,Q),a) > 0, a lies outside {\mathcal E}(q,Q); if {\bf dist}({\mathcal E}(q,Q),a) = 0, a is a boundary point of {\mathcal E}(q,Q); if {\bf dist}({\mathcal E}(q,Q),a) < 0, a is an internal point of {\mathcal E}(q,Q).

Given two ellipsoids, {\mathcal E}(q_1,Q_1) and {\mathcal E}(q_2,Q_2), the distance between them is

(7)\begin{aligned}
{\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) & = \max_{\langle l,l\rangle=1}
\left(-\rho(-l ~|~ {\mathcal E}(q_1,Q_1)) - \rho(l ~|~ {\mathcal E}(q_2,Q_2))\right) \\
& = \max_{\langle l,l\rangle=1}\left(\langle l,q_1\rangle -
\langle l,Q_1l\rangle^{1/2} - \langle l,q_2\rangle -
\langle l,Q_2l\rangle^{1/2}\right).
\end{aligned}

If {\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) > 0, the ellipsoids have no common points; if {\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) = 0, the ellipsoids have one common point - they touch; if {\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) < 0, the ellipsoids intersect.

Finding {\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) using QCQP is

d({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) = \min \langle (x-y), (x-y)\rangle

subject to:

\begin{aligned}
\langle (q_1-x), Q_1^{-1}(q_1-x)\rangle & \leqslant& 1,\\
\langle (q_2-x), Q_2^{-1}(q_2-y)\rangle & \leqslant& 1,\end{aligned}

where

d({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2))=\left\{\begin{array}{ll}
{\bf dist}^2({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2)) &
\mbox{ if } {\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2))>0, \\
0 & \mbox{ otherwise}. \end{array}\right.

Checking if k nondegenerate ellipsoids {\mathcal E}(q_1,Q_1),\cdots,{\mathcal E}(q_k,Q_k) have nonempty intersection, can be cast as a quadratically constrained quadratic programming (QCQP) problem:

\min 0

subject to:

\langle (x-q_i),Q_i^{-1}(x-q_i)\rangle - 1 \leqslant0, ~~~ i=1,\cdots,k.

If this problem is feasible, the intersection is nonempty.

Definition. Given compact convex set {\mathcal X}\subseteq{\bf R}^n, its polar set, denoted {\mathcal X}^\circ, is

{\mathcal X}^\circ = \{x\in{\bf R}^n ~|~ \langle x,y\rangle\leqslant1, ~ y\in{\mathcal X}\},

or, equivalently,

{\mathcal X}^\circ = \{l\in{\bf R}^n ~|~ \rho(l ~|~ {\mathcal X})\leqslant1\}.

The properties of the polar set are

  • If {\mathcal X} contains the origin, ({\mathcal X}^\circ)^\circ = {\mathcal X};
  • If {\mathcal X}_1\subseteq{\mathcal X}_2, {\mathcal X}_2^\circ\subseteq{\mathcal X}_1^\circ;
  • For any nonsingular matrix A\in{\bf R}^{n\times n}, (A{\mathcal X})^\circ = (A^T)^{-1}{\mathcal X}^\circ.

If a nondegenerate ellipsoid {\mathcal E}(q,Q) contains the origin, its polar set is also an ellipsoid:

\begin{aligned}
{\mathcal E}^\circ(q,Q) & = \{l\in{\bf R}^n ~|~ \langle l,q\rangle +
\langle l,Ql\rangle^{1/2}\leqslant1 \}\\
& = \{l\in{\bf R}^n ~|~ \langle l,(Q-qq^T)^{-1}l\rangle +
2\langle l,q\rangle\leqslant1 \}\\
& = \{l\in{\bf R}^n ~|~ \langle(l+(Q-qq^T)^{-1}q),
(Q-qq^T)(l+(Q-qq^T)^{-1}q)\rangle\leqslant1+\langle q,(Q-qq^T)^{-1}q\rangle \}.\end{aligned}

The special case is

{\mathcal E}^\circ(0,Q) = {\mathcal E}(0,Q^{-1}).

Definition. Given k compact sets {\mathcal X}_1, \cdots, {\mathcal X}_k\subseteq{\bf R}^n, their geometric (Minkowski) sum is

(8){\mathcal X}_1\oplus\cdots\oplus{\mathcal X}_k=\bigcup_{x_1\in{\mathcal X}_1}\cdots\bigcup_{x_k\in{\mathcal X}_k}
\{x_1 + \cdots + x_k\} .

Definition. Given two compact sets {\mathcal X}_1, {\mathcal X}_2 \subseteq{\bf R}^n, their geometric (Minkowski) difference is

(9){\mathcal X}_1\dot{-}{\mathcal X}_2 = \{x\in{\bf R}^n ~|~ x + {\mathcal X}_2 \subseteq {\mathcal X}_1 \}.

Ellipsoidal calculus concerns the following set of operations:

  • affine transformation of ellipsoid;
  • geometric sum of finite number of ellipsoids;
  • geometric difference of two ellipsoids;
  • intersection of finite number of ellipsoids.

These operations occur in reachability calculation and verification of piecewise affine dynamical systems. The result of all of these operations, except for the affine transformation, is not generally an ellipsoid but some convex set, for which we can compute external and internal ellipsoidal approximations.

Additional operations implemented in the Ellipsoidal Toolbox include external and internal approximations of intersections of ellipsoids with hyperplanes, halfspaces and polytopes.

Definition. Hyperplane H(c,\gamma) in {\bf R}^n is the set

(10)H = \{x\in{\bf R}^n ~|~ \langle c, x\rangle = \gamma\}

with c\in{\bf R}^n and \gamma\in{\bf R} fixed. The distance from ellipsoid {\mathcal E}(q,Q) to hyperplane H(c,\gamma) is

(11){\bf dist}({\mathcal E}(q,Q),H(c,\gamma)) =
\frac{\left|\gamma-\langle c,q\rangle\right| -
\langle c,Qc\rangle^{1/2}}{\langle c,c\rangle^{1/2}}.

If {\bf dist}({\mathcal E}(q,Q),H(c,\gamma))>0, the ellipsoid and the hyperplane do not intersect; if {\bf dist}({\mathcal E}(q,Q),H(c,\gamma))=0, the hyperplane is a supporting hyperplane for the ellipsoid; if {\bf dist}({\mathcal E}(q,Q),H(c,\gamma))<0, the ellipsoid intersects the hyperplane. The intersection of an ellipsoid with a hyperplane is always an ellipsoid and can be computed directly.

Checking if the intersection of k nondegenerate ellipsoids E(q_1,Q_1),\cdots,{\mathcal E}(q_k,Q_k) intersects hyperplane H(c,\gamma), is equivalent to the feasibility check of the QCQP problem:

\min 0

subject to:

\begin{aligned}
\langle (x-q_i),Q_i^{-1}(x-q_i)\rangle - 1 \leqslant0, & & i=1,\cdots,k,\\
\langle c, x\rangle - \gamma = 0. & &\end{aligned}

A hyperplane defines two (closed) halfspaces:

(12){\bf S}_1 = \{x\in{\bf R}^n ~|~ \langle c, x\rangle \leqslant\gamma\}

and

(13){\bf S}_2 = \{x\in{\bf R}^n ~|~ \langle c, x\rangle \geqslant\gamma\}.

To avoid confusion, however, we shall further assume that a hyperplane H(c,\gamma) specifies the halfspace in the sense (12). In order to refer to the other halfspace, the same hyperplane should be defined as H(-c,-\gamma).

The idea behind the calculation of intersection of an ellipsoid with a halfspace is to treat the halfspace as an unbounded ellipsoid, that is, as the ellipsoid with the shape matrix all but one of whose eigenvalues are \infty.

Definition. Polytope P(C,g) is the intersection of a finite number of closed halfspaces:

(14)P = \{x\in{\bf R}^n ~|~ Cx\leqslant g\},

wherein C=[c_1 ~ \cdots ~ c_m]^T\in{\bf R}^{m\times n} and g=[\gamma_1 ~ \cdots ~ \gamma_m]^T\in{\bf R}^m. The distance from ellipsoid {\mathcal E}(q,Q) to the polytope P(C,g) is

(15){\bf dist}({\mathcal E}(q,Q),P(C,g))=\min_{y\in P(C,g)}{\bf dist}({\mathcal E}(q,Q),y),

where {\bf dist}({\mathcal E}(q,Q),y) comes from ([dist:sub:point]). If {\bf dist}({\mathcal E}(q,Q),P(C,g))>0, the ellipsoid and the polytope do not intersect; if {\bf dist}({\mathcal E}(q,Q),P(C,g))=0, the ellipsoid touches the polytope; if {\bf dist}({\mathcal E}(q,Q),P(C,g))<0, the ellipsoid intersects the polytope.

Checking if the intersection of k nondegenerate ellipsoids E(q_1,Q_1),\cdots,{\mathcal E}(q_k,Q_k) intersects polytope P(C,g) is equivalent to the feasibility check of the QCQP problem:

\min 0

subject to:

\begin{aligned}
\langle (x-q_i),Q_i^{-1}(x-q_i)\rangle - 1 \leqslant0, & & i=1,\cdots,k,\\
\langle c_j, x\rangle - \gamma_j \leqslant0, & & j=1,\cdots,m.\end{aligned}

Operations with Ellipsoids

Affine Transformation

The simplest operation with ellipsoids is an affine transformation. Let ellipsoid {\mathcal E}(q,Q)\subseteq{\bf R}^n, matrix A\in{\bf R}^{m\times n} and vector b\in{\bf R}^m. Then

(16)A{\mathcal E}(q,Q) + b = {\mathcal E}(Aq+b, AQA^T) .

Thus, ellipsoids are preserved under affine transformation. If the rows of A are linearly independent (which implies m\leqslant n), and b=0, the affine transformation is called projection.

Geometric Sum

Consider the geometric sum (8) in which {\mathcal X}_1,\cdots,{\mathcal X}_k are nondegenerate ellipsoids {\mathcal E}(q_1,Q_1),\cdots, {\mathcal E}(q_k,Q_k)\subseteq{\bf R}^n. The resulting set is not generally an ellipsoid. However, it can be tightly approximated by the parametrized families of external and internal ellipsoids.

Let parameter l be some nonzero vector in {\bf R}^n. Then the external approximation {\mathcal E}(q,Q_l^+) and the internal approximation {\mathcal E}(q,Q_l^-) of the sum {\mathcal E}(q_1,Q_1)\oplus\cdots\oplus{\mathcal E}(q_k,Q_k) are tight along direction l, i.e.,

{\mathcal E}(q,Q_l^-)\subseteq{\mathcal E}(q_1,Q_1)\oplus\cdots\oplus{\mathcal E}(q_k,Q_k)
\subseteq{\mathcal E}(q,Q_l^+)

and

\rho(\pm l ~|~ {\mathcal E}(q,Q_l^-)) =
\rho(\pm l ~|~ {\mathcal E}(q_1,Q_1)\oplus\cdots\oplus{\mathcal E}(q_k,Q_k)) =
\rho(\pm l ~|~ {\mathcal E}(q,Q_l^+)).

Here the center q is

(17)q = q_1 + \cdots + q_k ,

the shape matrix of the external ellipsoid Q_l^+ is

(18)Q_l^+ = \left(\langle l,Q_1l\rangle^{1/2} + \cdots
+ \langle l,Q_kl\rangle^{1/2}\right)
\left(\frac{1}{\langle l,Q_1l\rangle^{1/2}}Q_1 + \cdots +
\frac{1}{\langle l,Q_kl\rangle^{1/2}}Q_k\right),

and the shape matrix of the internal ellipsoid Q_l^- is

(19)Q_l^- = \left(Q_1^{1/2} + S_2Q_2^{1/2} + \cdots + S_kQ_k^{1/2}\right)^T
\left(Q_1^{1/2} + S_2Q_2^{1/2} + \cdots + S_kQ_k^{1/2}\right),

with matrices S_i, i=2,\cdots,k, being orthogonal (S_iS_i^T=I) and such that vectors Q_1^{1/2}l, S_2Q_2^{1/2}l, \cdots, S_kQ_k^{1/2}l are parallel.

Varying vector l we get exact external and internal approximations,

\bigcup_{\langle l,l\rangle=1} {\mathcal E}(q,Q_l^-) =
{\mathcal E}(q_1,Q_1)\oplus\cdots\oplus{\mathcal E}(q_k,Q_k) =
\bigcap_{\langle l,l\rangle=1} {\mathcal E}(q,Q_l^+) .

For proofs of formulas given in this section, see [KUR1997], [KUR2000].

One last comment is about how to find orthogonal matrices S_2,\cdots,S_k that align vectors Q_2^{1/2}l, \cdots, Q_k^{1/2}l with Q_1^{1/2}l. Let v_1 and v_2 be some unit vectors in {\bf R}^n. We have to find matrix S such that Sv_2\uparrow\uparrow v_1. We suggest explicit formulas for the calculation of this matrix [DAR2012]:

(20)T = I + Q_1(S - I)Q_1^T,

(21)S = \begin{pmatrix}
     c & s\\
     -s & c
    \end{pmatrix},\quad c = \langle\hat{v_1},\ \hat{v_2}\rangle,\ \quad s = \sqrt{1 - c^2},\ \quad \hat{v_i} = \dfrac{v_i}{\|v_i\|}

(22)Q_1 = [q_1 \, q_2]\in \mathbb{R}^{n\times2},\ \quad q_1 = \hat{v_1},\ \quad q_2 = \begin{cases}
s^{-1}(\hat{v_2} - c\hat{v_1}),& s\ne 0\\
0,& s = 0.
\end{cases}

Geometric Difference

Consider the geometric difference (9) in which the sets {\mathcal X}_1 and {\mathcal X}_2 are nondegenerate ellipsoids {\mathcal E}(q_1,Q_1) and {\mathcal E}(q_2,Q_2). We say that ellipsoid {\mathcal E}(q_1,Q_1) is bigger than ellipsoid {\mathcal E}(q_2,Q_2) if

{\mathcal E}(0,Q_2) \subseteq {\mathcal E}(0,Q_1).

If this condition is not fulfilled, the geometric difference {\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2) is an empty set:

{\mathcal E}(0,Q_2) \not\subseteq {\mathcal E}(0,Q_1) ~~~ \Rightarrow ~~~
{\mathcal E}(q_1,Q_1) \dot{-}{\mathcal E}(q_2,Q_2) = \emptyset.

If {\mathcal E}(q_1,Q_1) is bigger than {\mathcal E}(q_2,Q_2) and {\mathcal E}(q_2,Q_2) is bigger than {\mathcal E}(q_1,Q_1), in other words, if Q_1=Q_2,

{\mathcal E}(q_1,Q_1) \dot{-}{\mathcal E}(q_2,Q_2) = \{q_1-q_2\} ~~~ \mbox{and} ~~~
{\mathcal E}(q_2,Q_2) \dot{-}{\mathcal E}(q_1,Q_1) = \{q_2-q_1\}.

To check if ellipsoid {\mathcal E}(q_1,Q_1) is bigger than ellipsoid {\mathcal E}(q_2,Q_2), we perform simultaneous diagonalization of matrices Q_1 and Q_2, that is, we find matrix T such that

TQ_1T^T = I ~~~ \mbox{and} ~~~ TQ_2T^T=D,

where D is some diagonal matrix. Simultaneous diagonalization of Q_1 and Q_2 is possible because both are symmetric positive definite (see [GANT1960]). To find such matrix T, we first do the SVD of Q_1:

(23)Q_1 = U_1\Sigma_1V_1^T .

Then the SVD of matrix \Sigma_1^{-1/2}U_1^TQ_2U_1\Sigma_1^{-1/2}:

(24)\Sigma_1^{-1/2}U_1^TQ_2U_1\Sigma_1^{-1/2} = U_2\Sigma_2V_2^T.

Now, T is defined as

(25)T = U_2^T \Sigma_1^{-1/2}U_1^T.

If the biggest diagonal element (eigenvalue) of matrix D=TQ_2T^T is less than or equal to 1, {\mathcal E}(0,Q_2)\subseteq{\mathcal E}(0,Q_1).

Once it is established that ellipsoid {\mathcal E}(q_1,Q_1) is bigger than ellipsoid {\mathcal E}(q_2,Q_2), we know that their geometric difference {\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2) is a nonempty convex compact set. Although it is not generally an ellipsoid, we can find tight external and internal approximations of this set parametrized by vector l\in{\bf R}^n. Unlike geometric sum, however, ellipsoidal approximations for the geometric difference do not exist for every direction l. Vectors for which the approximations do not exist are called bad directions.

Given two ellipsoids {\mathcal E}(q_1,Q_1) and {\mathcal E}(q_2,Q_2) with {\mathcal E}(0,Q_2)\subseteq{\mathcal E}(0,Q_1), l is a bad direction if

\frac{\langle l,Q_1l\rangle^{1/2}}{\langle l,Q_2l\rangle^{1/2}}>r,

in which r is a minimal root of the equation

{\bf det}(Q_1-rQ_2) = 0.

To find r, compute matrix T by (23)-(25) and define

r = \frac{1}{\max({\bf diag}(TQ_2T^T))}.

If l is not a bad direction, we can find tight external and internal ellipsoidal approximations {\mathcal E}(q,Q^+_l) and {\mathcal E}(q,Q^-_l) such that

{\mathcal E}(q,Q^-_l)\subseteq{\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2)\subseteq{\mathcal E}(q,Q^+_l)

and

\rho(\pm l ~|~ {\mathcal E}(q,Q_l^-)) =
\rho(\pm l ~|~ {\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2)) =
\rho(\pm l ~|~ {\mathcal E}(q,Q_l^+)).

The center q is

(26)q = q_1 - q_2;

the shape matrix of the internal ellipsoid Q^-_l is

\begin{aligned}
&& P = \frac{\sqrt{\langle l, Q_1 l\rangle}}{\sqrt{\langle l, Q_2 \rangle}};\nonumber\\
&& Q^-_l = \left(1 - \dfrac{1}{P}\right)Q_1 + \left(1 - P\right)Q_2.
\label{minkdiff_ia}\end{aligned}

and the shape matrix of the external ellipsoid Q^+_l is

(27)Q^+_l = \left(Q_1^{1/2} - SQ_2^{1/2}\right)^T
\left(Q_1^{1/2} - SQ_2^{1/2}\right).

Here S is an orthogonal matrix such that vectors Q_1^{1/2}l and SQ_2^{1/2}l are parallel. S is found from (20)-(22), with v_1=Q_2^{1/2}l and v_2=Q_1^{1/2}l.

Running l over all unit directions that are not bad, we get

\bigcup_{\langle l,l\rangle=1} {\mathcal E}(q,Q_l^-) =
{\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2) =
\bigcap_{\langle l,l\rangle=1} {\mathcal E}(q,Q_l^+) .

For proofs of formulas given in this section, see [KUR1997].

Geometric Difference-Sum

Given ellipsoids {\mathcal E}(q_1,Q_1), {\mathcal E}(q_2,Q_2) and {\mathcal E}(q_3,Q_3), it is possible to compute families of external and internal approximating ellipsoids for

(28){\mathcal E}(q_1,Q_1) \dot{-} {\mathcal E}(q_2,Q_2) \oplus {\mathcal E}(q_3,Q_3)

parametrized by direction l, if this set is nonempty ({\mathcal E}(0,Q_2)\subseteq{\mathcal E}(0,Q_1)).

First, using the result of the previous section, for any direction l that is not bad, we obtain tight external {\mathcal E}(q_1-q_2, Q_l^{0+}) and internal {\mathcal E}(q_1-q_2, Q_l^{0-}) approximations of the set {\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2).

The second and last step is, using the result of section 2.2.2, to find tight external ellipsoidal approximation {\mathcal E}(q_1-q_2+q_3,Q_l^+) of the sum {\mathcal E}(q_1-q_2,Q_l^{0+})\oplus{\mathcal E}(q_3,Q_3), and tight internal ellipsoidal approximation {\mathcal E}(q_1-q_2+q_3,Q_l^-) for the sum {\mathcal E}(q_1-q_2,Q_l^{0-})\oplus{\mathcal E}(q_3,Q_3).

As a result, we get

{\mathcal E}(q_1-q_2+q_3,Q_l^-) \subseteq
{\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2)\oplus{\mathcal E}(q_3,Q_3) \subseteq
{\mathcal E}(q_1-q_2+q_3,Q_l^+)

and

\rho(\pm l ~|~{\mathcal E}(q_1-q_2+q_3,Q_l^-)) =
\rho(\pm l ~|~ {\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2)\oplus{\mathcal E}(q_3,Q_3)) =
\rho(\pm l ~|~ {\mathcal E}(q_1-q_2+q_3,Q_l^+)).

Running l over all unit vectors that are not bad, this translates to

\bigcup_{\langle l,l\rangle=1} {\mathcal E}(q_1-q_2+q_3,Q_l^-) =
{\mathcal E}(q_1,Q_1)\dot{-}{\mathcal E}(q_2,Q_2)\oplus{\mathcal E}(q_3,Q_3) =
\bigcap_{\langle l,l\rangle=1} {\mathcal E}(q_1-q_2+q_3,Q_l^+) .

Geometric Sum-Difference

Given ellipsoids {\mathcal E}(q_1,Q1), {\mathcal E}(q_2,Q_2) and {\mathcal E}(q_3,Q_3), it is possible to compute families of external and internal approximating ellipsoids for

(29){\mathcal E}(q_1,Q_1) \oplus {\mathcal E}(q_2,Q_2) \dot{-} {\mathcal E}(q_3,Q_3)

parametrized by direction l, if this set is nonempty ({\mathcal E}(0,Q_3)\subseteq{\mathcal E}(0,Q_1)\oplus{\mathcal E}(0,Q_2)).

First, using the result of section 2.2.2, we obtain tight external {\mathcal E}(q_1+q_2,Q_l^{0+}) and internal {\mathcal E}(q_1+q_2,Q_l^{0-}) ellipsoidal approximations of the set {\mathcal E}(q_1,Q_1)\oplus{\mathcal E}(q_2,Q_2). In order for the set (29) to be nonempty, inclusion {\mathcal E}(0,Q_3)\subseteq{\mathcal E}(0,Q_l^{0+}) must be true for any l. Note, however, that even if (29) is nonempty, it may be that {\mathcal E}(0,Q_3)\not\subseteq{\mathcal E}(0,Q_l^{0-}), then internal approximation for this direction does not exist.

Assuming that (29) is nonempty and {\mathcal E}(0,Q_3)\subseteq{\mathcal E}(0,Q_l^{0-}), the second step would be, using the results of section 2.2.3, to compute tight external ellipsoidal approximation {\mathcal E}(q_1+q_2-q_3,Q_l^+) of the difference {\mathcal E}(q_1+q_2,Q_l^{0+})\dot{-}{\mathcal E}(q_3,Q_3), which exists only if l is not bad, and tight internal ellipsoidal approximation {\mathcal E}(q_1+q_2-q_3,Q_l^-) of the difference {\mathcal E}(q_1+q_2,Q_l^{0-})\dot{-}{\mathcal E}(q_3,Q_3), which exists only if l is not bad for this difference.

If approximation {\mathcal E}(q_1+q_2-q_3,Q_l^+) exists, then

{\mathcal E}(q_1,Q_1)\oplus{\mathcal E}(q_2,Q_2)\dot{-}{\mathcal E}(q_3,Q_3) \subseteq
{\mathcal E}(q_1+q_2-q_3,Q_l^+)

and

\rho(\pm l ~|~ {\mathcal E}(q_1,Q_1)\oplus{\mathcal E}(q_2,Q_2)\dot{-}{\mathcal E}(q_3,Q_3)) =
\rho(\pm l ~|~ {\mathcal E}(q_1+q_2-q_3,Q_l^+)).

If approximation {\mathcal E}(q_1+q_2-q_3,Q_l^-) exists, then

{\mathcal E}(q_1+q_2-q_3,Q_l^-) \subseteq
{\mathcal E}(q_1,Q_1)\oplus{\mathcal E}(q_2,Q_2)\dot{-}{\mathcal E}(q_3,Q_3)

and

\rho(\pm l ~|~{\mathcal E}(q_1+q_2-q_3,Q_l^-)) =
\rho(\pm l ~|~ {\mathcal E}(q_1,Q_1)\oplus{\mathcal E}(q_2,Q_2)\dot{-}{\mathcal E}(q_3,Q_3)) .

For any fixed direction l it may be the case that neither external nor internal tight ellipsoidal approximations exist.

Intersection of Ellipsoid and Hyperplane

Let nondegenerate ellipsoid {\mathcal E}(q,Q) and hyperplane H(c,\gamma) be such that {\bf dist}({\mathcal E}(q,Q),H(c,\gamma))<0. In other words,

{\mathcal E}_H(w,W) = {\mathcal E}(q,Q)\cap H(c,\gamma) \neq \emptyset .

The intersection of ellipsoid with hyperplane, if nonempty, is always an ellipsoid. Here we show how to find it.

First of all, we transform the hyperplane H(c,\gamma) into H([1~0~\cdots~0]^T, 0) by the affine transformation

y = Sx - \frac{\gamma}{\langle c,c\rangle^{1/2}}Sc,

where S is an orthogonal matrix found by (20)-(22) with v_1=c and v_2=[1~0~\cdots~0]^T. The ellipsoid in the new coordinates becomes {\mathcal E}(q',Q') with

\begin{aligned}
q' & = q-\frac{\gamma}{\langle c,c\rangle^{1/2}}Sc, \\
Q' & = SQS^T.\end{aligned}

Define matrix M=Q'^{-1}; m_{11} is its element in position (1,1), \bar{m} is the first column of M without the first element, and \bar{M} is the submatrix of M obtained by stripping M of its first row and first column:

M = \left[\begin{array}{c|cl}
m_{11} & & \bar{m}^T\\
 & \\
\hline
 & \\
\bar{m} & & \bar{M}\end{array}\right].

The ellipsoid resulting from the intersection is {\mathcal E}_H(w',W') with

\begin{aligned}
w' & = q' + q_1'\left[\begin{array}{c}
-1\\
\bar{M}^{-1}\bar{m}\end{array}\right],\\
W' & = \left(1-q_1'^2(m_{11}-
\langle\bar{m},\bar{M}^{-1}\bar{m}\rangle)\right)\left[\begin{array}{c|cl}
0 & & {\bf 0}\\
 & \\
\hline
 & \\
{\bf 0} & & \bar{M}^{-1}\end{array}\right],\end{aligned}

in which q_1' represents the first element of vector q'.

Finally, it remains to do the inverse transform of the coordinates to obtain ellipsoid {\mathcal E}_H(w,W):

\begin{aligned}
w & = S^Tw' + \frac{\gamma}{\langle c,c\rangle^{1/2}}c, \\
W & = S^TW'S.\end{aligned}

Intersection of Ellipsoid and Ellipsoid

Given two nondegenerate ellipsoids {\mathcal E}(q_1,Q_1) and {\mathcal E}(q_2,Q_2), {\bf dist}({\mathcal E}(q_1,Q_1),{\mathcal E}(q_2,Q_2))<0 implies that

{\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2)\neq\emptyset .

This intersection can be approximated by ellipsoids from the outside and from the inside. Trivially, both {\mathcal E}(q_1,Q_1) and {\mathcal E}(q_2,Q_2) are external approximations of this intersection. Here, however, we show how to find the external ellipsoidal approximation of minimal volume.

Define matrices

W_1 = Q_1^{-1}, ~~~~ W_2 = Q_2^{-1} .\label{wmatrices}

Minimal volume external ellipsoidal approximation {\mathcal E}(q+,Q^+) of the intersection {\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2) is determined from the set of equations:

(30)Q^+  = \alpha X^{-1}, \\

(31)X  =  \pi W_1 + (1-\pi)W_2,\\

(32)\alpha  =  1-\pi(1-\pi)\langle(q_2-q_1), W_2X^{-1}W_1(q_2-q_1)\rangle, \\

(33)q^+  = X^{-1}(\pi W_1q_1 + (1-\pi)W_2q_2), \\

(34)0 &=  \alpha({\bf det}(X))^2{\bf trace}(X^{-1}(W_1-W_2)) - {}\\
  &- n({\bf det}(X))^2 (2\langle q^+,W_1q_1-W_2q_2\rangle + \langle q^+,(W_2-W_1)q^+\rangle - {}\\
  &- \langle q_1,W_1q_1\rangle + \langle q_2,W_2q_2\rangle),

with 0\leqslant\pi\leqslant1. We substitute X, \alpha, q^+ defined in (31)-(33) into (34) and get a polynomial of degree 2n-1 with respect to \pi, which has only one root in the interval [0,1], \pi_0. Then, substituting \pi=\pi_0 into (30)-(33), we obtain q^+ and Q^+. Special cases are \pi_0=1, whence {\mathcal E}(q^+,Q^+)={\mathcal E}(q_1,Q_1), and \pi_0=0, whence {\mathcal E}(q^+,Q^+)={\mathcal E}(q_2,Q_2). These situations may occur if, for example, one ellipsoid is contained in the other:

{\mathcal E}(q_1,Q_1)\subseteq{\mathcal E}(q_2,Q_2) & \Rightarrow \pi_0 = 1,\\
{\mathcal E}(q_2,Q_2)\subseteq{\mathcal E}(q_1,Q_1) & \Rightarrow \pi_0 = 0.

The proof that the system of equations (30)-(34) correctly defines the minimal volume external ellipsoidal approximationi of the intersection {\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2) is given in [ROS2002].

To find the internal approximating ellipsoid {\mathcal E}(q^-,Q^-)\subseteq{\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2), define

(35)\beta_1 = \min_{\langle x,W_2x\rangle=1}\langle x,W_1x\rangle,

(36)\beta_2 = \min_{\langle x,W_1x\rangle=1}\langle x,W_2x\rangle,

Notice that (35) and (36) are QCQP problems. Parameters \beta_1 and \beta_2 are invariant with respect to affine coordinate transformation and describe the position of ellipsoids {\mathcal E}(q_1,Q_1), {\mathcal E}(q_2,Q_2) with respect to each other:

\beta_1\geqslant1,~\beta_2\geqslant1 & \Rightarrow
{\bf int}({\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2))=\emptyset, \\

\beta_1\geqslant1,~\beta_2\leqslant1 & \Rightarrow {\mathcal E}(q_1,Q_1)\subseteq{\mathcal E}(q_2,Q_2), \\

\beta_1\leqslant1,~\beta_2\geqslant1 & \Rightarrow {\mathcal E}(q_2,Q_2)\subseteq{\mathcal E}(q_1,Q_1), \\

\beta_1<1,~\beta_2<1 & \Rightarrow
{\bf int}({\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2))\neq\emptyset \\

&\mbox{and} ~ {\mathcal E}(q_1,Q_1)\not\subseteq{\mathcal E}(q_2,Q_2) \\

&\mbox{and} ~ {\mathcal E}(q_2,Q_2)\not\subseteq{\mathcal E}(q_1,Q_1).

Define parametrized family of internal ellipsoids {\mathcal E}(q^-_{\theta_1\theta_2},Q^-_{\theta_1\theta_2}) with

(37)q^-_{\theta_1\theta_2}  =  (\theta_1W_1 +
\theta_2W_2)^{-1}(\theta_1W_1q_1 + \theta_2W_2q_2),\\

(38)Q^-_{\theta_1\theta_2} =  (1 - \theta_1\langle q_1,W_1q_1\rangle -
\theta_2\langle q_2,W_2q_2\rangle +
\langle q^-_{\theta_1\theta_2},(Q^-)^{-1}q^-_{\theta_1\theta_2}\rangle)
(\theta_1W_1 + \theta_2W_2)^{-1} .

The best internal ellipsoid {\mathcal E}(q^-_{\hat{\theta}_1\hat{\theta}_2},Q^-_{\hat{\theta}_1\hat{\theta}_2}) in the class (37)-(38), namely, such that

{\mathcal E}(q^-_{{\theta}_1{\theta}_2},Q^-_{{\theta}_1{\theta}_2})\subseteq
{\mathcal E}(q^-_{\hat{\theta}_1\hat{\theta}_2},Q^-_{\hat{\theta}_1\hat{\theta}_2})
\subseteq {\mathcal E}(q_1,Q_1)\cap{\mathcal E}(q_2,Q_2)

for all 0\leqslant\theta_1,\theta_2\leqslant1, is specified by the parameters

(39)\hat{\theta}_1 = \frac{1-\hat{\beta}_2}{1-\hat{\beta}_1\hat{\beta}_2}, ~~~~
\hat{\theta}_2 = \frac{1-\hat{\beta}_1}{1-\hat{\beta}_1\hat{\beta}_2},

with

\hat{\beta}_1=\min(1,\beta_1), ~~~~ \hat{\beta}_2=\min(1,\beta_2).

It is the ellipsoid that we look for: {\mathcal E}(q^-,Q^-)={\mathcal E}(q^-_{\hat{\theta}_1\hat{\theta}_2},Q^-_{\hat{\theta}_1\hat{\theta}_2}). Two special cases are

\hat{\theta}_1=1, ~ \hat{\theta}_2=0 ~~~ \Rightarrow ~~~
{\mathcal E}(q_1,Q_1)\subseteq{\mathcal E}(q_2,Q_2) ~~~ \Rightarrow ~~~
{\mathcal E}(q^-,Q^-)={\mathcal E}(q_1,Q_1),

and

\hat{\theta}_1=0, ~ \hat{\theta}_2=1 ~~~ \Rightarrow ~~~
{\mathcal E}(q_2,Q_2)\subseteq{\mathcal E}(q_1,Q_1) ~~~ \Rightarrow ~~~
{\mathcal E}(q^-,Q^-)={\mathcal E}(q_2,Q_2).

The method of finding the internal ellipsoidal approximation of the intersection of two ellipsoids is described in [VAZ1999].

Intersection of Ellipsoid and Halfspace

Finding the intersection of ellipsoid and halfspace can be reduced to finding the intersection of two ellipsoids, one of which is unbounded. Let {\mathcal E}(q_1,Q_1) be a nondegenerate ellipsoid and let H(c,\gamma) define the halfspace

{\bf S}(c,\gamma) = \{x\in{\bf R}^n ~|~ \langle c,x\rangle\leqslant\gamma\}.

We have to determine if the intersection {\mathcal E}(q_1,Q_1)\cap{\bf S}(c,\gamma) is empty, and if not, find its external and internal ellipsoidal approximations, {\mathcal E}(q^+,Q^+) and {\mathcal E}(q^-,Q^-). Two trivial situations are:

  • {\bf dist}({\mathcal E}(q_1,Q_1),H(c,\gamma))>0 and \langle c, q_1\rangle>0, which implies that {\mathcal E}(q_1,Q_1)\cap{\bf S}(c,\gamma)=\emptyset;
  • {\bf dist}({\mathcal E}(q_1,Q_1),H(c,\gamma))>0 and \langle c, q_1\rangle<0, so that {\mathcal E}(q_1,Q_1)\subseteq{\bf S}(c,\gamma), and then {\mathcal E}(q^+,Q^+)={\mathcal E}(q^-,Q^-)={\mathcal E}(q_1,Q_1).

In case {\bf dist}({\mathcal E}(q_1,Q_1),H(c,\gamma)<0, i.e. the ellipsoid intersects the hyperplane,

{\mathcal E}(q_1,Q_1)\cap{\bf S}(c,\gamma) =
{\mathcal E}(q_1,Q_1)\cap\{x ~|~ \langle (x-q_2),W_2(x-q_2)\rangle\leqslant1\},

with

(40)q_2  =  (\gamma + 2\sqrt{\overline{\lambda}})c,\\

(41)W_2  =  \frac{1}{4\overline{\lambda}}cc^T,

\overline{\lambda} being the biggest eigenvalue of matrix Q_1. After defining W_1=Q_1^{-1}, we obtain {\mathcal E}(q^+,Q^+) from equations (30)-(34), and {\mathcal E}(q^-,Q^-) from (37)-(38), (39).

Remark. Notice that matrix W_2 has rank 1, which makes it singular for n>1. Nevertheless, expressions (30)-(31), (37)-(38) make sense because W_1 is nonsingular, \pi_0\neq0 and \hat{\theta}_1\neq0.

To find the ellipsoidal approximations {\mathcal E}(q^+,Q^+) and {\mathcal E}(q^-,Q^-) of the intersection of ellipsoid {\mathcal E}(q,Q) and polytope P(C,g), C\in{\bf R}^{m\times n}, b\in{\bf R}^m, such that

{\mathcal E}(q^-,Q^-)\subseteq{\mathcal E}(q,Q)\cap P(C,g)\subseteq{\mathcal E}(q^+,Q^+),

we first compute

{\mathcal E}(q^-_1,Q^-_1)\subseteq{\mathcal E}(q,Q)\cap{\bf S}(c_1,\gamma_1)\subseteq
{\mathcal E}(q^+_1,Q^+_1),

wherein {\bf S}(c_1,\gamma_1) is the halfspace defined by the first row of matrix C, c_1, and the first element of vector g, \gamma_1. Then, one by one, we get

\begin{aligned}
& & {\mathcal E}(q^-_2,Q^-_2)\subseteq{\mathcal E}(q^-_1,Q^-_1)\cap{\bf S}(c_2,\gamma_2), ~~~
{\mathcal E}(q^+_1,Q^+_1)\cap{\bf S}(c_2,\gamma_2)\subseteq{\mathcal E}(q^+_2,Q^+_2), \\
& & {\mathcal E}(q^-_3,Q^-_3)\subseteq{\mathcal E}(q^-_2,Q^-_2)\cap{\bf S}(c_3,\gamma_3), ~~~
{\mathcal E}(q^+_2,Q^+_2)\cap{\bf S}(c_3,\gamma_3)\subseteq{\mathcal E}(q^+_3,Q^+_3), \\
& & \cdots \\
& & {\mathcal E}(q^-_m,Q^-_m)\subseteq{\mathcal E}(q^-_{m-1},Q^-_{m-1})\cap{\bf S}(c_m,\gamma_m), ~~~
{\mathcal E}(q^+_{m-1},Q^+_{m-1})\cap{\bf S}(c_m,\gamma_m)\subseteq{\mathcal E}(q^+_m,Q^+_m), \\\end{aligned}

The resulting ellipsoidal approximations are

{\mathcal E}(q^+,Q^+)={\mathcal E}(q^+_m,Q^+_m), ~~~~ {\mathcal E}(q^-,Q^-)={\mathcal E}(q^-_m,Q^-_m) .

Checking if one ellipsoid contains another

Theorem of alternatives, also known as S-procedure [BOYD2004], states that the implication

\langle x, A_1x\rangle + 2\langle b_1,x\rangle + c_1 \leqslant0
~~ \Rightarrow ~~
\langle x, A_2x\rangle + 2\langle b_2,x\rangle + c_2 \leqslant0,

where A_i\in{\bf R}^{n\times n} are symmetric matrices, b_i\in{\bf R}^n, c_i\in{\bf R}, i=1,2, holds if and only if there exists \lambda>0 such that

\left[\begin{array}{cc}
A_2 & b_2\\
b_2^T & c_2\end{array}\right]
\preceq
\lambda\left[\begin{array}{cc}
A_1 & b_1\\
b_1^T & c_1\end{array}\right].

By S-procedure, {\mathcal E}(q_1,Q_1)\subseteq{\mathcal E}(q_2,Q_2) (both ellipsoids are assumed to be nondegenerate) if and only if the following SDP problem is feasible:

\min 0

subject to:

\begin{aligned}
\lambda & >  0, \\
\left[\begin{array}{cc}
Q_2^{-1} & -Q_2^{-1}q_2\\
(-Q_2^{-1}q_2)^T & q_2^TQ_2^{-1}q_2-1\end{array}\right]
& \preceq
\lambda \left[\begin{array}{cc}
Q_1^{-1} & -Q_1^{-1}q_1\\
(-Q_1^{-1}q_1)^T & q_1^TQ_1^{-1}q_1-1\end{array}\right]\end{aligned}

where \lambda\in{\bf R} is the variable.

Minimum Volume Ellipsoids

The minimum volume ellipsoid that contains set S is called Löwner-John ellipsoid of the set S. To characterize it we rewrite general ellipsoid {\mathcal E}(q,Q) as

{\mathcal E}(q,Q) = \{x ~|~ \langle (Ax + b), (Ax + b)\rangle \},

where

A = Q^{-1/2} ~~~ \mbox{ and } ~~~ b = -Aq .

For positive definite matrix A, the volume of {\mathcal E}(q,Q) is proportional to \det A^{-1}. So, finding the minimum volume ellipsoid containing S can be expressed as semidefinite programming (SDP) problem

\min \log \det A^{-1}

subject to:

\sup_{v\in S} \langle (Av + b), (Av + b)\rangle \leqslant1,

where the variables are A\in{\bf R}^{n\times n} and b\in{\bf R}^n, and there is an implicit constraint A\succ 0 (A is positive definite). The objective and constraint functions are both convex in A and b, so this problem is convex. Evaluating the constraint function, however, requires solving a convex maximization problem, and is tractable only in certain special cases.

For a finite set S=\{x_1,\cdots,x_m\}\subset{\bf R}^n, an ellipsoid covers S if and only if it covers its convex hull. So, finding the minimum volume ellipsoid covering S is the same as finding the minimum volume ellipsoid containing the polytope {\bf conv}\{x_1,\cdots,x_m\}. The SDP problem is

\min \log \det A^{-1}

subject to:

\begin{aligned}
A & \succ  0, \\
\langle (Ax_i + b), (Ax_i + b)\rangle & \leqslant 1, ~~~ i=1..m.\end{aligned}

We can find the minimum volume ellipsoid containing the union of ellipsoids \bigcup_{i=1}^m{\mathcal E}(q_i,Q_i). Using the fact that for i=1..m {\mathcal E}(q_i,Q_i)\subseteq{\mathcal E}(q,Q) if and only if there exists \lambda_i>0 such that

\left[\begin{array}{cc}
A^2 - \lambda_i Q_i^{-1} & Ab + \lambda_i Q_i^{-1}q_i\\
(Ab + \lambda_i Q_i^{-1}q_i)^T & b^Tb-1 - \lambda_i (q_i^TQ_i^{-1}q_i-1) \end{array}
\right] \preceq 0 .

Changing variable \tilde{b}=Ab, we get convex SDP in the variables A, \tilde{b}, and \lambda_1,\cdots,\lambda_m:

\min \log \det A^{-1}

subject to:

\begin{aligned}
\lambda_i & > 0,\\
\left[\begin{array}{ccc}
A^2-\lambda_iQ_i^{-1} & \tilde{b}+\lambda_iQ_i^{-1}q_i & 0 \\
(\tilde{b}+\lambda_iQ_i^{-1}q_i)^T & -1-\lambda_i(q_i^TQ_i^{-1}q_i-1) & \tilde{b}^T \\
0 & \tilde{b} & -A^2\end{array}\right] & \preceq 0, ~~~ i=1..m.\end{aligned}

After A and b are found,

q=-A^{-1}b ~~~ \mbox{ and } ~~~ Q=(A^TA)^{-1}.

The results on the minimum volume ellipsoids are explained and proven in [BOYD2004].

Maximum Volume Ellipsoids

Consider a problem of finding the maximum volume ellipsoid that lies inside a bounded convex set S with nonempty interior. To formulate this problem we rewrite general ellipsoid {\mathcal E}(q,Q) as

{\mathcal E}(q,Q) = \{Bx + q ~|~ \langle x,x\rangle\leqslant1\},

where B=Q^{1/2}, so the volume of {\mathcal E}(q,Q) is proportional to \det B.

The maximum volume ellipsoid that lies inside S can be found by solving the following SDP problem:

\max \log \det B

subject to:

\sup_{\langle v,v\rangle\leqslant1} I_S(Bv+q)\leqslant0 ,

in the variables B\in{\bf R}^{n\times n} - symmetric matrix, and q\in{\bf R}^n, with implicit constraint B\succ 0, where I_S is the indicator function:

I_S(x) = \left\{\begin{array}{ll}
0, & \mbox{ if } x\in S,\\
\infty, & \mbox{ otherwise.}\end{array}\right.

In case of polytope, S=P(C,g) with P(C,g) defined in (14), the SDP has the form

\min \log \det B^{-1}

subject to:

\begin{aligned}
B & \succ 0,\\
\langle c_i, Bc_i\rangle + \langle c_i, q\rangle & \leqslant \gamma_i,
~~~ i=1..m.\end{aligned}

We can find the maximum volume ellipsoid that lies inside the intersection of given ellipsoids \bigcap_{i=1}^m{\mathcal E}(q_i,Q_i). Using the fact that for i=1..m {\mathcal E}(q,Q)\subseteq{\mathcal E}(q_i,Q_i) if and only if there exists \lambda_i>0 such that

\left[\begin{array}{cc}
-\lambda_i - q^TQ_i^{-1}q + 2q_i^TQ_i^{-1}q - q_i^TQ_i^{-1}q_i + 1 & (Q_i^{-1}q-Q_i^{-1}q_i)^TB\\
B(Q_i^{-1}q-Q_i^{-1}q_i) & \lambda_iI-BQ_i^{-1}B\end{array}\right] \succeq 0.

To find the maximum volume ellipsoid, we solve convex SDP in variables B, q, and \lambda_1,\cdots,\lambda_m:

\min \log \det B^{-1}

subject to:

\begin{aligned}
\lambda_i & >  0, \\
\left[\begin{array}{ccc}
1-\lambda_i & 0 & (q - q_i)^T\\
0 & \lambda_iI & B\\
q - q_i & B & Q_i\end{array}\right] & \succeq  0, ~~~ i=1..m.\end{aligned}

After B and q are found,

Q = B^TB.

The results on the maximum volume ellipsoids are explained and proven in [BOYD2004].

References

[DAR2012]A. N. Dariyn and A. B. Kurzhanski. Method of invariant sets for linear systems of high dimensionality under uncertain disturbances. Doklady Akademii Nauk, 446(6):607–611, 2012.
[GANT1960]
    1. Gantmacher. Matrix Theory, I-II. Chelsea, 1960.
[ROS2002]F. Thomas L. Ros, A. Sabater. An Ellipsoidal Calculus Based on Propagation and Fusion. IEEE Transactions on Systems, Man and Cybernetics, Part B: Cybernetics, 32(4), 2002.
[VAZ1999]A. Yu. Vazhentsev. On Internal Ellipsoidal Approximations for Problems of Control and Synthesis with Bounded Coordinates. Izvestiya Rossiiskoi Akademii Nauk. Teoriya i Systemi Upravleniya., 1999.
[BOYD2004](1, 2, 3)
  1. Boyd and L. Vandenberghe. Convex Optimization. Cambridge University Press, 2004.