1 Introduction↩︎

Given a tuple \(\mathbf{x}=(x_1,\dots,x_n)\), we denote by \(\mathbb{R}[\mathbf{x}]\) the vector space of polynomials with real coefficients. Given \(r \in \mathbb{N}\) and a finite set of polynomials \(q, g_1,\dots,g_r \subset \mathbb{R}[\mathbf{x}]\), we consider the problem of minimizing \(q\) over the basic compact semialgebraic set \[\mathbf{X}:=\{\mathbf{x}\in \mathbb{R}[\mathbf{x}] : g_1(\mathbf{x}) \geq 0, \dots, g_r(\mathbf{x}) \geq 0 \},\] namely to compute \(q_{\min} := \inf_{\mathbf{x}\in \mathbf{X}} q(\mathbf{x})\). Note that this polynomial optimization problem is equivalent to \[\begin{align} \label{eq:pop} q_{\min} = \sup_{\tau \in \mathbb{R}} \{\tau : q - \tau \in \mathcal{P}_+(\mathbf{X}) \}, \end{align}\tag{1}\] where \(\mathcal{P}_+(\mathbf{X})\) is the set of polynomials being nonnegative on \(\mathbf{X}\). Polynomial optimization problems of the form 1 capture difficult instances such as the NP-hard MaxCut problem. We refer the interested reader to the recent monographs [1]–[3] that address several theoretical and applicative aspects of polynomial optimization.

Since the characterization of \(\mathcal{P}_+(\mathbf{X})\) is a highly non-trivial task, there have been substantial research efforts to provide tractable inner approximations, based on sums of squares decompositions. A polynomial \(p \in \mathbb{R}[\mathbf{x}]\) is called a sum of squares if it can be decomposed as \(p= \sum_{i} q_{i}^2\), for some \(q_{i} \in \mathbb{R}[\mathbf{x}]\). The set of sums of squares of polynomials is denoted by \(\Sigma[\mathbf{x}]\).

Let us define \(g_0:=1\), \([r] := \{1,\dots,r\}\) and \(g_J:= \displaystyle\prod_{j \in J} g_j\), for all \(J \subseteq [r]\).

For instance one can approximate \(\mathcal{P}_+(\mathbf{X})\) from the inside while relying either on the preordering \(\mathcal{T}(\mathbf{X}) \subseteq \mathcal{P}_+(\mathbf{X})\) or the \(\mathcal{Q}(\mathbf{X}) \subseteq \mathcal{P}_+(\mathbf{X})\), defined as \[\begin{align} \label{eq:preordering} \mathcal{T}(\mathbf{X}) & := \left\lbrace \sum_{J \in [r]} \sigma_J g_J : \sigma_J \in \Sigma[\mathbf{x}] \right\rbrace, \\ \mathcal{Q}(\mathbf{X}) & := \left\lbrace \sum_{j=0}^m \sigma_j g_j : \sigma_j \in \Sigma[\mathbf{x}] \right\rbrace. \end{align}\tag{2}\] When \(\mathbf{X}\) is compact, Schmüdgen’s Positivstellensatz states that every polynomial positive on \(\mathbf{X}\) belongs to \(\mathcal{T}(\mathbf{X})\):

Putinar’s Positivstellensatz allows one to replace \(\mathcal{T}(\mathbf{X})\) by \(\mathcal{Q}(\mathbf{X})\) in the above result under an assumption slightly stronger than compactness. Given a vector \(\mathbf{x}\in \mathbb{R}^n\), its Euclidean norm is denoted by \(\|\mathbf{x}\|:= \sqrt{\sum_{i=1}^n x_i^2}\). The quadratic module \(\mathcal{Q}(\mathbf{X})\) is called Archimedean if \(N - \|\mathbf{x}\|^2 \in \mathcal{Q}(\mathbf{X})\), for some \(N >0\).

1.1 Hierarchies of lower bounds for polynomial optimization↩︎

For a given integer \(t \geq \max_{1 \leq j \leq r} \{ \lceil \deg g_j / 2 \rceil \}\), if one considers the truncated preordering or quadratic module defined as \[\begin{align} \label{eq:truncated} \mathcal{T}(\mathbf{X})_{2t} & := \left\lbrace \sum_{J \in [r]} \sigma_J g_J : \sigma_J \in \Sigma[\mathbf{x}], \quad \deg (\sigma_J g_J) \leq 2t \right\rbrace,\\ \mathcal{Q}(\mathbf{X})_{2t} & := \left\lbrace \sum_{j=0}^m \sigma_j g_j : \sigma_j \in \Sigma[\mathbf{x}], \quad \deg (\sigma_j g_j) \leq 2t \right\rbrace, \end{align}\tag{3}\] then checking that a polynomial belongs to either \(\mathcal{T}(\mathbf{X})_{2t}\) or \(\mathcal{Q}(\mathbf{X})_{2t}\) can be done by solving a semidefinite program. Semidefinite programming is a particular class of convex programming, that consists of minimizing a linear function under linear matrix inequality constraints [6], and for which dedicated solvers can be readily used.

The hierarchies of Schmüdgen-type and Putinar-type lower bounds defined as \[\begin{align} \tag{4} \mathop{\mathrm{lb}}(q, \mathcal{T}(\mathbf{X}))_t := \sup_{\tau \in \mathbb{R}} \{\tau : q - \tau \in \mathcal{T}(\mathbf{X})_{2t} \}, \\ \mathop{\mathrm{lb}}(q, \mathcal{Q}(\mathbf{X}))_t := \sup_{\tau \in \mathbb{R}} \{\tau : q - \tau \in \mathcal{Q}(\mathbf{X})_{2t} \}, \tag{5} \end{align}\] both converge to \(q_{\min}\) when \(t\) goes to infinity, as a consequence of Schmüdgen’s and Putinar’s Positivstellensatz, respectively. The convergence proof for Putinar-type bounds is provided in [7]. Let us warn the reader that the lower bounds are indexed by the subscript \(t\) whenever the truncation order of the preordering or quadratic module is equal to \(2 t\). When providing convergence rates for the bi-sphere later on, the lower bounds will be indexed by \(2t\) and the truncation order of the preordering or quadratic module will be equal to \(4 t\). Note that when \(\mathbf{X}\) involves either equality constraints only or a single inequality constraint, then the preordering and quadratic module coincide. However, only the former remains true when one considers Cartesian products.

1.2 Related works on convergence rates↩︎

