Exact Global Optimization

Introduction

Exact global optimization problems can be solved exactly using Minimize and Maximize.

This computes the radius of the circle, centered at the origin, circumscribed about the set .
In[1]:=
Click for copyable input
Out[1]=
This computes the radius of the circle, centered at the origin, circumscribed about the set as a function of the parameters and .
In[2]:=
Click for copyable input
Out[2]=

Algorithms

Depending on the type of problem, several different algorithms can be used.

The most general method is based on the cylindrical algebraic decomposition (CAD) algorithm. It applies when the objective function and the constraints are real algebraic functions. The method can always compute global extrema (or extremal values, if the extrema are not attained). If parameters are present, the extrema can be computed as piecewise-algebraic functions of the parameters. A downside of the method is its high, doubly exponential complexity in the number of variables. In certain special cases not involving parameters, faster methods can be used.

When the objective function and all constraints are linear with rational number coefficients, global extrema can be computed exactly using the simplex algorithm.

For univariate problems, equation and inequality solving methods are used to find the constraint solution set and discontinuity points and zeros of the derivative of the objective function within the set.

Another approach to finding global extrema is to find all the local extrema, using the Lagrange multipliers or the Karush-Kuhn-Tucker (KKT) conditions, and pick the smallest (or the greatest). However, for a fully automatic method, there are additional complications. In addition to solving the necessary conditions for local extrema, it needs to check smoothness of the objective function and smoothness and nondegeneracy of the constraints. It also needs to check for extrema at the boundary of the set defined by the constraints and at infinity, if the set is unbounded. This in general requires exact solving of systems of equations and inequalities over the reals, which may lead to CAD computations that are harder than in the optimization by CAD algorithm. Currently Mathematica uses Lagrange multipliers only for equational constraints within a bounded box. The method also requires that the number of stationary points and the number of singular points of the constraints be finite. An advantage of this method over the CAD-based algorithm is that it can solve some transcendental problems, as long as they lead to systems of equations that Mathematica can solve.

Optimization by Cylindrical Algebraic Decomposition

Examples

This finds the point on the cubic curve which is closest to the origin.
In[3]:=
Click for copyable input
Out[3]=
This finds the point on the cubic curve which is closest to the origin, as a function of the parameter .
In[4]:=
Click for copyable input
Out[4]=
This visualization shows the point on the cubic curve which is closest to the origin, and the distance of the point from the origin. The value of parameter can be modified using the slider. The visualization uses the minimum computed by Minimize.
In[5]:=
Click for copyable input
Out[13]=

Customized CAD Algorithm for Optimization

The cylindrical algebraic decomposition (CAD) algorithm is available in Mathematica directly as CylindricalDecomposition. The algorithm is described in more detail in "Real Polynomial Systems". The following describes how to customize the CAD algorithm to solve the global optimization problem.

Reduction to Minimizing a Coordinate Function

Suppose it is required to minimize an algebraic function over the solution sets of algebraic constraints , where is a vector of variables and is a vector of parameters. Let be a new variable. The minimization of over the constraints is equivalent to the minimization of the coordinate function over the semialgebraic set given by .

If happens to be a monotonic function of one variable , a new variable is not needed, as can be minimized instead.

The Projection Phase of CAD

The variables are projected, with first, then the new variable , and then the parameters .

If algebraic functions are present, they are replaced with new variables; equations and inequalities satisfied by the new variables are added. The variables replacing algebraic functions are projected first. They also require special handling in the lifting phase of the algorithm.

Projection operator improvements by Hong, McCallum, and Brown can be used here, including the use of equational constraints. Note that if a new variable needs to be introduced, there is at least one equational constraint, namely .

The Lifting Phase of CAD

The parameters are lifted first, constructing the algebraic function equation and inequality description of the cells. If there are constraints that depend only on parameters and you can determine that is identically false over a parameter cell, you do not need to lift this cell further.

When lifting the minimization variable , you start with the smallest values of , and proceed (lifting the remaining variables in the depth-first manner) until you find the first cell for which the constraints are satisfied. If this cell corresponds to a root of a projection polynomial in , the root is the minimum value of , and the coordinates corresponding to of any point in the cell give a point at which the minimum is attained. If parameters are present, you get a parametric description of a point in the cell in terms of roots of polynomials bounding the cell. If there are no parameters, you can simply give the sample point used by the CAD algorithm. If the first cell satisfying the constraints corresponds to an interval , where is a root of a projection polynomial in , then is the infimum of values of , and the infimum value is not attained. Finally, if the first cell satisfying the constraints corresponds to an interval , is unbounded from below.

Strict Inequality Constraints

