Lagrange multipliers example with sympy - all minima but one maxima.

Question

Consider the following optimization problem:

Minimize $x^3+y^3$

Subject to: $x^2+y^2 \leq 1$

On the boundary of the constraint, we can consider $x=\cos\theta$ and $y=\sin\theta$.

Then, the objective function becomes $\cos^3\theta+\sin^3\theta$. Plotting it with $\theta$, we get the following graph:

It's clear that the local minima $\frac{\pi}{4}$, $\pi$ and $\frac{3 \pi}{2}$.

Now, I want to get the local minima of the constrained optimization problem above using Lagrange multipliers.

The Lagrangian becomes: $L(x,y,\lambda) = x^3+y^3-\lambda(-x^2-y^2+1)$. And the KKT conditions (equations 12.30 in the book by Nocedal and Wright) yield:

$$3x^2+2\lambda x = 0 \tag{1}$$ $$3y^2+2\lambda y = 0 \tag{2}$$ $$x^2+y^2 \leq 1 \tag{3}$$ $$\lambda(x^2+y^2-1)=0 \tag{4}$$ $$\lambda \geq 0 \tag{5}$$

Now, we convert these into a system of polynomial equations. First, we replace $\lambda$ by $\lambda^2$ so we don't have to worry about $\lambda \geq 0$ and eliminate equation (5). Next, we convert equation (3) into an equality by introducing a dummy variable.

$$x^2+y^2-1=-\kappa^2 \tag{6}$$

Since $\kappa^2 \geq 0$ for real $\kappa$, this is equivalent to (3). Now, I plug equations (1), (2), (6) and (4) into sympy's polynomial equation solver:

from sympy import *
x, y, z, l, m, k = symbols('x y z l m k')
solve([Eq(3*x**2+2*x*l**2,0),
   Eq(3*y**2+2*y*l**2,0),
   Eq(x**2+y**2+k**2,1),
   Eq(x**2*l+y**2*l-l,0)], [x,y,l,k])

This produces the following solutions to this system:

[(-1, 0, -sqrt(6)/2, 0),
 (-1, 0, sqrt(6)/2, 0),
 (0, -1, -sqrt(6)/2, 0),
 (0, -1, sqrt(6)/2, 0),
 (0, 0, 0, -1),
 (0, 0, 0, 1),
 (0, 1, -sqrt(6)*I/2, 0),
 (0, 1, sqrt(6)*I/2, 0),
 (1, 0, -sqrt(6)*I/2, 0),
 (1, 0, sqrt(6)*I/2, 0),
 (-sqrt(2)/2, -sqrt(2)/2, -2**(1/4)*sqrt(3)/2, 0),
 (-sqrt(2)/2, -sqrt(2)/2, 2**(1/4)*sqrt(3)/2, 0),
 (sqrt(2)/2, sqrt(2)/2, -2**(1/4)*sqrt(3)*I/2, 0),
 (sqrt(2)/2, sqrt(2)/2, 2**(1/4)*sqrt(3)*I/2, 0)]

Any solution that involves imaginary numbers for $\lambda$ or $\kappa$ should be ignored since we require their squares to be $\geq 0$.

This gives us $(-1,0)$ which corresponds to $\pi$ in the graph above and this is a local minima, $(0,-1)$ which corresponds to $\frac{3 \pi}{2}$ in the graph above and this corresponds to a local minima as well. Next, we have $(0,0)$ and this corresponds to an inflection point at the very center of the feasible region. And finally, $(-\frac{1}{\sqrt{2}},-\frac{1}{\sqrt{2}})$ which corresponds to $\frac{5 \pi}{4}$. Herein lies the problem. All the points before this, clearly corresponded to local minima. But this one actually corresponds to a local maxima as can be seen in the graph above.

So it would seem the KKT conditions are picking all the local minima, but somehow swapping one of them with a local maxima. What am I missing here?

Answer 1

In this case, since the only critical point in the interior of the feasible set is a saddle point, all local maxima and minima will occur on the boundary. In this sense it is enough to use Lagrange multipliers handling only the equality constraint. However, if you really want to use the KKT conditions, keep in mind that:

1. If $a \in D$ is a local minimiser (some regularity assumptions are require), it satisfies

$$ \begin{cases} 3x^2 -2 \lambda x = 0\\ 3y^2 - 2 \lambda y = 0\\ x^2+y^2 \leq 1 \\\lambda \leq 0 \\ \lambda (x^2+y^2-1) = 0 \end{cases} $$

2. If $a \in D$ is a local maximiser (again, some regularity assumptions are require), it satisfies

$$ \begin{cases} 3x^2 -2 \lambda x = 0\\ 3y^2 - 2 \lambda y = 0\\ x^2+y^2 \leq 1 \\\lambda \ge 0 \\ \lambda (x^2+y^2-1) = 0 \end{cases} $$

Solving these two systems you will get all the candidates to local minimiser and maximiser. The point is that when you fix the sign of the multiplier you a restricting yourself to find either minima or maxima.

Keep in mind that these are only necessary conditions.

EDIT:

After solving both systems, you will see that the candidates to minimisers are $$ (0,0), (-1,0), (0,-1), (-\sqrt{2}/2, -\sqrt{2}/2) $$

and the candidates to maximisers are