Recently, there has been a lot of dedicated research efforts in analyzing the asymptotic behavior of the hierarchies of Schmüdgen-type and Putinar-type bounds defined in 4 5 , respectively. For a polynomial \(q\) achieving its minimal value \(q_{\min}\) over a general compact semialgebraic set \(\mathbf{X}\), it was proved in [8] that the Schmüdgen-type bounds \(\mathop{\mathrm{lb}}(q, \mathcal{T}(\mathbf{X}))_t\) converge to \(q_{\min}\) at a rate in \(O(1/t^c)\), where \(c\) is a constant depending on \(\mathbf{X}\). If \(\mathcal{Q}(\mathbf{X})\) is Archimedean, an initial result provided in [9] shows that the Putinar-type bounds \(\mathop{\mathrm{lb}}(q, \mathcal{Q}(\mathbf{X}))_t\) converge to \(q_{\min}\) at a rate in \(O(1/\log(t)^c)\). Recently, this rate has been exponentially improved in [10], where the authors prove that the Putinar-type bounds converge to \(q_{\min}\) at a rate in \(O(1/t^c)\), thus matching the convergence rate of Schmüdgen-type bounds. Without any compactness assumption on the feasible set, one can rely on Putinar-type hierarchies based on sums of squares of rational fractions with (prescribed) uniform denominators. For these hierarchies, a similar rate in \(O(1/t^c)\) has been obtained in [11].

For special cases of feasible sets including the unit sphere/ball, the standard simplex and the hypercube, the value of the constant \(c\) can be explicitly given. An improved convergence rate in \(O(1/t^2)\) for Schmüdgen-type bounds has been obtained in [12] for the unit sphere, in [13] for the unit ball and the standard simplex, and in [14] for the hypercube. The case of the binary hypercube \(\{0,1\}^n\) is covered in [15]. For Putinar-type bounds, a rate in \(O(1/t)\) has been recently obtained in [16] for the hypercube. For polynomial optimization with correlative sparsity, stronger rates have been obtained in [17] under the assumption that the input data are (sufficiently) sparse.

In addition, convergence rates are available for the generalized moment problem (GMP), which consists of minimizing a linear function over the cone of positive Borel measures supported on \(\mathbb{R}^n\). Note that polynomial optimization can be cast as a special instance of the GMP, and that Schmüdgen-type and Putinar-type hierarchies can be naturally derived for other GMP instances. When the GMP instance involves finitely many constraints, one can obtain a rate similar to the one for polynomial optimization on the measure support [18]. Specialized rates have been obtained for certain GMP instances involving infinitely (countably) many constraints, in particular for optimal control [19], volume estimation [20], and exit location estimation of stochastic processes [21].

We refer to the recent survey [22] for a more detailed overview of convergence rates and related aspects, in particular the tightness of the performance analysis, the exponential convergence under local optimality conditions, and finite convergence.

1.3 Our contributions↩︎

In this work, we mainly focus on hierarchies of Schmüdgen-type bounds for polynomial optimization over Cartesian products of compact semialgebraic sets.

We first provide in Section 3.1 a convergence rate of Schmüdgen-type, or equivalently Putinar-type, bounds when minimizing a polynomial on the product of spheres \(\mathbf{X}= S^{n-1} \times S^{n-1}\). Given two \(n\)-dimensional tuples of real variables \(\mathbf{x}=(x_1,\dots,x_n)\) and \(\mathbf{y}=(y_1,\dots,y_n)\), and a polynomial \(q \in \mathbb{R}[\mathbf{x},\mathbf{y}]\), we define \(q_{\min}:= \min_{(\mathbf{x},\mathbf{y}) \in \mathbf{X}} q(\mathbf{x},\mathbf{y}) = \min_{\mathbf{x},\mathbf{y}\in S^{n-1}} q(\mathbf{x},\mathbf{y})\) and \(q_{\max} := \max_{(\mathbf{x},\mathbf{y}) \in \mathbf{X}} q(\mathbf{x},\mathbf{y})\).

The proof of Theorem 3 is given in Section 3.1.

We extend this result in Section 3.2 (Theorem 7) to the case of arbitrary sphere products, and in Section 3.3 (Theorem 8) to the case of products of distinct sets. As proved in [8], the hierarchy of Schmüdgen-type relaxations has the convergence rate \(O(1/t^c)\) for a general compact set \(\mathbf{X}\), where the constant \(c\) depends on \(\mathbf{X}\). We emphasize that the result from [8] only yields the existence of such a constant \(c\) but cannot directly provide the value, e.g., \(c=2\) for specific constraint sets. The drawback of [8] is that \(c\) depends on the description of the constraint set in an unspecified way. For any concrete situation, such as the sphere, the unit ball, the unit cube or the standard simplex, one could hope to extract a suitable \(c\) from the proof of Section 3 from [8], see [8] for more explanation. It would in particular require to estimate intermediate constants satisfying the Łojasiewicz inequality used in [8], which appears to be delicate. Theorem 7 implies that one can select \(c=2\) when \(\mathbf{X}\) is the product of two or more spheres.

Eventually we analyze in Section 4 the convergence rate of a problem arising from quantum information theory, called the quantum Wasserstein distance. Computing this distance boils down to solving a specific instance of the generalized moment problem on the complex bi-sphere. We provide a rate in \(O(1/t^2)\), stated in Theorem 9, by first changing variables from complex to real, then by combining Theorem 3 with the recent convergence rate results from [18], obtained for the generalized moment problem.

After the initial submission of this work in May 2025, we became aware of the related work [23], that was submitted as an arXiv preprint two months later in July 2025. In [23] (see also Corollary 3.4), the authors provide similar convergence rates with a proof also based on polynomial kernel products. The authors of [23] also provide convergence rate results for the Hausdorff distance between (truncated) pseudo-moment and moment sequences; see Section 3.2 and Section 4 therein. We emphasize that the two works were carried out independently.

1.4 Outline of the proof technique↩︎

