Projected Gradient

In this notebook, we are going to use projected gradient to solve the following optimization problem.

$$\begin{align} \min \quad 2(x_1 - 1)^2 &+ (x_2 - 2)^2 \\ \mathrm{s.t.} \quad x_1 + 4x_2 - 2 &\le 0 \\ \quad -x_1 + x_2 &\le 0 \\ \quad x_1, x_2 &\ge 0 \end{align} $$

Iteration 1

In [1]:
import numpy as np
from numpy.linalg import inv
from scipy.optimize import minimize_scalar
In [2]:
x1 = 0
x2 = 0
g1 = 4*(x1 - 1)
g2 = 2*(x2 - 2)
gradF = np.array([g1,g2])

Note that we have constraints -x1 + x2 <= 0, x1 >= 0, x2 >= 0.

We only use -x1 + x2 <= 0 and x2 >= 0 as it automatically implies x1 >= 0. If we use all three, then we cannot invert the matrix.

In [3]:
E = np.array([[-1,1],[0,-1]])
I = np.identity(2)
ET = np.transpose(E)
P = I - np.dot(np.dot(ET,inv(np.dot(E,ET))),E)
print(P)

[[ 0.  0.]
 [ 0.  0.]]
In [4]:
d = -np.dot(P, gradF)
print(d)

[-0. -0.]

Since d=0, we calculate Q.

In [5]:
Q = -np.dot(np.dot(inv(np.dot(E,ET)),E),gradF)
print(Q)

[-4. -8.]

Since the second row is the most negative, we delete the second row of the matrix E.

In [6]:
E = np.array([[-1,1]])
I = np.identity(2)
ET = np.transpose(E)
P = I - np.dot(np.dot(ET,inv(np.dot(E,ET))),E)
print(P)

[[ 0.5  0.5]
 [ 0.5  0.5]]
In [7]:
d = -np.dot(P, gradF)
print(d)

[ 4.  4.]

x1 = 0 + lam*4

x2 = 0 + lam*4

lam <= 2/20

In [8]:
def fun(l):
    return 2*(0 + 4*l - 1)**2 + (0 + 4*l - 2)**2

lam = minimize_scalar(fun, bounds=(0, 0.1), method='bounded') 
print(lam.x)

0.0999933893039
In [9]:
x1 = 0 + 4*lam.x
x2 = 0 + 4*lam.x

print("x1= ", x1)
print("x2= ", x2)

x1=  0.399973557215
x2=  0.399973557215

Iteration 2

In [10]:
x1 = 0.4
x2 = 0.4
g1 = 4*(x1 - 1)
g2 = 2*(x2 - 2)
gradF = np.array([g1,g2])
In [11]:
E = np.array([[1,4],[-1,1]])
I = np.identity(2)
ET = np.transpose(E)
P = I - np.dot(np.dot(ET,inv(np.dot(E,ET))),E)
print(P)

[[  1.11022302e-16   1.11022302e-16]
 [  0.00000000e+00   2.22044605e-16]]
In [12]:
d = -np.dot(P, gradF)
print(d)

[  6.21724894e-16   7.10542736e-16]
In [13]:
Q = -np.dot(np.dot(inv(np.dot(E,ET)),E),gradF)
print(Q)

[ 1.12 -1.28]

Since the second row has negative values, delete the second row from the E matrix.

In [14]:
E = np.array([[1,4]])
I = np.identity(2)
ET = np.transpose(E)
P = I - np.dot(np.dot(ET,inv(np.dot(E,ET))),E)
print(P)

[[ 0.94117647 -0.23529412]
 [-0.23529412  0.05882353]]
In [15]:
d = -np.dot(P, gradF)
print(d)

[ 1.50588235 -0.37647059]
In [16]:
def fun(l):
    return 2*(0.4 + 1.505*l - 1)**2 + (0.4 - 0.3764*l - 2)**2

lam = minimize_scalar(fun, bounds=(0, 1.0627), method='bounded') 
print(lam.x)

0.257669168234
In [17]:
x1 = 0.4 + 1.505*lam.x
x2 = 0.4 - 0.3764*lam.x

print("x1= ", x1)
print("x2= ", x2)

x1=  0.787792098192
x2=  0.303013325077