小白自学数学建模丨M001 - 线性规划(LP)
M001 - 线性规划 (LP)
线性规划的实例与定义
例1.1 某机床厂生产甲、乙两种机床,每台销售后的利润分别为4千元和3千元。生产甲机床需要使用A,B机器加工,加工时间分别为每台2小时和1小时;生产乙机床需要A,B,C三种机器加工,加工时间为每台各一小时。若每天可用于加工的机器时数分别为A机器10小时、B机器8小时和C机器7小时,问该厂应生产甲、乙机床各几台,才能使总利润最大?
上述问题的数学模型: 设该厂生产\(x_1\)台甲机床和\(x_2\)乙机床时总利润z最大,则\(x_1, x_2\)应满足:
变量\(x_1,x_2\)称为决策变量,(1.1)式被称为问题的目标函数,(1.2)中的几个不等式是问题的约束条件,记为s.t.(即subject to)。
线性规划问题的解的概念
其中c和x为n维列向量
A、Aeq为适当维数的矩阵
b、beq为适当维数的列向量
一般线性规划问题的(数学)标准型为
可行解: 满足约束条件(1.4)的解\(x=[x_1, L, x_n]^T\),称为线性规划问题的可行解,而使目标函数(1.3)达到最大值的可行解叫最优解。
可行域: 所有可行解构成的集合称为问题的可行域,记为R。
线性规划的Matlab标准形式及软件求解
Matlab规定线性规划的标准形式为
其中c、x、b、beq、lb、ub为列向量
f称为价值向量
b称为资源向量
A、Aeq为矩阵
Matlab中求解线性规划的命令如下
[x,fval] = linprog(c,A,b)
[x,fval] = linprog(c,A,b,Aeq,beq)
[x,fval] = linprog(c,A,b,Aeq,beq,lb,ub)
其中
x返回的是决策向量的取值,
fval返回的是目标函数的最优解,
c为价值向量,
A,b对应的是线性不等式约束,
Aeq,beq对应的是线性等式约束,
lb和ub分别对应的是决策向量的下界向量和上界向量
例1.2 求解下列线性规划问题
使用Matlab的解法
f = [-2; -3; 5];
a = [-2, 5, -1; 1, 3, 1];
b = [-10; 12];
Aeq = [1, 1, 1];
beq = 7;
[x, y] = linprog(f, a, b, Aeq, beq, zeros(3 , 1));
x, y = -y
使用Python的解法
from scipy import optimize
import numpy as np
f = np.array([-2, -3, 5])
a = np.array([[-2, 5, -1], [1, 3, 1]])
b = np.array([-10, 12])
Aeq = np.array([[1, 1, 1]])
beq = np.array(7)
res = optimize.linprog(f, a, b, Aeq, beq)
print(res)
Python结果如下:
con: array([1.80714554e-09])
fun: -14.571428565645032
message: 'Optimization terminated successfully.'
nit: 5
slack: array([-2.24602559e-10, 3.85714286e+00])
status: 0
success: True
x: array([6.42857143e+00, 5.71428571e-01, 2.35900788e-10])
例1.3 求解下列线性规划问题
使用Python的解法一
from scipy import optimize
import numpy as np
f = np.array([-2, -3, -1])
a = np.array([[-2, -4, -2], [-3, -2, 0]])
b = np.array([-8, -6])
Aeq = np.array([[1, 2, 4]])
beq = np.array(101)
res = optimize.linprog(f, a, b, Aeq, beq)
print(res)
解法一结果如下所示 :
con: array([1.46057431e-07])
fun: -201.99999956122394
message: 'Optimization terminated successfully.'
nit: 7
slack: array([193.99999969, 296.99999902])
status: 0
success: True
x: array([1.01000000e+02, 1.24887119e-07, 3.11058324e-09])
使用Python的解法二
import pulp as pp
z = [2, 3, 1]
a = [[1, 4, 2], [3, 2, 0]]
b = [8, 6]
Aeq = [[1, 2, 4]]
beq = [101]
#确定最大化还是最小化问题
m = pp.LpProblem(sense = pp.LpMaximize)
#定义三个变量放到列表中如x1,x2,x3
x = [pp.LpVariable(f'x{i}', lowBound = 0) for i in[1, 2, 3]]
#定义一个目标函数并将目标函数加入求解的问题中
m += pp.lpDot(z, x)#LpDot用于计算点积
#设置比较条件
for i in range(len(a)):
m += (pp.lpDot(a[i], x) >= b[i])
#设置等式条件
for i in range(len(Aeq)):
m += (pp.lpDot(Aeq[i], x) == beq[i])
#求解问题
m.solve()
#输出结果
print(f'优化结果: {pp.value(m.objective)}')
print(f'参数取值: {[pp.value(var) for var in x]}')
解法二结果如下所示 :
Welcome to the CBC MILP Solver
Version: 2.10.3
Build Date: Dec 15 2019
command line - E:\PyChram Workspace\Learning\venv\lib\site-packages\pulp\apis\..\solverdir\cbc\win\64\cbc.exe
C:\Users\Lucifer\AppData\Local\Temp\···pulp.mps max timeMode elapsed branch printingOptions all solution
C:\Users\Lucifer\AppData\Local\Temp\···pulp.sol (default strategy 1)
At line 2 NAME MODEL
At line 3 ROWS
At line 8 COLUMNS
At line 20 RHS
At line 24 BOUNDS
At line 25 ENDATA
Problem MODEL has 3 rows, 3 columns and 8 elements
Coin0008I MODEL read with 0 errors
Option for timeMode changed from cpu to elapsed
Presolve 3 (0) rows, 3 (0) columns and 8 (0) elements
0 Obj -0 Primal inf 29.25 (3) Dual inf 5.9999997 (3)
0 Obj -0 Primal inf 29.25 (3) Dual inf 5.1666667e+10 (3)
2 Obj 202
Optimal - objective value 202
Optimal objective 202 - 2 iterations time 0.002
Option for printingOptions changed from normal to all
Total time (CPU seconds): 0.00 (Wallclock seconds): 0.00
优化结果: 202.0
参数取值: [101.0, 0.0, 0.0]