PolynomialOptimize

PolynomialOptimize
Bounds the optimal value of a homogeneous polynomial on the unit sphere
Other toolboxes required	none
Related functions	CopositivePolynomial PolynomialAsMatrix PolynomialSOS
Function category	Polynomial optimization
Usable within CVX?	yes

PolynomialOptimize is a function that computes upper and lower bounds on the optimum (i.e., maximum or minimum) value of an even-degree homogeneous polynomial on the unit sphere. That is, given an even-degree n-variable homogeneous polynomial p, this function bounds the (very difficult to compute in general) quantity

\( \max\big\{ p(x_1,x_2,\ldots,x_n) : x_1^2 + x_2^2 + \cdots + x_n^2 = 1 \big\}, \)

or the corresponding minimization problem.

Syntax

• OB = PolynomialOptimize(P,N,D,K)
• OB = PolynomialOptimize(P,N,D,K,OPTTYPE)
• [OB,IB] = PolynomialOptimize(P,N,D,K)
• [OB,IB] = PolynomialOptimize(P,N,D,K,OPTTYPE)

Argument descriptions

Input arguments

• P: The polynomial, represented as a vector of coefficients of its monomials in lexicographic order.
• N: The number of variables in the polynomial.
• D: Half the degree of the polynomial.
• K: A non-negative integer that specifies the level of the hierarchy that is used to compute the bound(s) on the optimal value. Larger values of K result in better bounds at the expense of increases memory usage and longer computation times.
• OPTTYPE (optional, default 'max'): Either 'max' or 'min', indicating whether the polynomial should be maximized or minimized.

Output arguments

• OB: An "outer" bound on the optimal value. If maximizing, this is an upper bound on the maximum; if minimizing, this is a lower bound on the minimum.
• IB: An "inner" bound on the optimal value. If maximizing, this is a lower bound on the maximum; if minimizing, this is an upper bound on the minimum.

Examples

The Motzkin polynomial

The Motzkin polynomial is the following degree-6 homogeneous polynomial\[p(x,y,z) = x^4y^2 + x^2y^4 - 3x^2y^2z^2 + z^6\]. We can load this polynomial into a list-of-coefficients form that PolynomialOptimize understands via the exp2ind helper function:

>> n = 3;% number of variables
>> d = 3;% half the degree of the polynomial
>> p = zeros(nchoosek(n+2*d-1,2*d),1);% initialize p to have the correct size
>> p(exp2ind([4 2 0])) = 1;% the "4 2 0" here are the exponents in the term x^4y^2z^0, while the "1" is the coefficient of this term
>> p(exp2ind([2 4 0])) = 1;
>> p(exp2ind([2 2 2])) = -3;
>> p(exp2ind([0 0 6])) = 1;% now p represents the Motzkin polynomial

The Motzkin polynomial, despite the fact that it cannot be written as a sum of squares, is non-negative. In fact, the minimum value of this polynomial over the unit sphere is exactly 0. We can get lower bounds for the optimal value as follows:

>> PolynomialOptimize(p,n,d,0,'min')

ans =

   -0.2000

>> PolynomialOptimize(p,n,d,1,'min')

ans =

   -0.1104

>> PolynomialOptimize(p,n,d,2,'min')

ans =

   -0.0780

We can see above that increasing the value of K (i.e., the level of the hierarchy) results in a better lower bound (i.e., a bound that is closer to the true minimum value of 0). This function is fast enough that we can go to a much higher level of the hierarchy, and we can also get upper bounds on this minimum value:

>> [ob,ib] = PolynomialOptimize(p,n,d,10,'min')

ob =

   -0.0158


ib =

    1.0333e-06

>> [ob,ib] = PolynomialOptimize(p,n,d,100,'min')

ob =

   -0.0017


ib =

    1.5581e-09

>> [ob,ib] = PolynomialOptimize(p,n,d,500,'min')

ob =

  -3.3898e-04


ib =

    4.6048e-13

In other words, by the 500-th level of the hierarchy we have learned that the minimum value of the Motzkin polynomial is between -3.3898e-04 and 4.6048e-13. Note that the inner bounds are computed by an algorithm that includes some randomness, so when you run the code above you might get different inner bounds.

Source code

Click here to view this function's source code on github.

Notes

The inner bounds, even though they are somewhat random, do (deterministically) converge to the true optimal value of the polynomial as K approaches infinity.

References

Coming soon.