03. 非线性规划
定义
如果目标函数或约束条件中包含非线性函数，就称这种规划问题为非线性规划问题。
非线性规划目前还没有适于各种问题的一般算法。
非线性规划模型描述如：
min
?
f
(
x
)
,
\min f(x),
minf(x),
s
.
t
.
{
g
i
(
x
)
≤
0
,
i
=
1
,
2
,
…
,
m
,
h
j
(
x
)
=
0
,
j
=
1
,
2
,
…
,
l
,
s.t.\left\{\begin{aligned} &g_i(x)\le0,i=1,2,\dots,m,\\ &h_j(x)=0,j=1,2,\dots,l, \end{aligned}\right.
s.t.{?gi?(x)≤0,i=1,2,…,m,hj?(x)=0,j=1,2,…,l,?
其中
x
=
[
x
1
,
x
2
,
…
,
x
n
]
T
x=[x_1,x_2,\dots,x_n]^T
x=[x1?,x2?,…,xn?]T，而
f
,
g
i
,
h
j
f,g_i,h_j
f,gi?,hj? 都是实值函数。
一般非线性规划只能得到局部最优解，不能保证是全局最优解。
无约束非线性规划求解
设
f
(
x
)
f(x)
f(x) 具有连续的一阶偏导数，且
x
?
x^*
x?? 是无约束问题的局部极小点，则
?
f
(
x
?
)
=
0
\nabla f(x^*)=0
?f(x?)=0
其中
?
f
(
x
)
\nabla f(x)
?f(x) 表示
f
(
x
)
f(x)
f(x) 的梯度。
设
f
(
x
)
f(x)
f(x) 具有对各个变量的二阶偏导数，称矩阵
(
?
2
f
?
x
1
2
?
2
f
?
x
1
?
x
2
…
?
2
f
?
x
1
?
x
n
?
2
f
?
x
2
?
x
1
?
2
f
?
x
2
2
…
?
2
f
?
x
2
?
x
n
?
?
?
?
?
2
f
?
x
n
?
x
1
?
2
f
?
x
n
?
x
2
…
?
2
f
?
x
n
2
)
\begin{pmatrix} \frac{\partial^2f}{\partial x_1^2} & \frac{\partial^2f}{\partial x_1\partial x_2} & \dots & \frac{\partial^2f}{\partial x_1\partial x_n} \\ \frac{\partial^2f}{\partial x_2\partial x_1} & \frac{\partial^2f}{\partial x_2^2} & \dots & \frac{\partial^2f}{\partial x_2\partial x_n} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial^2f}{\partial x_n\partial x_1} & \frac{\partial^2f}{\partial x_n\partial x_2} & \dots & \frac{\partial^2f}{\partial x_n^2} \\ \end{pmatrix}
?????????x12??2f??x2??x1??2f???xn??x1??2f???x1??x2??2f??x22??2f???xn??x2??2f??……?…??x1??xn??2f??x2??xn??2f???xn2??2f??????????
为函数的黑塞矩阵，记
?
2
f
(
x
)
\nabla^2f(x)
?2f(x)。
因此，只要找到
x
?
x^*
x?，满足
?
f
(
x
?
)
=
0
\nabla f(x^*)=0
?f(x?)=0，且
?
2
f
(
x
?
)
\nabla^2f(x^*)
?2f(x?) 为正定矩阵，则
x
?
x^*
x? 为无约束优化问题的局部最优解。
找到
x
?
x^*
x? 的具体方法有最速降线法、牛顿法等。
有约束非线性规划求解
有等式约束非线性规划的 Lagrange 乘数法
对于只有等式约束的非线性规划问题：
min
?
f
(
x
)
,
\min f(x),
minf(x),
s
.
t
.
{
h
j
(
x
)
=
0
,
j
=
1
,
2
,
…
,
l
,
x
∈
R
n
.
s.t.\left\{\begin{aligned} & h_j(x)=0,j=1,2,\dots,l,\\ & x\in R^n. \end{aligned}\right.
s.t.{?hj?(x)=0,j=1,2,…,l,x∈Rn.?
设函数
f
,
h
j
(
j
=
1
,
2
,
…
,
l
)
f,h_j(j=1,2,\dots,l)
f,hj?(j=1,2,…,l) 在可行点
x
?
x^*
x? 的某个邻域
N
(
x
?
,
?
)
N(x^*,\epsilon)
N(x?,?)? 内可微，向量组
?
h
j
(
x
?
)
\nabla h_j(x^*)
?hj?(x?) 线性无关，令
L
(
x
,
λ
)
=
f
(
x
)
?
λ
T
H
(
x
)
,
L(x,\lambda)=f(x)-\lambda^TH(x),
L(x,λ)=f(x)?λTH(x),
其中
λ
=
[
λ
1
,
λ
2
,
…
,
λ
l
]
T
\lambda=[\lambda_1,\lambda_2,\dots,\lambda_l]^T
λ=[λ1?,λ2?,…,λl?]T
H
(
x
)
=
[
h
1
(
x
)
,
h
2
(
x
)
,
…
,
h
l
(
x
)
]
T
H(x)=[h_1(x),h_2(x),\dots,h_l(x)]^T
H(x)=[h1?(x),h2?(x),…,hl?(x)]T
若
x
?
x^*
x? 是局部最优解，则有
λ
?
=
[
λ
1
?
,
λ
2
?
,
…
,
λ
l
?
]
T
\lambda^*=[\lambda_1^*,\lambda_2^*,\dots,\lambda_l^*]^T
λ?=[λ1??,λ2??,…,λl??]T
使得
?
L
(
x
?
,
λ
?
)
=
0
\nabla L(x^*,\lambda^*)=0
?L(x?,λ?)=0，即
?
f
(
x
?
)
?
∑
j
=
1
l
λ
j
?
?
h
j
(
x
?
)
=
0
\nabla f(x^*)-\sum_{j=1}^l\lambda_j^*\nabla h_j(x^*)=0
?f(x?)?j=1∑l?λj???hj?(x?)=0
从而将有约束条件转化为无约束条件。
有约束非线性规划的罚函数法
构造带参数的所谓增广目标函数，从而把有约束非线性规划问题转化为一系列无约束非线性规划问题。
增广目标函数通常由两部分组成：
对于外点法，增广目标函数形式如下：
T
(
x
,
M
)
=
f
(
x
)
+
M
∑
i
=
1
m
[
max
?
{
0
,
g
i
(
x
)
}
]
+
M
∑
j
=
1
l
[
h
j
(
x
)
]
2
T(x,M)=f(x)+M\sum_{i=1}^m[\max\{0,g_i(x)\}]+M\sum_{j=1}^l[h_j(x)]^2
T(x,M)=f(x)+Mi=1∑m?[max{0,gi?(x)}]+Mj=1∑l?[hj?(x)]2
其中 M 是一个较大的正数。
从而转化为无约束问题：
min
?
T
(
x
,
M
)
,
x
∈
R
n
\min T(x,M),x\in R^n
minT(x,M),x∈Rn
罚函数法的计算精度可能较差。
Python 代码
利用 scipy、cvxopt、cvxpy 包，可以实现非线性规划求解。
scipy 求解
对于问题：
min
?
2
+
x
1
1
+
x
2
?
3
x
1
+
4
x
3
,
\min \frac{2+x_1}{1+x_2}-3x_1+4x_3,
min1+x2?2+x1???3x1?+4x3?,
s
.
t
.
0.1
≤
x
i
≤
0.9
,
i
=
1
,
2
,
3
s.t.\quad0.1\le x_i\le0.9,i=1,2,3
s.t.0.1≤xi?≤0.9,i=1,2,3
使用 scipy 求解，代码如下：
import numpy as np
from scipy.optimize import minimize
def obj(x):
x1,x2,x3 = x
return (2+x1)/(1+x2)-3*x1+4*x3
bound = [(.1, .9) for _ in range(3)]
res = minimize(obj, np.ones(3), bounds=bound)
print('best obj =', res.fun)
print('success =', res.success)
print('best x =', res.x)
输出如下：
best obj = -0.7736842105263159
success = True
best x = [0.9 0.9 0.1]
cvxpy 求解
对于问题：
min
?
z
=
x
1
2
+
x
2
2
+
3
x
3
2
+
4
x
4
2
+
2
x
5
2
?
8
x
1
?
2
x
2
?
3
x
3
?
x
4
?
2
x
5
,
\min z=x_1^2+x_2^2+3x_3^2+4x_4^2+2x_5^2-8x_1-2x_2-3x_3-x_4-2x_5,
minz=x12?+x22?+3x32?+4x42?+2x52??8x1??2x2??3x3??x4??2x5?,
s
.
t
.
{
0
≤
x
i
≤
99
,
i
=
1
,
2
,
…
,
5
x
i
∈
Z
+
,
i
=
1
,
2
,
…
,
5
x
1
+
x
2
+
x
3
+
x
4
+
x
5
≤
400
,
x
1
+
2
x
2
+
2
x
3
+
x
4
+
6
x
5
≤
800
,
2
x
1
+
x
2
+
6
x
3
≤
200
,
x
3
+
x
4
+
5
x
5
≤
200.
s.t.\left\{\begin{aligned} & 0\le x_i\le99,i=1,2,\dots,5\\ & x_i\in Z^+,i=1,2,\dots,5\\ & x_1+x_2+x_3+x_4+x_5\le400,\\ & x_1+2x_2+2x_3+x_4+6x_5\le800,\\ & 2x_1+x_2+6x_3\le200,\\ & x_3+x_4+5x_5\le200. \end{aligned}\right.
s.t.???????????????????????0≤xi?≤99,i=1,2,…,5xi?∈Z+,i=1,2,…,5x1?+x2?+x3?+x4?+x5?≤400,x1?+2x2?+2x3?+x4?+6x5?≤800,2x1?+x2?+6x3?≤200,x3?+x4?+5x5?≤200.?
代码如下：
import numpy as np
import cvxpy as cp
c1 = np.array([1, 1, 3, 4, 2])
c2 = np.array([-8, -2, -3, -1, -2])
a = np.array([[1, 1, 1, 1, 1],
[1, 2, 2, 1, 6],
[2, 1, 6, 0, 0],
[0, 0, 1, 1, 5]])
b = np.array([400, 800, 200, 200])
x = cp.Variable(5, integer=True)
obj = cp.Minimize(c1 @ x**2 + c2 @ x)
con = [0 <= x, x <= 99, a@x <= b]
prob = cp.Problem(obj, con)
prob.solve(solver='CPLEX')
print('best z =', prob.value)
print('best x =', x.value)
输出如下：
best z = -17.0
best x = [4. 1. 0. 0. 0.]
|