I am trying to solve a linear programming problem with Python. Ling program fails to find a solution. But quad program works. I do not understand why, and I am not sure my formulation of the program in linprog and quad program are equivalent.
The liner programming problem, my code, and the error message from linprog are below.
Code
k = 6
n = 3
indexes = [1, 2, 5, 6, 7, 9]
V = np.zeros((1, k))
count = 0
for ID in xrange(1, 4):
ind = count * n + ID
p = indexes.index(ind)
V[0, p] = 1
count += 1
bounds = []
for i in xrange(6):
bounds.append((0, 1))
bounds = tuple(bounds)
W2 = np.identity(k)
W2 = np.multiply(W2, -1.0)
b2 = np.transpose(np.zeros(k))
V = [[1.0, 0.0, 1.0, 0.0, 0.0, 1.0]]
W1 = [[10., 0., 0., 0., 30., 0.],
[ 0., 10., 20., 0., 0., 0.],
[ 0., 0., 0., 20., 0., 30.]]
W1 = np.array(W1)
W3 = [[1., 1., 0., 0., 0., 0.],
[0., 0., 1., 1., 0., 0.],
[0., 0., 0., 0., 1., 1.]]
W2 = np.array(W2)
b1 = [10., 20., 30.]
b3 = [1., 1., 1.,]
EQ = np.vstack([W1, W3]).T
Eb = np.vstack([b1, b3]).T
M = np.identity(k)
P = dot(M.T, M)
q = np.zeros(k)
def quadprog_solve_qp(P, q, W2, b2, W1, b1, W3, b3):
qp_G = .5 * (P + P.T) # make sure P is symmetric
qp_a = -q
qp_C = -np.vstack([W3, W1, W2]).T
qp_b = -np.hstack([b3, b1, b2])
meq = W1.shape[0]
return quadprog.solve_qp(qp_G, qp_a, qp_C, qp_b, meq)[0]
a = quadprog_solve_qp(P=P, q=q, W2=W2, b2=b2, W1=W1, b1=b1, W3=W3, b3=b3)
print a
V = [1] * k
res = opt.linprog(c=V, A_eq=EQ, b_eq=Eb, bounds=bounds, options={"disp": True})
Error message from linprog failing
Optimization failed. Unable to find a feasible starting point.
You can install quad program with
sudo pip install quadprog
For examples using quadprog see
https://scaron.info/blog/quadratic-programming-in-python.html
Your problem has no solution. You have the constraints [0, 0, 0.4, 0, 0, 0] * x = 0.8
and x[2] <= 1
but x[2]
has to be 2 to fulfil the first constraint.
Here is an example:
'''
min [1, 1, 1, 1, 1, 1] * x
with [ 0 0 0 0 0 0 ] [ 0 ]
[ 0 0 0.4 0 0 0 ] * x = [ 0.4 ]
[ 0 0 0 0.5 0 0 ] [ 0.25 ]
and [ 0 0 0 0 0 0 ] [ 0 ]
[ 0 0 0 0 0.4 0 ] * x = [ 0.04 ]
[ 0 0 0 0 0 0.5 ] [ 0.125 ]
and 0 <= x[i] <= 1, for 0 <= i < 6
'''
import scipy.optimize as opt
import numpy as np
k = 6
n = 3
W1 = np.zeros((n, k))
W1[1, 2] = 0.4
W1[2, 3] = 0.5
b1 = np.zeros([n, 1])
b1[1] = 0.4
b1[2] = 0.25
W3 = np.zeros((n, k))
W3[1, 4] = 0.4
W3[2, 5] = 0.5
b3 = np.zeros([n, 1])
b3[1] = 0.04
b3[2] = 0.125
c = np.array([1] * k)
A_eq = np.vstack([W1, W3])
b_eq = np.vstack([b1, b3])
bounds = np.array([(0, 1)] * k)
options = {"disp": True}
result = opt.linprog(c=c, A_eq=A_eq, b_eq=b_eq, bounds=bounds, options=options)
print(result.x)
I fixed the dimensions of the vectors and matrices. Also I changed the entries of matrices and vectors, so that a solution can be found. Here is the output:
Optimization terminated successfully.
Current function value: 1.850000
Iterations: 4
[0. 0. 1. 0.5 0.1 0.25]
I didn't try quadprog but if your quadprog function finds a solution, then I'm sure that the problems are not equivalent as there is no solution for the original problem.
Here is another example with other values:
import scipy.optimize as opt
import numpy as np
k = 6
n = 3
W1 = np.array([[0.35855621, 0, 0, 0, 0.85993925, 0 ],
[0, 0.35855621, 0.05538291, 0, 0, 0 ],
[0, 0, 0, 0.05538291, 0, 0.85993925]])
b1 = np.array([[0.35855621],
[0.05538291],
[0.85993925]])
c = np.array([1] * k)
A_eq = W1
b_eq = b1
bounds = np.array([(0, 1)] * k)
options = {"disp": True}
result = opt.linprog(c=c, A_eq=A_eq, b_eq=b_eq, bounds=bounds, options=options)
print(result.x)
The output is:
Optimization terminated successfully.
Current function value: 1.571416
Iterations: 3
[0. 0.15446089 0. 0. 0.41695528 1. ]