Maximizing a quadratic form involving a trace-bounded positive definite matrix?

Question

$\newcommand{\tr}{\mathrm{tr}}$Suppose $P, Q$ are two real, symmetric positive definite matrices and $v$ a nonzero unit vector.

Consider $$ f(X) = v^T(P + X^{-1})^{-1} v + v^T(Q + X^{-1})^{-1} v. $$ Question: What is the value of$$ F^\star = \sup\{ f(X) : X \succ 0, \tr(X) \leq 1\}? $$

---

In the case $P = Q$: Since $X \preceq I$, clearly $f(X) \leq 2 v^T (P + I)^{-1} v$. Approximating $xx^T$ by positive definite matrices gives $$F^\star \geq 2\, \sup_{x : \|x\|_2 = 1} \frac{(v^T x)^2}{x^T(P+I)x} = 2 v^T (P + I)^{-1} v,$$ by Cauchy-Schwarz. Thus, if $P = Q$,$$ F^\star = 2 v^T (P + I)^{-1} v. $$ I am unsure how to proceed to $P \neq Q$.

Answer 1

We assume that the dimension of the matrices is $n\geq 2$. Let $Sym$ (resp. $Skew$) denote the sets of symmetric (resp. skew-symmetric) matrices.

i) Up to an orthonormal change of basis, we may assume that $v=[1,0,\cdots,0]^T$; indeed, the matrices $P,Q,X$ remain in $S_{>0}$ the space of symmetric $>0$ matrices.

Thus $f(X)=[(P+X^{-1})^{-1}+(Q+X^{-1})^{-1}]_{1,1}$.

ii) We search the critical points of $f$ -$Df_X(H)=0$-

$Df_X:H\in Sym\mapsto v^T(P+X^{-1})^{-1}(X^{-1}HX^{-1})(P+X^{-1})^{-1}v+(\cdots)$.

$DF_X(H)=tr([X^{-1}(P+X^{-1})^{-1}vv^T(P+X^{-1})^{-1}X^{-1}+(\cdots)]H)=$

$tr((C_P{C_P}^T+C_Q{C_Q}^T)H)$, where $C_P$ is the column $X^{-1}(P+X^{-1})^{-1}v$.

Then $X$ is critical iff $C_P{C_P}^T+C_Q{C_Q}^T\in Skew$, that is, $C_P{C_P}^T+C_Q{C_Q}^T=0$.

Thus $C_P=C_Q=0$, that is impossible. Then $f$ admits no free extremum.

iii) Note that, for every $X$, $f(X)\leq ((P+I)^{-1}+(Q+I)^{-1})_{1,1}$.

Then $f$ may have a maximum for $tr(X)=1,X>0$ or have an upper bound on the edge of $\{X\geq 0,tr(X)=1\}$ associated to a sequence $(X_k)$ s.t. $\det(X_k)\rightarrow 0$.

iv) We add the condition (1) $tr(X)=1$. The necessary condition for an extremum becomes

there is $\lambda\in \mathbb{R}$ s.t., for every $H\in Sym$, $Df_X(H)-\lambda tr(H)=0$ or $tr((C_P{C_p}^T+C_Q{C_Q}^T-\lambda I)H)=0$.

Then $C_P{C_P}^T+C_Q{C_Q}^T-\lambda I\in Skew$, that implies (2) $C_P{C_p}^T+C_Q{C_Q}^T=\lambda I$.

Since $1\leq rank(C_P{C_P}^T+C_Q{C_Q}^T)\leq 2$, the above equality is possible only when $n=2,\lambda >0$.

v) $\textbf{The case $n=2$.}$ Let $X=\begin{pmatrix}a&b\\b&1-a\end{pmatrix}$. We must solve the system (2) of $2$ equations with $2$ unknowns. This algebraic system, when it is generic, can have around $20$ complex solutions. We only keep the real solutions $X$ s.t. $a\in(0,1),ac-b^2>0$.

We give $2$ examples. In the first case, there is a maximum. In the second case, there is a non-reached upper bound.

a) $P=\begin{pmatrix}14&9\\9&7\end{pmatrix},Q=\begin{pmatrix}19&-9\\-9&10\end{pmatrix}$.

$\bullet$ We find an equation in $b$ of degree $16$ which has $14$ real roots. Only one solution satisfies the conditions $a\in(0,1),ac-b^2>0$.

We find $a\approx 0.47831,b\approx -0.33816, \max\approx 0.2391$.

$\bullet$ Using a Maple software of optimization, under the conditions $a\in [0,1],ac-b^2= 0$, we obtain -numerically- an approximation of the upper bound on the edge

$a\approx 0.451226,b\approx -0.497615,\sup\approx 0.22375<\max$.

b) $P=diag(1,2),Q=\begin{pmatrix}2&1\\1&3\end{pmatrix}$.

$\bullet$ We find an equation in $b$ of degree $15$ which has $15$ real solutions. No solutions satisfy the conditions $a\in (0,1),ac-b^2>0$.

$\bullet$ Using a Maple software of optimization, under the conditions $a\in [0,1],ac-b^2\geq 0$, we obtain -numerically- an approximation of the upper bound on the edge.

$a\approx 0.9898182,b\approx -0.1003899,\sup\approx 0.8447258$.

Note that $ac-b^2\approx 0$ but $a\not=1,b\not=0$; when $a\rightarrow 1^-,b=0$, $f(X)$ tends to $\dfrac{5}{6}\approx 0.83333$.

$\textbf{Remark.}$ When $n=2$ it suffices to randomly choose several initial points $X_0\in\{Y>0,Trace(Y)=1\}$ and that works after some tests. Since the calculations are unstable, we must remove the obtained solutions s.t. $ac-b^2<0$, even if $ac-b^2=-10^{-20}$.

vi) $\textbf{The case $n\geq 3$.}$ We only consider the case $n=3$; of course, when $n$ increases, the complexity of the calculations increases very quickly. There is a non-reached upper bound on the edge. Here, unfortunately, $X$ -a priori- doesn't satisfy some standard matricial equation.

Let $\Delta_i$ be the principal minor of dimension $i$ of $X$. Using a Maple software of optimization, under the conditions

$(x_{i,i})_{1\leq i\leq n}\geq 0,(\Delta_i)_{2\leq i\leq n}\geq 0,tr(X)=1$, we obtain -numerically- an approximation of the upper bound on the edge. The above remark also applies in these cases: a solution with $\Delta_2\approx -10^{-20}$ is forbidden. Even when $n=3$ we need quite a few tests to get an acceptable approximation of the $\sup$.

An example. $P=\begin{pmatrix}43&47&7\\47&57&7\\7&7&29\end{pmatrix},Q=\begin{pmatrix}50&0&-10\\0&11&5\\-10&5&19\end{pmatrix}$.

EDIT. We find $\sup\approx 0.180$ with $X$ neighbor of $\begin{pmatrix}0.571&*&*\\-0.420&0.374&*\\0.0686&0.00117&0.0548\end{pmatrix}$.