数学建模：非线性规划的 Python 求解_python非线性规划求解

如果目标函数或约束条件中包含非线性函数，就称这种规划问题为非线性规划问题。一般来说，求解非线性规划要比线性规划困难得多，而且，不像线性规划有通用的方法。非线性规划目前还没有适用于各种问题的一般算法，各个方法都有自己特定的适用范围。

一般形式

非线性规划的一般形式为
min ⁡ f ( x ) ,  s.t.  { G ( x ) ≤ 0 , H ( x ) = 0 ,

$$\begin{gathered} \min \quad f(x), \\ \text { s.t. }\left\{\begin{array}{l} G(x) \leq 0, \\ H(x)=0, \end{array}\right. \\ \end{gathered}$$ minf(x), s.t. {G(x)≤0,H(x)=0,​​

其中 $G(x)={\left[g_1(x),g_2(x),\cdots ,g_m(x)\right]}^T$, $H(x)={\left[h_1(x),h_2(x),\cdots ,h_l(x)\right]}^T$，$f,g_i,h_j$ 都是定义在 ${\mathbb{R}}^n$ 上的实值函数

Python 求解

在使用 Python 求解非线性规划时要注意求解器的适用范围，各求解器的适用范围如下

1. scipy.optimize.minimize 函数

调用方式:

import scipy
res = scipy.optimize.minimize(fun, x0, args=(), method=None, jac=None, hess=None, hessp=None, bounds=None, constraints=(), tol=None, callback=None, options=None)
1
2

参数:

fun: 目标函数: fun(x, *args) x 是形状为 (n,) 的一维数组 (n 是自变量的数量)

x0: 初值, 形状 (n，) 的 ndarray

method: 求解器的类型，如果未给出，有约束时则选择 'SLSQP'。

constraints: 求解器 SLSQP 的约束被定义为字典的列表，字典用到的字段主要有：

• 'type': str: 约束类型：“eq” 表示等式，“ineq” 表示不等式 (SLSQP 的不等式约束是 $f(x)\ge 0$ 的形式)
• 'fun': 可调用的函数或匿名函数：定义约束的函数

bounds: 变量界限

返回:

print(res.fun) #最优值
print(res.success) #求解状态
print(res.x)  #最优解
1
2
3

e.g.

from scipy.optimize import minimize
import numpy as np
def obj(x):
    # 定义目标函数
    return f(x)
LB=[LB]; UB=[UB] # n 个元素的列表
cons=[{'type':'ineq','fun':lambda x: G(x)},{'type':'eq','fun':lambda x: H(x)}]
bd=tuple(zip(LB, UB)) # 生成决策向量界限的元组
res=minimize(obj,np.ones(n),constraints=cons,bounds=bd) #第2个参数为初值
print(res.fun) #最优值
print(res.success) #求解状态
print(res.x)  #最优解
1
2
3
4
5
6
7
8
9
10
11
12

例 1

$$\min \frac{2+x_1}{1+x_2}-3x_1+4x_3\text{, }$$

$$s.t.0.1\le x_i\le 0.9,i=1,2,3$$

from scipy.optimize import minimize
from numpy import ones
def obj(x):
    x1,x2,x3=x
    return (2+x1)/(1+x2)-3*x1+4*x3
LB=[0.1]*3; UB=[0.9]*3
bound=tuple(zip(LB, UB))   
res=minimize(obj,ones(3),bounds=bound) 
print(res.fun,'\n',res.success,'\n',res.x)  #输出最优值、求解状态、最优解
1
2
3
4
5
6
7
8
9

-0.7736842105263159 
 True 
 [0.9 0.9 0.1]
1
2
3

例 2

max ⁡ 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 ,  s.t.  { 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 0 ≤ x i ≤ 99 , i = 1 , 2 , ⋯   , 5.

$$\begin{gathered} \max 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}, \\ \text { s.t. }\left\{\begin{array}{l} x_{1}+x_{2}+x_{3}+x_{4}+x_{5} \leq 400, \\ x_{1}+2 x_{2}+2 x_{3}+x_{4}+6 x_{5} \leq 800 \\ 2 x_{1}+x_{2}+6 x_{3} \leq 200 \\ x_{3}+x_{4}+5 x_{5} \leq 200 \\ 0 \leq x_{i} \leq 99, \quad i=1,2, \cdots, 5 . \end{array}\right. \end{gathered}$$ maxz=x12​+x22​+3x32​+4x42​+2x52​−8x1​−2x2​−3x3​−x4​−2x5​, s.t. ⎩ ⎨ ⎧​x1​+x2​+x3​+x4​+x5​≤400,x1​+2x2​+2x3​+x4​+6x5​≤8002x1​+x2​+6x3​≤200x3​+x4​+5x5​≤2000≤xi​≤99,i=1,2,⋯,5.​​

from scipy.optimize import minimize
import numpy as np
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])
obj=lambda x: np.dot(-c1,x**2)+np.dot(-c2,x)
# np.dot()向量点乘或矩阵乘法
cons={'type':'ineq','fun':lambda x:b-A@x}
bd=[(0,99) for i in range(A.shape[1])]
#res=minimize(obj,np.ones(5)*90,constraints=cons,bounds=bd)
res=minimize(obj,np.ones(5),constraints=cons,bounds=bd)
print(res.fun) 
print(res.success) 
print(res.x)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

-42677.56889347886
True
[0.00000000e+00 1.43988405e-12 3.33333334e+01 9.90000000e+01
 1.35333333e+01]
1
2
3
4

本题只能求得局部最优解，可使用不同初值尝试求解

2. cvxopt.solvers 模块求解二次规划

标准型：

min ⁡ 1 2 x T P x + q T x ,  s.t.  { A x ≤ b , A e q ⋅ x = b e q .

$$\begin{gathered} \min \quad \frac{1}{2} x^{T} P x+q^{T} x, \\ \text { s.t. }\left\{\begin{array}{l} A x \leq b, \\ A e q \cdot x=b e q . \end{array}\right. \end{gathered}$$ min21​xTPx+qTx, s.t. {Ax≤b,Aeq⋅x=beq.​​

调用方式:

import numpy as np
from cvxopt import matrix,solvers
sol=solvers.qp(P,q,A,b,Aeq,beq) #二次规划
'''
输入到cvxopt.solvers的矩阵 P,q,A,b,Aeq,beq 都是 cvxopt.matrix 形式
在程序中虽然没有直接使用 NumPy 库中的函数，也要加载
数据如果全部为整型数据，也必须写成浮点型数据
cvxopt.matrix(array,dims) 把array按照dims重新排成矩阵
    省略dims: 如果array为np.array, 则为其原本形式; 
        如果array为list, 则认为 list 中用逗号分隔开的为一"列"
'''
print(sol['x']) #最优解
print(sol['primal objective']) #最优值
1
2
3
4
5
6
7
8
9
10
11
12
13

例 3

min ⁡ z = 1.5 x 1 2 + x 2 2 + 0.85 x 3 2 + 3 x 1 − 8.2 x 2 − 1.95 x 3 ,  s.t.  { x 1 + x 3 ≤ 2 − x 1 + 2 x 2 ≤ 2 x 2 + 2 x 3 ≤ 3 x 1 + x 2 + x 3 = 3

$$\begin{gathered} \min \quad z=1.5 x_{1}^{2}+x_{2}^{2}+0.85 x_{3}^{2}+3 x_{1}-8.2 x_{2}-1.95 x_{3}, \\ \text { s.t. }\left\{\begin{array}{l} x_{1}+x_{3} \leq 2 \\ -x_{1}+2 x_{2} \leq 2 \\ x_{2}+2 x_{3} \leq 3 \\ x_{1}+x_{2}+x_{3}=3 \end{array}\right. \end{gathered}$$ minz=1.5x12​+x22​+0.85x32​+3x1​−8.2x2​−1.95x3​, s.t. ⎩ ⎨ ⎧​x1​+x3​≤2−x1​+2x2​≤2x2​+2x3​≤3x1​+x2​+x3​=3​​

解:

P = [ 3 0 0 0 2 0 0 0 1.7 ] , q = [ 3 − 8.2 − 1.95 ] , A = [ 1 0 1 − 1 2 0 0 1 2 ] , b = [ 2 2 3 ] ,

$$\begin{gathered} P=\left[\begin{array}{ccc} 3 & 0 & 0 \\ 0 & 2 & 0 \\ 0 & 0 & 1.7 \end{array}\right], q=\left[\begin{array}{c} 3 \\ -8.2 \\ -1.95 \end{array}\right], A=\left[\begin{array}{ccc} 1 & 0 & 1 \\ -1 & 2 & 0 \\ 0 & 1 & 2 \end{array}\right], b=\left[\begin{array}{l} 2 \\ 2 \\ 3 \end{array}\right], \\ \end{gathered}$$ P=⎣ ⎡​300​020​001.7​⎦ ⎤​,q=⎣ ⎡​3−8.2−1.95​⎦ ⎤​,A=⎣ ⎡​1−10​021​102​⎦ ⎤​,b=⎣ ⎡​223​⎦ ⎤​,​

$Aeq=[1,1,1],beq=[3].$

import numpy as np
from cvxopt import matrix,solvers
n=3; P=matrix(0.,(n,n)); #零矩阵
P[::n+1]=[3,2,1.7]; #对角元变为3, 2, 1.7
q=matrix([3,-8.2,-1.95])
A=matrix([[1.,0,1],[-1,2,0],[0,1,2]]).T
b=matrix([2.,2,3])
Aeq=matrix(1.,(1,n)); beq=matrix(3.)
s=solvers.qp(P,q,A,b,Aeq,beq)
print("最优解为：",s['x'])
print("最优值为：",s['primal objective'])
1
2
3
4
5
6
7
8
9
10
11

     pcost       dcost       gap    pres   dres
 0: -1.3148e+01 -4.4315e+00  1e+01  1e+00  9e-01
 1: -6.4674e+00 -7.5675e+00  1e+00  2e-16  1e-16
 2: -7.1538e+00 -7.1854e+00  3e-02  1e-16  4e-16
 3: -7.1758e+00 -7.1761e+00  3e-04  1e-16  2e-15
 4: -7.1760e+00 -7.1760e+00  3e-06  5e-17  2e-15
Optimal solution found.
最优解为： [ 8.00e-01]
[ 1.40e+00]
[ 8.00e-01]

最优值为： -7.1759977687772745
1
2
3
4
5
6
7
8
9
10
11
12

3. cvxpy 库

使用方式与这里相似，只需选择支持的求解器即可

例 4

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 ,  s.t.  { 0 ≤ x i ≤ 99 ,  且  x i  为整数  ( i = 1 , ⋯   , 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.

$$\begin{gathered} \min \quad 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}, \\ \text { s.t. }\left\{\begin{array}{l} 0 \leq x_{i} \leq 99, \text { 且 } x_{i} \text { 为整数 }(i=1, \cdots, 5), \\ x_{1}+x_{2}+x_{3}+x_{4}+x_{5} \leq 400, \\ x_{1}+2 x_{2}+2 x_{3}+x_{4}+6 x_{5} \leq 800 \\ 2 x_{1}+x_{2}+6 x_{3} \leq 200 \\ x_{3}+x_{4}+5 x_{5} \leq 200 . \end{array}\right. \end{gathered}$$ minz=x12​+x22​+3x32​+4x42​+2x52​−8x1​−2x2​−3x3​−x4​−2x5​, s.t. ⎩ ⎨ ⎧​0≤xi​≤99, 且 xi​ 为整数 (i=1,⋯,5),x1​+x2​+x3​+x4​+x5​≤400,x1​+2x2​+2x3​+x4​+6x5​≤8002x1​+x2​+6x3​≤200x3​+x4​+5x5​≤200.​​

import cvxpy as cp
import numpy as np
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()
print("最优值为:",prob.value)
print("最优解为:",x.value)
1
2
3
4
5
6
7
8
9
10
11
12
13

最优值为: -17.0
最优解为: [4. 1. 0. 0. 0.]
1
2