For a given compact semialgebraic set \(\mathbf{X}\subset \mathbb{R}^{n}\) and a polynomial \(p \in \mathbb{R}[\mathbf{x}]\), we denote the supremum-norm of \(p\) on \(\mathbf{X}\) by \(\|p\|_{\mathbf{X}}:= \max_{\mathbf{x}\in \mathbf{X}} |p(\mathbf{x})|\). To obtain the convergence rate of Theorem 3, we follow the exact same technical strategy pursued in [12], [13] for the case of the unit sphere \(S^{n-1}\), the unit ball and the standard simplex, as well as in [14] for the case of the hypercube. Let us denote by \(\mathbb{R}[\mathbf{x},\mathbf{y}]_d\) the vector space of polynomials of degree at most \(d\). Similarly let \(\mathcal{P}_+( \mathbf{X})_d\) be the subset of polynomials from \(\mathbb{R}[\mathbf{x},\mathbf{y}]_d\) that are nonnegative on \(\mathbf{X}\). Let \(q \in \mathbb{R}[\mathbf{x},\mathbf{y}]_d\) and \(\mathbf{X}= S^{n-1} \times S^{n-1}\). The goal is to prove that \[{q - q_{\min} + \varepsilon \in \mathcal{Q}(\mathbf{X})_{4t} (=\mathcal{T}(\mathbf{X})_{4 t}) }\] for some small \(\varepsilon\). Up to translation and scaling, we assume that \(q_{\min}=0\) and \(q_{\max}=1\). We will construct an invertible linear operator \(\mathop{\mathrm{\mathbf{K}}}: \mathbb{R}[\mathbf{x},\mathbf{y}]_d \to \mathbb{R}[\mathbf{x},\mathbf{y}]_d\) which satisfies the following properties: \[\begin{align} \tag{6} \mathop{\mathrm{\mathbf{K}}}(1)=1, \\ \tag{7} \mathop{\mathrm{\mathbf{K}}}p \in {\mathcal{Q}(\mathbf{X})}_{4 t}, \quad \text{for all } p \in \mathcal{P}_+( \mathbf{X})_d,\\ \tag{8} \| \mathop{\mathrm{\mathbf{K}}}^{-1} q - q \|_{\mathbf{X}} \leq \varepsilon. \end{align}\] Then it follows from [13] that \(q + \varepsilon \in {\mathcal{Q}(\mathbf{X})}_{4 t}\) and thus \(q_{\min} - \mathop{\mathrm{lb}}(q, {\mathcal{Q}(\mathbf{X})})_{2t} \leq \varepsilon\). As in [12]–[14], the statement of Theorem 3 is proven by showing the existence, for each large enough \(t\), of an operator \(\mathop{\mathrm{\mathbf{K}}}\) satisfying 6 , 7 and 8 with \(\varepsilon = O(1/t^2)\).

The way to construct such an operator is to follow the so-called polynomial kernel method. In Section 2, we explain how such a kernel is found in [12] for the case of the unit sphere. Our rate for the bi-sphere is obtained by multiplying the two kernels associated to each unit sphere. This idea of kernel multiplication is the same as the one previously used for the hypercube in [14] (see also [17]).

2 Preliminaries: convergence rate on the unit sphere↩︎

For a multi-index \(\alpha = (\alpha_1,\dots,\alpha_n) \in \mathbb{N}^n\), we use the notation \(|\alpha|=\alpha_1 + \dots + \alpha_n\). Given \(d \in \mathbb{N}\), let \(\mathbb{N}_d^n := \{\alpha \in \mathbb{N}^n : |\alpha| \leq d \}\). The dimension of \(\mathbb{N}_d^n\) is equal to the number of monomials of degree at most \(d\), i.e., \(\binom{n+d}{d}\).
A map \(\mathop{\mathrm{C}}: S^{n-1} \times S^{n-1} : \mathbb{R}\) on \(S^{n-1}\) is called a polynomial kernel when \(\mathop{\mathrm{C}}(\mathbf{x},\mathbf{x}')\) is a polynomial in the variables \(\mathbf{x}\), \(\mathbf{x}'\). For each \(k\), let \(H_k[\mathbf{x}] := \mathop{\mathrm{span}}\{ P_{\alpha} : \alpha \in \mathbb{N}^n, |\alpha|=k \}\) be the spanning set of orthonormal spherical harmonics of degree \(k\), i.e., homogeneous polynomials of degree \(k\) lying in the kernel of the Laplace operator. Let us equip \(S^{n-1}\) with the uniform Haar measure \(\mu\), and consider the Christoffel-Darboux kernel \(\mathop{\mathrm{C}}_{2 t}\), which is defined in terms of the orthonormal basis \(\{P_{\alpha} : \alpha \in \mathbb{N}^n\}\) for \(\mathbb{R}[\mathbf{x}]\) with respect to \((S^{n-1},\mu)\) as: \[\begin{align} \label{eq:cdksphere} \mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}') := \sum_{k=0}^{2t} \mathop{\mathrm{C}}^{(k)} (\mathbf{x},\mathbf{x}'), \quad \text{where } \mathop{\mathrm{C}}^{(k)} (\mathbf{x},\mathbf{x}') := \sum_{|\alpha|=k} P_{\alpha}(\mathbf{x}) P_{\alpha}(\mathbf{x}'). \end{align}\tag{9}\] Then one associates to \(\mathop{\mathrm{C}}_{2t}\) an operator \(\mathop{\mathrm{\mathbf{C}}}_{2t}: \mathbb{R}[\mathbf{x}] \to \mathbb{R}[\mathbf{x}]\), defined as follows: \[\begin{align} \label{eq:ksphere} \mathop{\mathrm{\mathbf{C}}}_{2t} p(\mathbf{x}) := \int_{S^{n-1}} \mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}') p(\mathbf{x}') \mathop{\mathrm{d}}\mu(\mathbf{x}') . \end{align}\tag{10}\] This operator \(\mathop{\mathrm{\mathbf{C}}}_{2t}\) is called reproducing as it satisfies \(\mathop{\mathrm{\mathbf{C}}}_{2t} p = p\), for all \(p \in \mathbb{R}[\mathbf{x}]_{2t}\). Therefore, all its eigenvalues are equal to \(1\).

One has \(p(\mathbf{x}) = \sum_{k=0}^d p_k(\mathbf{x})\), for all \(\mathbf{x}\in S^{n-1}\), with \(p_k \in H_k\) given by \[p_k = \int_{S^{n-1}} \mathop{\mathrm{C}}^{(k)} (\mathbf{x},\mathbf{x}') p(\mathbf{x}') \mathop{\mathrm{d}}\mu(\mathbf{x}') .\]

In this case one has the following closed form of the Christoffel-Darboux kernel \[\begin{align} \label{eq:funkhecke} \mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}') = \sum_{k=0}^{2t} \mathcal{G}_k^{(n-1)} (\mathbf{x}\cdot \mathbf{x}'), \quad \text{for all } \mathbf{x},\mathbf{x}' \in S^{n-1}. \end{align}\tag{11}\] In the above formula, the set \(\{\mathcal{G}_k^{(n)} : k \in \mathbb{N}\}\) of Gegenbauer polynomials are defined, for all \(n \geq 2\), as the set of orthogonal polynomials with respect to the weight function \(w_n(x):= \hat{w}_n (1-x^2)^{\frac{n-2}{2}}\), where \(\hat{w}_n\) is a normalization constant ensuring that \(\int_{-1}^1 w_n(x) \mathop{\mathrm{d}}x = 1\). We refer the interested reader to [24] for a detailed survey about univariate orthogonal polynomials. The formula 11 is also known as the Funk-Hecke formula, stated for instance in [25].

The next step is to perturb the eigenvalues of the reproducing operator \(\mathop{\mathrm{\mathbf{C}}}_{2t}\) by considering the following kernel instead: \[\begin{align} \label{eq:cdksphereperturb} \mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}'; \lambda) := \sum_{k=0}^{2t} \lambda_k \mathop{\mathrm{C}}^{(k)} (\mathbf{x},\mathbf{x}') = \sum_{k=0}^{2t} \lambda_k \mathcal{G}_k^{(n-1)} (\mathbf{x}\cdot \mathbf{x}'), \quad \text{for all } \mathbf{x},\mathbf{x}' \in S^{n-1}. \end{align}\tag{12}\] Then \(\lambda = (\lambda_k)\) is the set of eigenvalues of the operator obtained in 10 by considering \(\mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}'; \lambda)\) instead of \(\mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}')\). The following result provides one way to prove that \(\mathop{\mathrm{\mathbf{C}}}_{2t}\) satisfies 7 .

