Eigenvectors that are tensor products?

Question

Consider a fixed $N\times N$ positive definite symmetric matrix $A$. Assume $N=d^r$ for some $d,r\geq 1$.

I wonder if one can find a closed formula for the maximizer/maximum of the function $$f(x):=\langle x^{\otimes r}, A x^{\otimes r}\rangle$$ over $x\in\mathbb R^d$ under the constraint $\|x\|=1$ (unit sphere).

Of course when $r=1$ the maximum is reached at the eigenvector associated with the maximal eigenvalue of $A$, which is the maximal value of $f$. But what if $r\geq 2$?

When $r\geq 2$, this is still true when $A$ is the $r$-fold product of a $d\times d$ matrix $B$, namely $A=B^{\otimes r}$.

So my question is: What can we say when $r\geq 2$ and $A$ not the $r$-fold product of a $d\times d$ matrix? Somehow we are looking at eigenvectors of $A$ that are tensor products, hence the title.

Answer 1

Your $f(x) := \langle x^{\otimes r}, A x^{\otimes r}\rangle$ is a homogeneous polynomial of degree $2r$ of $d$ variables, and $q(x) := ||x||^2$ is a homogeneous quadratic polynomial, so you are looking for extrema of a rational function $g(x) := \frac{f(x)}{q(x)^r}$ under the assumption that $q$ is everywhere positive.

In fact, the linear map $S^2 \big(V^{\otimes r}\big)^* \to S^{2r} V^*$ is surjective, and polynomials in the image of the positive-definite cone are sums of squares of degree $r$ polynomials, known as polynomials SOS.

The critical points of $g$ are given by $d\log(f) = r\cdot d\log(q)$, i.e. $d$ equations of $d$ variables of degree $2r+1$ each, $\frac{\partial{f}}{\partial x_i} \cdot (x_1^2 + \dots + x_d^2) = 2r \cdot x_i \cdot f(x)$, and I'm afraid no formula is much "closer" than this.

Geometrically you can also think of the function $g(x)$ as a pencil of degree $2r$ hypersurfaces in $(d-1)$-dimensional projective space $\mathbf{P}(V)$, such that the fiber over $\infty$ is non-reduced, a quadric $q(x)=0$ with multiplicity $r$, and the fiber over $0$ if $f(x)=0$, and then you are looking for the biggest $R$ such that fiber over $R$ is singular and has at least one real singular point. For generic $f(x)$ the function $g(x)$ is Morse, so the number $D(d,r)$ of its complex critical values can be computed topologically and explicitly, from Euler characteristics. Your number $R = max g(x)$ is a root of the discriminantal polynomial $P(z) \in \mathbf{R}[z]$ of one variable and degree $D(d,r)$. For generic $f(x)$ (and $d,r$) this polynomial $P$ will be irreducible, so in some sense the closest formula is saying that $R$ is the biggest real root of $P$.

Maybe there exists some more practical non-algebraic formula, like a continued fraction, based on some iterative approximations.

Answer 2

I strongly doubt any efficient algorithm exists for even approximately finding the maximizer in general. With slightly different notation, your problem is the same as the problem of "tensor principal component analysis" (tensor PCA). See the definition on page 2 of https://arxiv.org/pdf/1411.1076.pdf There, a reference is given that an exact solution is NP-hard. Further, that paper considers a statistical model in which the tensor is a sum of a rank-1 symmetric tensor (i.e., a tensor which is of the form $v^{\otimes k}$ for some vector $v$ where $k$ is the number of indices of the tensor) plus a random tensor with entries that are i.i.d. Gaussian up to the requirement of symmetry. In that case, at low noise (i.e., the Gaussian random tensor has small variance), there are polynomial algorithms which (with high probability) allow one to recover the vector $v$ up to some error, given the tensor as input. The simplest of these algorithms uses an eigendecomposition of a matrix, and the eigenvector with leading eigenvalue is close to $v^{\otimes k/2}$. As the noise increases, it is conjectured that there is a regime where no polynomial time algorithm can approximately recover $v$, even though approximate recovery is information-theoretically possible. At even larger noise, recovery becomes information-theoretically impossible.

Answer 3

This is NP-hard already when $r = 2$. To see this, I will consider the problem of minimizing your function $f$ instead of maximizing it, but it's not too much work to flip things around to see that maximization is hard too. I'll also ignore the positive semidefiniteness constraint and instead just assume that $A$ is symmetric, since we can freely add multiples of the identity to $A$, thus making $A$ PSD if we like without changing the difficulty of maximizing/minimizing.

With that out of the way, let $C$ be a $d \times d$ symmetric matrix, and define $A$ (which is $d^2 \times d^2$) by $$ A = \sum_{i,j=1}^d c_{i,j}\mathbf{e_ie_j^T} \otimes \mathbf{e_ie_j^T}, $$ where $\mathbf{e_i}$ is the $i$-th standard basis vector of $\mathbb{R}^d$. A calculation shows that $$ f(x) = \langle \mathbf{y}, C\mathbf{y}\rangle, $$ where $\mathbf{y} = (x_1^2,\ldots,x_d^2)$. It follows that the minimum value of $f$ is non-negative if and only if $C$ is copositive. Since determining copositivity of a matrix is NP-hard, so is your problem.