$$ (0,0), (1,0), (0,1), (\sqrt{2}/2, \sqrt{2}/2) $$

Looking at the behaviour of the objective function on the boundary and on the straight line $y=x$, you can see that $(0,0)$, $(\sqrt{2}/2, \sqrt{2}/2)$ and $(-\sqrt{2}/2, -\sqrt{2}/2)$ are saddle points. The local maximisers are $(0,1),(1,0)$ and the local minimisers are $(0,-1),(-1,0)$. Compactness arguments and the inexistence of irregular points actually allows to show that these are the global maximisers/minimisers.

Answer 2

Calling

$$ f(x,y) = x^3+y^3\\ g(x,y) = x^2+y^2-1 $$

we have the Lagrangian

$$ L(x,y,\lambda,\epsilon) = f+\lambda(g+\epsilon^2) $$

So the stationary points are obtained by solving

$$ \nabla L = 0 = \left\{ \begin{array}{l} \epsilon ^2+x^2+y^2-1 \\ 3 x^2+2 \lambda x \\ 3 y^2+2 \lambda y \\ 2 \epsilon \lambda \\ \end{array} \right. $$

with results

$$ \left[ \begin{array}{cccccc} \lambda & x & y & \epsilon & f(x,y) & \det(H_g)\\ \frac{3}{2} & -1 & 0 & 0 & -1 & 0 \\ \frac{3}{2} & 0 & -1 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 \\ -\frac{3}{2} & 0 & 1 & 0 & 1 & 0 \\ -\frac{3}{2} & 1 & 0 & 0 & 1 & 0 \\ \frac{3}{2 \sqrt{2}} & -\frac{1}{\sqrt{2}} & -\frac{1}{\sqrt{2}} & 0 & -\frac{1}{\sqrt{2}} & 12 \sqrt{2} \\ -\frac{3}{2 \sqrt{2}} & \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} & 0 & \frac{1}{\sqrt{2}} & -12 \sqrt{2} \\ \end{array} \right] $$

where $\det(H_g)$ represents the bordered Hessian determinant

$$ H_g = \left( \begin{array}{ccc} 0 & 2 x & 2 y \\ 2 x & 6 x & 0 \\ 2 y & 0 & 6 y \\ \end{array} \right) $$

needed to qualify the stationary points. So as we can verify, the optimal points are at $(0,-1)$ and $(-1,0)$ as global minima and $(0,1), \ (1,0)$ as global maxima. the last two solutions cover local extrema. The first a local maximum and the last a local minimum

NOTE

Attached a plot showing the objective function level surface lines into the feasible region, showing in red the constraint gradient and in black the objective function gradient at each stationary point. At the maximum points the combination $\nabla f = \lambda \nabla g$ has $\lambda > 0$ as well as the minimum we have $-\nabla f = \lambda f $ also with $\lambda > 0$

Attached the MATHEMATICA script to produce the shown graphics.

usols = {{x -> -1, y -> 0, s1 -> 0, l1 -> 3/2}, {x -> 0, y -> -1, 
s1 -> 0, l1 -> 3/2}, {x -> 0, y -> 0, s1 -> -1, l1 -> 0}, {x -> 0,
 y -> 1, s1 -> 0, l1 -> -(3/2)}, {x -> 1, y -> 0, s1 -> 0, 
l1 -> -(3/2)}, {x -> -(1/Sqrt[2]), y -> -(1/Sqrt[2]), s1 -> 0, 
l1 -> 3/(2 Sqrt[2])}, {x -> 1/Sqrt[2], y -> 1/Sqrt[2], s1 -> 0, 
l1 -> -(3/(2 Sqrt[2]))}};

toler = 0.000001;

\[Lambda] = Sqrt[2]/5;

xinf = -1.6;
yinf = -1.6;
xsup = 1.6;
ysup = 1.6;

f[x_, y_] := x^3 + y^3
g1[x_, y_] := x^2 + y^2 - 1
G = {g1[x, y]};
For [i = 1; ListaSols = {}; ListaVecs = {}, i <= Length[usols], i++,
  s = {s1} /. usols[[i]];
  For[j = 1, j <= Length[s], j++, If[Abs[s[[j]]] < toler,
    v = Grad[-G[[j]], {x, y}] /. usols[[i]];
    nv = Norm[v];
    If[nv > 0,
      p1 = {x, y} /. usols[[i]];
      p2 = p1 + \[Lambda] v/nv;
    AppendTo[ListaVecs, Graphics[{Red, Arrow[{p1, p2}]}]]]
 ]
];
 v = Grad[f[x, y], {x, y}] /. usols[[i]];
 nv = Norm[v];
 p1 = {x, y} /. usols[[i]];
 If[nv > 0,
   p2 = p1 + \[Lambda] v/nv;
   AppendTo[ListaVecs, Graphics[{Arrow[{p1, p2}]}]]];
   AppendTo[ListaSols, {{Red, PointSize[0.015], Point[p1]}}]
 ]

grp = ContourPlot[f[x, y], {x, xinf, xsup}, {y, yinf, ysup}, 
RegionFunction -> ((g1[#1, #2] < 0) &), Contours -> 125, 
Epilog -> ListaSols, PlotPoints -> 50];
Show[grp, ListaVecs, AspectRatio -> 1]