The proof of Lemma 1 relies on Tchakaloff’s Theorem [26] and the exisence of cubature rules. Now coming back to the case of the unit sphere, assume that \(\lambda = (\lambda_k)\) is the set of coefficients of a univariate sum of squares in the basis of Gegenbauer polynomials. Then at each fixed \(\mathbf{x}'\), the polynomial \(\mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}'; \lambda)\), defined in 12 , is equal to a sum of squares of degree at most \(2t\) on \(\mathbf{X}\), so \(\mathop{\mathrm{C}}_{2t}(\mathbf{x},\mathbf{x}'; \lambda) \in {\mathcal{Q}(\mathbf{X})_{2t}}\) and Lemma 1 implies that \(\mathop{\mathrm{\mathbf{C}}}_{2t}\) satisfies 7 .

Next, if the harmonic decomposition of \(q\) is \(q(\mathbf{x})=\sum_{k=0}^d q_k(\mathbf{x})\) for all \(\mathbf{x}\in S^{n-1}\), with \(q_k \in H_k[\mathbf{x}]\), then one has \(\mathop{\mathrm{\mathbf{C}}}_{2t} q = \sum_{k=0}^d \lambda_k q_k\). If in addition, \(\lambda_0=1\) and \(\lambda_k \neq 0\) for all \(k\) then \(\mathop{\mathrm{\mathbf{C}}}_{2t}^{-1} q = \sum_{k=0}^d 1/\lambda_k q_k\). By [13], the resulting perturbed operator satisfies 6 if \(\lambda_0=1\), and it satisfies 8 with \[\varepsilon = \max_{1\leq k \leq d} \max_{\|\mathbf{x}\|=1} |q_k(\mathbf{x})| \cdot \left|1 - 1/\lambda_k \right|.\] It remains to give an upper bound of the two multiplied quantities. For the first quantity, one can define the so-called harmonic constant: \[\begin{align} \label{eq:harmonic} \gamma(S^{n-1})_d := \max_{p \in \mathbb{R}[\mathbf{x}]_d} \max_{0\leq k \leq d} \frac{\|p_k\|_{S^{n-1}}}{\|p\|_{S^{n-1}}}. \end{align}\tag{13}\] In case of minimizing a homogeneous polynomial of degree \(d\), one can define a similar constant denoted by \(\gamma(S^{n-1})_{=d}\). This latter constant can be bounded independently of the dimension \(n\) as a consequence of [12]. More generally, one has the following result.

To finally obtain the desired convergence rate, the second quantity can be bounded as follows.

3 Convergence rate for polynomial optimization on set products↩︎

3.1 Products of two spheres↩︎

This section is dedicated to the proof of Theorem 3.

3.2 Extension to arbitrary sphere products↩︎

One can easily extend the rate given in Theorem 3 when optimizing a polynomial over the Cartesian product of \(m\) spheres \(\mathbf{X}= S^{n-1} \underbrace{\times \dots \times}_{m-1} S^{n-1}\).

3.3 Extension to products of distinct sets↩︎

Given a compact set \(\mathbf{X}\subset \mathbb{R}^n\), we denote by \(\mathcal{M}(\mathbf{X})\) the space of positive Borel measures supported on \(\mathbf{X}\). In this section, we assume that convergence rates have been obtained for polynomial optimization over two given sets \(\mathbf{X}_1\) and \(\mathbf{X}_2\), thanks to perturbed kernels and associated operators satisfying 6 , 7 and 8 . We show how to derive a convergence rate for polynomial optimization over the Cartesian product \(\mathbf{X}= \mathbf{X}_1 \times \mathbf{X}_2\).

Before stating our result, we briefly describe the quantities that will be involved.
For each set \(\mathbf{X}_i\), one associates an integer \(m_i\) that is the number of indices needed when decomposing polynomials into orthonormal components. As explained after Theorem 8, this number is equal to \(1\) when \(\mathbf{X}_i\) is either the unit sphere, the unit ball or the standard simplex, and is equal to \(n\) when \(\mathbf{X}_i\) is the hypercube \([-1,1]^n\). Its value differs if \(\mathbf{X}_i\) is itself a product set.
Next, one associates to each set \(\mathbf{X}_i\) a function \(d_i : \mathbb{N}\to \mathbb{N}\) that takes as input, an integer \(t\), and outputs the relaxation order \(d_i(t)\) used for the lower bound \(\mathop{\mathrm{lb}}(q, \mathcal{T}(\mathbf{X}_i))_{d_i(t)}\). In the case of the unit sphere, \(d_i\) is the identity, but for the bi-sphere, it is twice the identity.
Eventually, one associates to each set \(\mathbf{X}_i\) a function \(\eta^{(i)} : \mathbb{N}\to \mathbb{N}\to \mathbb{R}_{\geq 0}\) converging to 0 at infinity. This function is used to quantify the convergence rate. Given a relaxation order \(t\) as input, the output is proportional to \(1/t^2\) for the unit sphere, the unit ball, the standard simplex and the hypercube, tough the constant of proportionality depends on the set.

We now give examples of sets where the assumptions of Theorem 8 are fulfilled.

3.3.0.1 The unit sphere [12], the unit ball [13] or the standard simplex[13]

 
If \(\mathbf{X}_i\) is either \(S^{n-1} := \{\mathbf{x}\in \mathbb{R}^n : \|\mathbf{x}\|^2 = 1\}\), or \(\mathbb{B}^n := \{\mathbf{x}\in \mathbb{R}^n : \|\mathbf{x}\|^2 \leq 1\}\) or \(\Delta^n : \{\mathbf{x}\in \mathbb{R}^n : x_1 \geq 0, \dots, x_n \geq 0, \;1 - \sum_{i=1}^n x_i \geq 0\}\) then one can select

• \(m_i = 1\), \(t_i=2 n d \sqrt{d}\);

• \(d_i : = t \mapsto 2 t\);

• \(\mu_i\) to be either the uniform Haar measure for \(S^{n-1}\), or the measure with density \(c_n (1 - \|\mathbf{x}\|^2)^{-\frac{1}{2}}\) (with \(c_n> 0\) being a normalization constant) for \(\mathbb{B}^n\), or the measure with density \(c_n x_1^{-1/2} \cdots x_n^{-1/2} (1 - \sum_{i=1}^n x_i)^{-\frac{1}{2}}\) (with \(c_n> 0\) being a normalization constant) for \(\Delta^n\);

• \(\eta^{(i)} : t \mapsto \frac{n^2 d^3}{t^2}\);

• \(\mathcal{P}(S^{n-1}) = \displaystyle\bigoplus_{\kappa^{(i)}=0}^{\infty} H_{\kappa^{(i)}}\), where \(H_{\kappa^{(i)}}\) is spanned by the set of orthonormal polynomials with respect to \((S^{n-1},\mu_i)\) of exact degree \(\kappa^{(i)}\);

• the vector \(\lambda^{(i)}=\left(\lambda_{\kappa^{(i)}}^{(i)}\right)\), indexed by \(\mathbb{N}_{2t}\), for all \(t \geq t_i\), as in Lemma 3.

3.3.0.2 The hypercube [14]

 
If \(\mathbf{X}_i = [-1, 1]^n\), then one can select

• \(m_i = n\), \(t_i= \pi d \sqrt{2n}\);

• \(d_i : t \mapsto n (t+1)\);

• \(\mu_i\) to be the normalized Chebyshev measure with density \(\frac{1}{\pi \sqrt{1-x_1^2}} \cdots \frac{1}{\pi \sqrt{1-x_n^2}}\);

• \(\eta^{(i)} : t \mapsto \frac{n \pi^2 d^2}{t^2}\),

• \(\mathcal{P}([-1,1]^n) = \displaystyle\bigoplus_{\kappa^{(i)} \in \mathbb{N}^n} H_{\kappa^{(i)}}\), where \(H_{\kappa^{(i)}} = \mathbb{R}T_{\kappa^{(i)}}\) and \(T_{\kappa^{(i)}}\) is the multivariate Chebyshev polynomial defined as the product of univariate Chebyshev polynomials: \[T_{\kappa^{(i)}}(\mathbf{x}) := \prod_{j=1}^n T_{\kappa_j^{(i)}}(x_j), \quad T_k(\cos \phi) := \cos (k \phi) \quad (\phi \in \mathbb{R}, k \in \mathbb{N}) ;\]

• the vector \(\lambda^{(i)}=\left(\lambda_{\kappa^{(i)}}^{(i)}\right)\), indexed by \(\mathbb{N}^n_{2t}\), for all \(t \geq t_i\), as in [14].

From the above, we obtain a direct corollary of Theorem 8.

4 Application to quantum Wasserstein distance↩︎

In this section we analyze the convergence rate of the hierarchy derived in [28] to approximate as closely as desired the so-called order 2 quantum Wasserstein distance. Roughly speaking, the classical Wasserstein distance quantifies the cost of transporting one probability measure to another. In order to extend this notion to quantum states, several notions of quantum Wasserstein distances have been proposed; see, e.g., [29] and [30]. In [28], the authors specifically focus on the order 2 quantum Wasserstein distance that has been recently proposed in [31].

4.1 Problem statement and main result↩︎

4.1.0.1 Additional complex notation.

Let \(\mathbf{i}= \sqrt{-1} \in \mathbb{C}\). Given \(x \in \mathbb{C}\), one can write \(x = a + \mathbf{i}\, b\), where \(a = \mathop{\mathrm{Re}}(x)\) is the real part of \(x\) and \(b=\mathop{\mathrm{Im}}(x)\) is the imaginary part of \(x\). The complex conjugate of \(x\) is \(\overline{x} = a - \mathbf{i}\, b\). We extend this to vectors and matrices by applying the corresponding operations entrywise. The modulus of \(x\) is denoted by \(|x|:=\sqrt{a^2+b^2}\). The set of polynomials in variables \(\mathbf{x}=(x_1,\dots,x_n)\) and \(\overline{\mathbf{x}}=(\overline{x}_1,\dots,\overline{x}_n)\) with complex coefficients is denoted by \(\mathbb{C}[\mathbf{x},\overline{\mathbf{x}}]\). A polynomial \(p \in \mathbb{C}[\mathbf{x},\overline{\mathbf{x}}]\) satisfying \(p=\overline{p}\) is called Hermitian. Given a vector \(u\in \mathbb{C}^n\), \(u^*\) is the row vector containing the conjugate entries of \(u\). We will use the same notation as for the real case for the Euclidean norm \(\|u\| = \sqrt{u^* u}\). The \(\ell_1\)-norm \(\|u\|_1\) is the sum of moduli of the entries of \(u\). A linear functional \(L : \mathbb{C}[\mathbf{x},\overline{\mathbf{x}}] \to \mathbb{C}\) is called Hermitian if \(\overline{L(p)} = L(\overline{p})\) for all \(p \in \mathbb{C}[\mathbf{x},\overline{\mathbf{x}}]\). Given a compact set \(\mathbf{X}\subset \mathbb{C}^n\), we denote by \(\mathcal{M}(\mathbf{X})\) the space of complex positive Borel measures supported on \(\mathbf{X}\).

Given \(n \in \mathbb{N}\), let \(\rho, \nu \in \mathbb{C}^{n \times n}\) be two normalized quantum states, i.e., Hermitian matrices with nonnegative eigenvalues and unit trace. A quantum transport plan between \(\rho\) and \(\nu\) is a finite set of triples \(\{(\omega_{\ell}, u_{\ell}, v_{\ell})\}_{\ell}\) such that \[\begin{align} \sum_{\ell} \omega_{\ell} \, u_{\ell} u_{\ell}^* = \rho \quad \text{ and} \quad\sum_{\ell} \omega_{\ell} \, v_{\ell} v_{\ell}^* = \nu, \end{align}\] where \(\omega_{\ell} > 0\), \(\sum_{\ell} \omega_{\ell} =1\), \(u_{\ell}, v_{\ell} \in \mathbb{C}^n\), \(\|u_{\ell}\|=\| v_{\ell}\| = 1\). The set of all such quantum transport plans between two states \(\rho\) and \(\nu\) is denoted by \(\mathcal{Q}(\rho,\nu)\). Given a plan \(e=\{(\omega_{\ell}, u_{\ell}, v_{\ell})\}_{\ell}\), the associated (order \(2\)) quantum transport cost is defined as \[\begin{align} T_2(e) := \sum_{\ell} \omega_{\ell} \mathop{\mathrm{Tr}}\left[ (u_{\ell} u_{\ell}^* - v_{\ell} v_{\ell}^*) \overline{(u_{\ell} u_{\ell}^* - v_{\ell} v_{\ell}^*)} \right], \end{align}\] and the order \(2\) quantum Wasserstein distance is defined as \[\begin{align} \label{eq:qwtwo} W_2(\rho,\nu) := \Big(\inf_{e\in \mathcal{Q}(\rho,\nu)} T_2(e)\Big)^{1/2}. \end{align}\tag{27}\] In [32], the authors propose a moment reformulation of the DPS hierarchy [33], initially designed to distinguish separable and entangled states. Inspired by this reformulation, it is proved in [28] that \(W_2^2(\rho,\nu)\) is the optimal value of a particular instance of the generalized moment problem (GMP), i.e., an infinite-dimensional linear problem over positive Borel measures: \[\label{eq:measC} \begin{align} W^2_2(\rho,\sigma) = \inf_{\mu \in \mathcal{M}(\mathbf{X})} \quad & \int_{\mathbf{X}} f \mathop{\mathrm{d}}\mu \\ \mathop{\mathrm{s.t.}} \quad & \int_{\mathbf{X}} \mathbf{x}\mathbf{x}^* \mathop{\mathrm{d}}\mu = \rho, \\ \quad & \int_{\mathbf{X}} \mathbf{y}\mathbf{y}^* \mathop{\mathrm{d}}\mu = \nu, \end{align}\tag{28}\] where \(f(\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}})= \mathop{\mathrm{Tr}}[(\mathbf{x}\mathbf{x}^* - \mathbf{y}\mathbf{y}^*) \overline{(\mathbf{x}\mathbf{x}^* - \mathbf{y}\mathbf{y}^*)}]\) and \(\mathbf{X}= \{(\mathbf{x},\mathbf{y}) \in \mathbb{C}^{2n}: \|\mathbf{x}\|=\|\mathbf{y}\|=1 \}\) is the Cartesian product of two complex unit spheres. The dual of problem 28 is given by \[\label{eq:dualC} \begin{align} \sup_{\Lambda,\Gamma \in \mathbb{C}^{n \times n}} \quad & \mathop{\mathrm{Tr}}(\rho \Lambda + \nu \Gamma) \\ \text{s.t.} \quad & f - \mathbf{x}^* \Lambda \mathbf{x}- \mathbf{y}^* \Gamma \mathbf{y}\geq 0 \text{ on } \mathbf{X}, \\ \quad & \Lambda^*=\Lambda, \quad \Gamma^*=\Gamma. \end{align}\tag{29}\] The Putinar-type (or equivalently Schmüdgen-type) moment hierarchy for 28 is given by \[\label{eq:momC} \begin{align} W^2_2(\rho,\sigma)_t := \inf_{L} \quad & L(f) \\ \text{s.t.} \quad & L: \mathbb{C}[\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}]_{2t} \to \mathbb{C}\text{ Hermitian linear}, \\ \quad & L(\mathbf{x}\mathbf{x}^*) = \rho, \\ \quad & L(\mathbf{y}\mathbf{y}^*) = \nu, \\ \quad & L \geq 0 \text{ on } {\mathcal{Q}(\mathbf{X})}_{2t}. \end{align}\tag{30}\] It is proved in [28] that the corresponding sequence of lower bounds given by this hierarchy converges to the order \(2\) quantum Wasserstein distance, i.e., \(W^2_2(\rho,\sigma)_t \to W^2_2(\rho,\sigma)\) when \(t\) goes to infinity. In the recent contribution [18], the authors analyze convergence rates of hierarchies of lower bounds for the GMP. We recall this result in the case where the GMP instance involves a single measure and equality constraints only.

