QP definition
参考文献见:https://courses.csail.mit.edu/6.867/wiki/images/a/a7/Qp-cvxopt.pdf
二次型规划可划归为以下模型
minx12xTPx+qTxsubjecttoGx≤hAx=b \begin{aligned} min_{x} \quad \frac{1}{2}x^TPx+q^Tx \\ subject \quad to \quad Gx \le h \\ \quad \quad \quad Ax=b \\ \end{aligned} minx21xTPx+qTxsubjecttoGx≤hAx=b
一个例子
- 原问题
minx,y12x2+3x+4ysubjecttox,y≥0x+3y≥152x+5y≤1003x+4y≤80 \begin{aligned} min_{x,y} \quad \frac{1}{2}x^2+3x+4y \\ subject \quad to\quad x,y\ge0 \\ \quad \quad \quad x+3y\ge15 \\ \quad \quad \quad 2x+5y\le100 \\ \quad \quad \quad 3x+4y\le80 \\ \end{aligned} minx,y21x2+3x+4ysubjecttox,y≥0x+3y≥152x+5y≤1003x+4y≤80 - 矩阵表达,公式
minx,y12[xy]T[1000][xy]+[34][xy][−100−1−1−32534][xy]≤[00−1510080] \begin{aligned}min_{x,y} \quad \frac{1}{2} \begin{bmatrix}x \\y \end{bmatrix}^T \begin{bmatrix}1 &0 \\0&0 \end{bmatrix}\begin{bmatrix}x \\y \end{bmatrix}+ \begin{bmatrix}3 \\4 \end{bmatrix}\begin{bmatrix}x \\y \end{bmatrix} \\\begin{bmatrix}-1 & 0 \\ 0 & -1 \\ -1 & -3 \\2 &5 \\3 & 4\end{bmatrix}\begin{bmatrix}x \\y \end{bmatrix} \le\begin{bmatrix} 0 \\ 0 \\ -15 \\ 100 \\80\end{bmatrix} \end{aligned} minx,y21[xy]T[1000][xy]+[34][xy]⎣⎢⎢⎢⎢⎡−10−1230−1−354⎦⎥⎥⎥⎥⎤[xy]≤⎣⎢⎢⎢⎢⎡00−1510080⎦⎥⎥⎥⎥⎤
求解
import numpy as np
import cvxopt
from cvxopt import matrix,solvers
使用numpy 和 matrix结合
P = matrix(np.diag([1,0]), tc='d')
q = matrix(np.array([3,4]), tc='d')
G = matrix(np.array([[-1,0],[0,-1],[-1,-3],[2,5],[3,4]]),tc='d')
h = matrix(np.array([0,0,-15,100,80]), tc='d')
sol = solvers.qp(P,q,G,h)
#sol = solvers.qp(P,q,G,h,A,b)
print(sol['x'])
print(sol['primal objective'])
pcost dcost gap pres dres
0: 1.0780e+02 -7.6366e+02 9e+02 4e-17 4e+01
1: 9.3245e+01 9.7637e+00 8e+01 8e-17 3e+00
2: 6.7311e+01 3.2553e+01 3e+01 8e-17 1e+00
3: 2.6071e+01 1.5068e+01 1e+01 7e-17 7e-01
4: 3.7092e+01 2.3152e+01 1e+01 1e-16 4e-01
5: 2.5352e+01 1.8652e+01 7e+00 9e-17 4e-16
6: 2.0062e+01 1.9974e+01 9e-02 7e-17 2e-16
7: 2.0001e+01 2.0000e+01 9e-04 8e-17 2e-16
8: 2.0000e+01 2.0000e+01 9e-06 1e-16 2e-16
Optimal solution found.
[ 7.13e-07]
[ 5.00e+00]
20.00000617311241
使用 matrix直接求解
这里matrix的输入是转置的,比如G
Glist = [[-1.0,0.0,-1.0,2.0,3.0],[0.0,-1.0,-3.0,5.0,4.0]]
G = matrix([[-1.0,0.0,-1.0,2.0,3.0],[0.0,-1.0,-3.0,5.0,4.0]])
print(Glist)
print(G)
[[-1.0, 0.0, -1.0, 2.0, 3.0], [0.0, -1.0, -3.0, 5.0, 4.0]]
[-1.00e+00 0.00e+00]
[ 0.00e+00 -1.00e+00]
[-1.00e+00 -3.00e+00]
[ 2.00e+00 5.00e+00]
[ 3.00e+00 4.00e+00]
P = matrix(np.diag([1,0]), tc='d')
q = matrix(np.array([3,4]), tc='d')
G = matrix(np.array([[-1,0],[0,-1],[-1,-3],[2,5],[3,4]]),tc='d')
h = matrix(np.array([0,0,-15,100,80]), tc='d')
sol = solvers.qp(P,q,G,h)
#sol = solvers.qp(P,q,G,h,A,b)
print(sol['x'])
print(sol['primal objective'])
pcost dcost gap pres dres
0: 1.0780e+02 -7.6366e+02 9e+02 4e-17 4e+01
1: 9.3245e+01 9.7637e+00 8e+01 8e-17 3e+00
2: 6.7311e+01 3.2553e+01 3e+01 8e-17 1e+00
3: 2.6071e+01 1.5068e+01 1e+01 7e-17 7e-01
4: 3.7092e+01 2.3152e+01 1e+01 1e-16 4e-01
5: 2.5352e+01 1.8652e+01 7e+00 9e-17 4e-16
6: 2.0062e+01 1.9974e+01 9e-02 7e-17 2e-16
7: 2.0001e+01 2.0000e+01 9e-04 8e-17 2e-16
8: 2.0000e+01 2.0000e+01 9e-06 1e-16 2e-16
Optimal solution found.
[ 7.13e-07]
[ 5.00e+00]
20.00000617311241