If there are no parameters, all constraints are strict inequalities, and you only need the extremum value, then a significantly simpler version of the algorithm can be used. (You can safely make inequality constraints strict if you know that , where is the solution set of the constraints.) In this case many lower-dimensional cells can be disregarded; hence, the projection may only consist of the leading coefficients, the resultants, and the discriminants. In the lifting phase, only full-dimensional cells need be constructed; hence, there is no need for algebraic number computations.

Experimental`Infimum[{f,cons},{x,y,...}]
find the infimum of values of f on the set of points satisfying the constraints cons.
Experimental`Supremum[{f,cons},{x,y,...}]
find the supremum of values of f on the set of points satisfying the constraints cons.

Finding extremum values.

This finds the smallest ball centered at the origin which contains the set bounded by the surface . A full Maximize call with the same input did not finish in 10 minutes.
In[14]:=
Click for copyable input
Out[14]=

Linear Optimization

When the objective function and all constraints are linear, global extrema can be computed exactly using the simplex algorithm.

This solves a random linear minimization problem with ten variables.
In[15]:=
Click for copyable input
Out[22]=
Optimization problems where the objective is a fraction of linear functions and the constraints are linear (linear fractional programs) reduce to linear optimization problems. This solves a random linear fractional minimization problem with ten variables.
In[23]:=
Click for copyable input
Out[31]=

Univariate Optimization

For univariate problems, equation and inequality solving methods are used to find the constraint solution set and discontinuity points and zeros of the derivative of the objective function within the set.

This solves a univariate optimization problem with a transcendental objective function.
In[32]:=
Click for copyable input
Out[32]=
Here is a visualization of the minimum found.
In[33]:=
Click for copyable input
Out[33]=
Here Mathematica recognizes that the objective functions and the constraints are periodic.
In[34]:=
Click for copyable input
Out[34]=

Optimization by Finding Stationary and Singular Points

Here is an example where the minimum is attained at a singular point of the constraints.
In[35]:=
Click for copyable input
Out[35]=
In[36]:=
Click for copyable input
Out[36]=
The maximum of the same objective function is attained on the boundary of the set defined by the constraints.
In[37]:=
Click for copyable input
Out[37]=
In[38]:=
Click for copyable input
Out[38]=
There are no stationary points in this example.
In[39]:=
Click for copyable input
Out[39]=

Here is a set of 3-dimensional examples with the same constraints. Depending on the objective function, the maximum is attained at a stationary point of the objective function on the solution set of the constraints, at a stationary point of the restriction of the objective function to the boundary of the solution set of the constraints, and at the boundary of the boundary of the solution set of the constraints.

Here the maximum is attained at a stationary point of the objective function on the solution set of the constraints.
In[40]:=
Click for copyable input
Out[40]=
In[41]:=
Click for copyable input
Out[41]=
Here the maximum is attained at a stationary point of the restriction of the objective function to the boundary of the solution set of the constraints.
In[42]:=
Click for copyable input
Out[42]=
In[43]:=
Click for copyable input
Out[43]=
Here the maximum is attained at the boundary of the boundary of the solution set of the constraints.
In[44]:=
Click for copyable input
Out[44]=
In[45]:=
Click for copyable input
Out[45]=
The Lagrange multiplier method works for some optimization problems involving transcendental functions.
In[46]:=
Click for copyable input
Out[46]=
In[47]:=
Click for copyable input
Out[47]=

Optimization over the Integers

Integer Linear Programming

An integer linear programming problem is an optimization problem in which the objective function is linear, the constraints are linear and bounded, and the variables range over the integers.

To solve an integer linear programming problem Mathematica first solves the equational constraints, reducing the problem to one containing inequality constraints only. Then it uses lattice reduction techniques to put the inequality system in a simpler form. Finally, it solves the simplified optimization problem using a branch-and-bound method.

This solves a randomly generated integer linear programming problem with 7 variables.
In[48]:=
Click for copyable input
Out[58]=

Optimization over the Reals Combined with Integer Solution Finding

Suppose a function needs to be minimized over the integer solution set of constraints . Let be the minimum of over the real solution set of . If there exists an integer point satisfying , then is the minimum of over the integer solution set of . Otherwise you try to find an integer solution of , and so on. This heuristic works if you can solve the real optimization problem and all the integer solution finding problems, and you can find an integer solution within a predefined number of attempts. (By default Mathematica tries 10 candidate optimum values. This can be changed using the system option.)

This finds a point with integer coordinates maximizing among the points lying below the cubic .
In[59]:=
Click for copyable input
Out[59]=
In[60]:=
Click for copyable input
Out[60]=
New to Mathematica? Find your learning path »
Have a question? Ask support »