Our final result, given in Theorem 9 below, provides the convergence rate of the moment hierarchy 30 for the order \(2\) quantum Wasserstein distance, assuming that there is a bounded dual optimal solution of 29 . This result is obtained by combining Lemma 4 with our analysis derived in Theorem 3.

The proof of Theorem 9 is postponed in Section 4.3. As shown in [28], when both \(\rho\) and \(\nu\) are positive definite quantum states, then the dual 29 is attained at a bounded tuple. In this case, the constant \(\kappa\) depends on \(n\) and the reciprocal of the states’ minimal eigenvalues.

4.2 Changing variables from complex to real↩︎

In order to apply Lemma 4, we provide an equivalent formulation of the moment hierarchy 30 involving linear functionals acting on real polynomials. This reformulation is obtained thanks to [32], where the authors describe how to transform polynomials in \(\mathbb{C}[\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}]\) into real polynomials in \(\mathbb{R}[\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}]\) via the change of variables \(\mathbf{x}= \mathbf{a}+ \mathbf{i}\, \mathbf{b}\), \(\mathbf{y}= \mathbf{c}+ \mathbf{i}\, \mathbf{d}\). In this way, any polynomial \(p \in \mathbb{C}[\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}]\) corresponds to a unique pair of real polynomials \[\begin{align} p_{\mathop{\mathrm{Re}}} (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}) = \mathop{\mathrm{Re}}(p(\mathbf{a}+\mathbf{i}\, \mathbf{b}, \mathbf{a}-\mathbf{i}\, \mathbf{b}, \mathbf{c}+\mathbf{i}\, \mathbf{d}, \mathbf{c}-\mathbf{i}\, \mathbf{d})), \\ p_{\mathop{\mathrm{Im}}} (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}) = \mathop{\mathrm{Im}}(p(\mathbf{a}+\mathbf{i}\mathbf{b}, \mathbf{a}-\mathbf{i}\, \mathbf{b}, \mathbf{c}+\mathbf{i}\, \mathbf{d}, \mathbf{c}-\mathbf{i}\, \mathbf{d})), \end{align}\] satisfying \(p(\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}})= p_{\mathop{\mathrm{Re}}} (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}) + \mathbf{i}\, p_{\mathop{\mathrm{Im}}} (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d})\). When \(p\) is Hermitian, one has \(p(\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}) = p_{\mathop{\mathrm{Re}}} (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d})\), so in particular one has \[f(\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}})= \mathop{\mathrm{Tr}}[(\mathbf{x}\mathbf{x}^* - \mathbf{y}\mathbf{y}^*) \overline{(\mathbf{x}\mathbf{x}^* - \mathbf{y}\mathbf{y}^*)}] = f_{\mathop{\mathrm{Re}}} (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}) .\] Given a set of Hermitian polynomials, we define its real analog by applying the map \(\mathop{\mathrm{Re}}(\cdot)\) entrywise to this set. In particular, the real analog of the complex bi-sphere \(\mathbf{X}\) is the set \(\mathbf{X}_{\mathop{\mathrm{Re}}} = \{ (\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}) \in \mathbb{R}^{4n} : 1 - \|\mathbf{a}\|^2 - \|\mathbf{b}\|^2 = 1 - \|\mathbf{c}\|^2 - \|\mathbf{d}\|^2 = 0 \}\). The associated truncated real quadratic module is \[\begin{align} {\mathcal{Q}}^{\mathbb{R}}(\mathbf{X}_{\mathop{\mathrm{Re}}})_{2t} := \left\{ \sigma + q_1 (1 - \|\mathbf{a}\|^2 - \|\mathbf{b}\|^2) + q_2 (1 - \|\mathbf{c}\|^2 - \|\mathbf{d}\|^2) : \right. & \sigma \in \Sigma[\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}]_{2t}, \\ & \left. q_1,q_2 \in \mathbb{R}[\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}]_{2t-2} \right\} . \end{align}\] Similarly, the maps \(\mathop{\mathrm{Re}}(\cdot)\) and \(\mathop{\mathrm{Im}}(\cdot)\) act entrywise for vectors and matrices with polynomial entries in \(\mathbb{C}[\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}]\), e.g., \(\mathop{\mathrm{Re}}(\mathbf{x}\mathbf{x}^*) = \mathbf{a}\mathbf{a}^T + \mathbf{b}\mathbf{b}^T\) and \(\mathop{\mathrm{Im}}(\mathbf{x}\mathbf{x}^*) = \mathbf{b}\mathbf{a}^T - \mathbf{a}\mathbf{b}^T\). Eventually, for a linear functional \(L : \mathbb{C}[\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}] \to \mathbb{C}\), one has \(L(p) = \mathop{\mathrm{Re}}(L(p))+ \mathbf{i}\, \mathop{\mathrm{Im}}(L(p))\). If in addition \(L\) is Hermitian, that is, satisfies \(\overline{L(p)}=L(\overline{p})\), then one can associate to \(L\) a real linear functional \(L^{\mathbb{R}} : \mathbb{R}[\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}] \to \mathbb{R}\), defined by \[L^{\mathbb{R}} (q) = L\left( q \left( \frac{\mathbf{x}+\overline{\mathbf{x}}}{2}, \frac{\mathbf{x}-\overline{\mathbf{x}}}{2 \mathbf{i}}, \frac{\mathbf{y}+\overline{\mathbf{y}}}{2}, \frac{\mathbf{y}-\overline{\mathbf{y}}}{2 \mathbf{i}} \right) \right) .\] Then for any \(p \in \mathbb{C}[\mathbf{x},\overline{\mathbf{x}},\mathbf{y},\overline{\mathbf{y}}]\), one has \(L(p) = L^{\mathbb{R}} (p_{\mathop{\mathrm{Re}}}) + \mathbf{i}\, L^{\mathbb{R}} (p_{\mathop{\mathrm{Im}}})\). For any \(t \geq 2\), we can now provide the following relaxation over real linear functionals: \[\label{eq:momR} \begin{align} \inf_{L^{\mathbb{R}}} \quad & L^{\mathbb{R}}(f^{\mathop{\mathrm{Re}}}) \\ \text{s.t.} \quad & L^{\mathbb{R}}: \mathbb{R}[\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}]_{2t} \to \mathbb{R}\text{ linear}, \\ \quad & L^{\mathbb{R}}(\mathbf{a}\mathbf{a}^T + \mathbf{b}\mathbf{b}^T) = \mathop{\mathrm{Re}}(\rho), \\ \quad & L^{\mathbb{R}}(\mathbf{b}\mathbf{a}^T - \mathbf{a}\mathbf{b}^T) = \mathop{\mathrm{Im}}(\rho), \\ \quad & L^{\mathbb{R}}(\mathbf{c}\mathbf{c}^T + \mathbf{d}\mathbf{d}^T) = \mathop{\mathrm{Re}}(\nu), \\ \quad & L^{\mathbb{R}}(\mathbf{d}\mathbf{c}^T - \mathbf{c}\mathbf{d}^T) = \mathop{\mathrm{Im}}(\nu), \\ \quad & L^{\mathbb{R}} \geq 0 \text{ on } \mathcal{T}^{\mathbb{R}}(\mathbf{X}_{\mathop{\mathrm{Re}}})_{2t}. \end{align}\tag{31}\] As a consequence of the above defined transformations, the relaxation 31 is equivalent to the moment relaxation 30 over Hermitian linear functionals. In particular the equivalence \[L \geq 0 \text{ on } {\mathcal{Q}}(\mathbf{X})_{2t} \iff L^{\mathbb{R}} \geq 0 \text{ on } {\mathcal{Q}}^{\mathbb{R}}(\mathbf{X}_{\mathop{\mathrm{Re}}})_{2t}\] is proved in [32]. As in the complex setting, one can show that the hierarchy of relaxations converges to \(W^2_2(\rho,\nu)\) when \(t\) goes to infinity, that is \[\label{eq:infmomR} \begin{align} W^2_2(\rho,\nu) = \inf_{L^{\mathbb{R}}} \quad & L^{\mathbb{R}}(f_{\mathop{\mathrm{Re}}}) \\ \text{s.t.} \quad & L^{\mathbb{R}}: \mathbb{R}[\mathbf{a},\mathbf{b},\mathbf{c},\mathbf{d}] \to \mathbb{R}\text{ linear}, \\ \quad & L^{\mathbb{R}}(\mathbf{a}\mathbf{a}^T + \mathbf{b}\mathbf{b}^T) = \mathop{\mathrm{Re}}(\rho), \\ \quad & L^{\mathbb{R}}(\mathbf{b}\mathbf{a}^T - \mathbf{a}\mathbf{b}^T) = \mathop{\mathrm{Im}}(\rho), \\ \quad & L^{\mathbb{R}}(\mathbf{c}\mathbf{c}^T + \mathbf{d}\mathbf{d}^T) = \mathop{\mathrm{Re}}(\nu), \\ \quad & L^{\mathbb{R}}(\mathbf{d}\mathbf{c}^T - \mathbf{c}\mathbf{d}^T) = \mathop{\mathrm{Im}}(\nu), \\ \quad & L^{\mathbb{R}} \geq 0 \text{ on } {\mathcal{Q}}^{\mathbb{R}}(\mathbf{X}_{\mathop{\mathrm{Re}}}). \end{align}\tag{32}\] The latter formulation can be written equivalently as an infinite-dimensional linear problem over real positive Borel measures supported on \(\mathbf{X}_{\mathop{\mathrm{Re}}}\): \[\label{eq:measR} \begin{align} W^2_2(\rho,\nu) = \inf_{\mu^{\mathbb{R}}\in \mathcal{M}(\mathbf{X}_{\mathop{\mathrm{Re}}})} \quad & \int f_{\mathop{\mathrm{Re}}} \mathop{\mathrm{d}}\mu^{\mathbb{R}} \\ \text{s.t.} \quad & \int (\mathbf{a}\mathbf{a}^T + \mathbf{b}\mathbf{b}^T) \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \mathop{\mathrm{Re}}(\rho), \\ \quad & \int (\mathbf{b}\mathbf{a}^T - \mathbf{a}\mathbf{b}^T) \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \mathop{\mathrm{Im}}(\rho), \\ \quad & \int (\mathbf{c}\mathbf{c}^T + \mathbf{d}\mathbf{d}^T) \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \mathop{\mathrm{Re}}(\nu), \\ \quad & \int (\mathbf{d}\mathbf{c}^T - \mathbf{c}\mathbf{d}^T) \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \mathop{\mathrm{Im}}(\nu). \\ \end{align}\tag{33}\] Note that every feasible solution \(\mu_{\mathbb{R}}\) of 33 must be a probability measure. Indeed, combining the first equality constraint together with the fact that \(\mu_{\mathbb{R}}\) is supported on \(\mathbf{X}_{\mathop{\mathrm{Re}}}\) yields \[\int \mathop{\mathrm{Tr}}(\mathbf{a}\mathbf{a}^T + \mathbf{b}\mathbf{b}^T) \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \int (\|\mathbf{a}\|^2+\|\mathbf{b}\|^2) \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \int \mathop{\mathrm{d}}\mu^{\mathbb{R}} = \mathop{\mathrm{Tr}}(\mathop{\mathrm{Re}}(\rho))=1.\] The dual of 33 is given by \[\label{eq:dualR} \begin{align} \sup_{\Lambda_i,\Gamma_i \in \mathbb{R}^{n \times n}} \quad & \mathop{\mathrm{Tr}}(\mathop{\mathrm{Re}}(\rho) \Lambda_1 + \mathop{\mathrm{Re}}(\nu) \Gamma_1 + \mathop{\mathrm{Im}}(\rho) \Lambda_2 + \mathop{\mathrm{Im}}(\nu) \Gamma_2) \\ \text{s.t.} \quad & f_{\mathop{\mathrm{Re}}} \geq \mathop{\mathrm{Tr}}((\mathbf{a}\mathbf{a}^T + \mathbf{b}\mathbf{b}^T)\Lambda_1 + (\mathbf{c}\mathbf{c}^T + \mathbf{d}\mathbf{d}^T)\Gamma_1 + (\mathbf{b}\mathbf{a}^T - \mathbf{a}\mathbf{b}^T)\Lambda_2 \\ \quad & \quad \quad \quad + (\mathbf{d}\mathbf{c}^T - \mathbf{c}\mathbf{d}^T)\Gamma_2) \text{ on } \mathbf{X}_{\mathop{\mathrm{Re}}}, \\ \quad & \Lambda_1^T=\Lambda_1, \quad \Lambda_2^T=-\Lambda_2, \quad \Gamma_1^T=\Gamma_1, \quad \Gamma_2^T=-\Gamma_2. \end{align}\tag{34}\]

4.3 Proof of Theorem 9↩︎

We are now in the appropriate setup to apply Lemma 4.

5 Concluding remarks↩︎

We have proved a convergence rate in \(O(1/t^2)\) for the Schmüdgen-type (or equivalently Putinar-type) hierarchies of lower bounds for the minimization of a polynomial over the bi-sphere. Thanks to this result, we could provide a similar rate for a hierarchy of semidefinite programs approximating the \(2\) order quantum Wasserstein distance. We also extended this result to arbitrary sphere products, and products of distinct sets while assuming that a rate is available for each set. Our proof technique heavily relies on the polynomial kernel method.

We would like to emphasize that all obtained convergence rates for hierarchies of lower bounds also hold for the hierarchies of upper bounds developed in [34], by following the same reasoning as the one from [13]. In addition, our derived rates can be combined with the general result from [18] to obtain rates for GMP instances with measures supported on set products.

Recently, the convergence rate of the Putinar-type hierarchy for polynomial optimization over the hypercube has been revisited in [35]. The idea is also to rely on the polynomial kernel method with a suitably chosen Gaussian density. An interesting open question would be to investigate whether this idea could be extended for set products.

For the case of the unit sphere, we would like also to mention the recent works [36], [37], where the authors provide convergence rates in \(O(1/t)\) for a specialized hierarchy of spectral bounds, using quantum and real de Finetti theorems, respectively. A natural research track would be to extend this analysis to the case of the bi-sphere.

Eventually, it has been proved in [38] that the Putinar-type hierarchy converges in finitely many steps, under genericity assumptions on the objective function and the polynomials describing the feasible set. When optimizing a generic polynomial over the unit sphere, it has been shown in [39] that one also has finite convergence. Here again, it would be interesting to obtain similar finite convergence results for specific set products such as the bi-sphere.

5.0.0.1 Acknowledgments.

We thank Jonas Britz, Saroj Prasad Chhatoi, Monique Laurent and Lucas Slot for fruitful discussions. This work was supported by the European Union’s HORIZON–MSCA-2023-DN-JD programme under the Horizon Europe (HORIZON) Marie Sklodowska-Curie Actions, grant agreement 101120296 (TENORS).

References↩︎