数学建模1·线性规划

最新推荐文章于 2023-07-10 22:24:21 发布

weixin_45965693

最新推荐文章于 2023-07-10 22:24:21 发布

阅读量845

点赞数 1

分类专栏： LaTex 文章标签：数学建模线性规划

本文链接：https://blog.csdn.net/weixin_45965693/article/details/107841336

版权

LaTex 专栏收录该内容

4 篇文章 0 订阅

订阅专栏

线性规划的标准形式
$\min _{x} c^{T} x\\ \ \\ \text { s.t. }\left\{\begin{array}{l} A x \leq b \\ A e q \cdot x=b e q \\ l b \leq x \leq u b \end{array}\right.$

其中c和x为n维列向量，A、Aeq为适当维数的矩阵，b、beq为适当维数的列向量

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分别对应的是决策向量的下界向量和上界向量

python scipy库求解

from scipy import optimize 
import numpy as np 
#求解函数 
res = optimize.linprog(c,A,b,Aeq,beq,LB,UB,X0,OPTI ONS) 
#目标函数最小值 
print(res.fun) 
#最优解 
print(res.x)

python pulp库求解

例1.求解
$\begin{array}{ll} \max z=2 x_{1}+3 x_{2}-5 x_{3} \\ \text { s.t. } \quad x_{1}+x_{2}+x_{3}=7 \\ 2 x_{1}-5 x_{2}+x_{3} \geq 10 \\ x_{1}+3 x_{2}+x_{3} \leq 12 \\ x_{1}, x_{2}, x_{3} \geq 0 \end{array}$

%注意这里求的是最大值,想把最大变成最小有三处需要添负号
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

#导入包
from scipy import optimize
import numpy as np
#确定c,A,b,Aeq,beq 
c = 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(-c,A,b,Aeq,beq) 
print(res)

$\begin{array}{ll} \max z=2 x_{1}+3 x_{2}+x_{3} \\ \text { s.t. } \quad x_{1}+2x_{2}+4x_{3}=101 \\ x_{1}+4x_{2}+2x_{3} \geq 8 \\ 3x_{1}+2x_{2} \geq 6 \\ x_{1}, x_{2}, x_{3} \geq 0 \end{array}$

import pulp
#目标函数系数
z = [2,3,1]
#约束 
a = [[1,4,2],[3,2,0]] 
b = [8, 6] 
#确定最大化最小化问题，最大化只要把Min改成Max即可 
m = pulp.LpProblem(sense=pulp.LpMinimize) 
#定义三个变量放到列表中 
x = [pulp.LpVariable(f’x{i}’, lowBound=0) for i in [1,2,3]] 
#定义目标函数，lpDot可以将两个列表（必须是1维的！）的对应位相乘再加和 
#相当于z[0]*x[0]+z[1]*x[1]+z[2]*x[2] 
m += pulp.lpDot(z, x)
#设置约束条件 
for i in range(len(a)): 
	m += (pulp.lpDot(a[i], x) >= b[i]) 
#求解 
m.solve() 
#输出结果 
print(f’优化结果：{pulp.value(m.objective)}’)#目标值 
print(f’参数取值：{[pulp.value(var) for var in x]})#最优解

非线性的可以转化为线性的规划问题
例2.
$\begin{array}{l} \min \left|x_{1}\right|+\left|x_{2}\right|+\cdots+\left|x_{n}\right| \\ \text { s.t. } \quad A x \leq b \end{array}$

其中 $x=\left[\begin{array}{lll} x_{1} & \cdots & x_{n} \end{array}\right]^{T}$ ，A和b为相应维数的矩阵和向量
转化：对任意的 $x_{i}$ ，存在 $x_{i},v_{i}>0$ 满足
$x_{i}=u_{i}-v_{i}, \quad\left|x_{i}\right|=u_{i}+v_{i}$
记 $u=\left[\begin{array}{lllll} u_{1} & \cdots & u_{n} \end{array}\right]^{T}, \quad v=\left[\begin{array}{lll} v_{1} & \cdots & v_{n} \end{array}\right]^{T}$

上面问题变成：
$\begin{array}{l} \min \sum_{i=1}^{n}\left(u_{i}+v_{i}\right) \\ \text { s.t. }\left\{\begin{array}{l} A(u-v) \leq b \\ u, v \geq 0 \end{array}\right. \end{array}$

案例1——投资的收益和风险

1.1问题提出

市场上有n种资产 $S_{i}(i=1,2,L,n)$ 可以选择，现用数额为M的相当大的资金作为一个时期的投资。
这n种资产在这一时期内购买 $S_{i}$ 的平均收益率为 $r_{i}$ ，风险损失率为 $q_{i}$ ，投资越分散，总的风险越少，总体风险可用投资的 $S_{i}$ 中最大的一个风险来度量。
购买 $S_{i}$ 时要付交易费，费率为 $p_{i}$ ，当购买额不超过给定值 $u_{i}$ 时，交易费按购买 $u_{i}$ 计算。
另外，假定同期银行存款利率是 $r_{0}$ ，既无交易费又无风险（ $r_{0}=5%$ )
已知 $n = 4$ 时相关数据如表

$S_{i}$	$r_{i}(\%)$	$q_{i}(\%)$	$p_{i}(\%)$	$u_{i}(元)$
$S_{1}$	28	2.5	1	103
$S_{2}$	21	1.5	2	198
$S_{3}$	23	5.5	4.5	52
$S_{4}$	25	2.6	6.5	40

试给该公司设计一种投资组合方案，即用给定资金M，有选择性地购买若干种资产或存银行生息，使净收益尽可能大，使总体风险尽可能小。

1.2符号规定和基本假设

符号	含义
$S_{i}$	表示第i种投资项目，其中 $S_{0}$ 指存入银行
$r_{i},p_{i},q_{i}$	分别表示平均收益率，交易费率，风险损失率，其中 $p_{0}=0,q_{0}=0$
$u_{i}$	表示 $s_{i}$ 的交易定额
$x_{i}$	表示投资项目 $s_{i}$ 的资金
$a$	表示投资风险度
$Q$	表示总体收益

基本假设

投资数额M相当大，为了便于计算，假设M=1;、
投资越分散，总的风险越小；
总体风险用投资项目 $s_{i}$ 中最大的一个风险来度量；
n+1中资产 $s_{i}$ 之间是相互独立的；
在投资的这一时期内， $r_{i},p_{i},q_{i}$ 为定值，不受意外因素影响；
净收益和总体风险只受 $r_{i},p_{i},q_{i}$ 影响，不受其它因素干扰。

1.3模型的分析与建立

总体风险用所投资的 $s_{i}$ 中最大的一个风险来衡量，即 $\max \left\{\boldsymbol{q}_{i} \boldsymbol{x}_{i} \mid \boldsymbol{i}=\mathbf{1}, \boldsymbol{2}, \mathbf{L}, \boldsymbol{n}\right\}$
购买 $s_{i}$ 所付交易费是一个分段函数，即 $\text { 交易费 }=\left\{\begin{array}{ll} \boldsymbol{p}_{i} \boldsymbol{x}_{i}, & \boldsymbol{x}_{i}>\boldsymbol{u}_{i} \\ \boldsymbol{p}_{i} \boldsymbol{u}_{i}, & \boldsymbol{x}_{i} \leq \boldsymbol{u}_{i} \end{array}\right.$
而题目所给的定值 $u_{i}$ 相对总投资M很少， $p_{i}u_{i}$ 更小，这样购买 $s_{i}$ 的净收益可以简化为 $r_{i}-p_{i})x_{i}$
要使净收益尽可能大，总体风险尽可能小，这是一个多目标规划模型

目标函数为
$\left\{\begin{array}{l} \max \sum_{i=0}^{n}\left(r_{i}-p_{i}\right) x_{i} \\ \min \max \left\{q_{i} x_{i}\right\} \end{array}\right.$
约束条件为
$\left\{\begin{array}{l} \sum_{i=0}^{n}\left(1+p_{i}\right) x_{i}=M \\ x_{i} \geq 0, \quad i=0,1, \cdots, n \end{array}\right.$

a）在实际投资中，投资者承受风险的程度不一样，若给定风险一个界限a，使最大的风险 $\frac{q_{i} x_{i}}{M} \leq a$ ，可找到相应的投资方案。这样把多目标规划变成一个目标的线性规划。

模型一固定风险水平，优化收益
$\max \sum_{i=0}^{n}\left(r_{i}-p_{i}\right) x_{i}\\ \ \\ \text { s.t. }\left\{\begin{array}{l} \frac{q_{i} x_{i}}{M} \leq a \\ \sum_{i=0}^{n}\left(1+p_{i}\right) x_{i}=M, \quad x_{i} \geq 0, \quad i=0,1, \cdots, n \end{array}\right.$
模型二固定盈利水平，极小化风险
$\min \left\{\max \left\{q_{i} x_{i}\right\}\right\}\\ \ \\ \text { s.t. }\left\{\begin{array}{l} \sum_{i=0}^{n}\left(r_{i}-p_{i}\right) x_{i} \geq k \\ \sum_{i=0}^{n}\left(1+p_{i}\right) x_{i}=M, \quad x_{i} \geq 0, \quad i=0,1, \cdots, n \end{array}\right.$

b)投资者在权衡资产风险和预期收益两方面时，希望选择一个令自己满意的投资组合。因此对风险、收益分别赋予权重 $\leq 1$ 和 $(1 - s)$ ，s称为投资偏好系数。

模型三
$\min s\left\{\max \left\{q_{i} x_{i}\right\}\right\}-(1-s) \sum_{i=0}^{n}\left(r_{i}-p_{i}\right) x_{i}\\ \text { s.t. } \sum_{i=0}^{n}\left(1+p_{i}\right) x_{i}=M, \quad x_{i} \geq 0, \quad i=0,1,2, \cdots, n$

1.4模型一的求解

$\begin{array}{l} \text { min } f=(-0.05,-0.27,-0.19,-0.185,-0.185)\left(x_{0}, x_{1}, x_{2}, x_{3}, x_{4}\right)^{T} \\ \text { s.t. }\left\{\begin{array}{l} x_{0}+1.01 x_{1}+1.02 x_{2}+1.045 x_{3}+1.065 x_{4}=1 \\ 0.025 x_{1} \leq a \\ 0.015 x_{2} \leq a \\ 0.055 x_{3} \leq a \\ 0.026 x_{4} \leq a \\ x_{i} \geq 0(i=0,1, \cdots, 4) \end{array}\right. \end{array}$

由于a是任意给定的风险度，到底怎样没有一个准则，不同的投资者有不同的风险度。我们从a=0开始，以步长△a=0.001进行循环搜索，编程如下

clc,clear
a=0;hold on
while a<0.05
	c=[-0.05,-0.27,-0.19,-0.185,-0.185];
	A=[zeros(4,1),diag([0.025,0.015,0.055,0.026])];
	b=a*ones(4,1);
	Aeq=[1,1.01,1.02,1.045,1.065];
	beq=1;LB=zeros(5,1);
	[x,Q]=linprog(c,A,b,Aeq,beq,LB);
	Q=-Q;plot(a,Q,'*k');
	a=a+0.001;
end
xlabel('a'),ylabel('Q')

可以看出

（1）风险大，收益也打
（2）当投资越分散式，投资者承担风险越小，这与题意一致。冒险的投资者会出现集中投资的情况，保守的投资者则尽量分散投资
（3）在 $a = 0.006$ 附件有一个转折点。在这一点左边，风险增加很少时，利润增长很快；在这一点右边，风险大增加很大时，利润增长很缓慢。所以对于风险和收益没有特殊偏好的投资者涞水，应该选择曲线的转折点作为最优投资组合，大约是 $a=0.6\%,Q=20\%$ ,所对应的投资方案为 $\begin{aligned} &\text { 风险度 } a=0.006, \text { 收益 } Q=0.2019,\\ &\quad x_{0}=0, \quad x_{1}=0.24\\ &x_{2}=0.4, x_{3}=0.1091, x_{4}=0.2212 \end{aligned}$

案例2——运输问题

某商品有m个产地、n个销地，各产地的产量分别为 $a_{1},\cdots,a_{m}$ ，各销地的需求量分别为 $b_{1},\cdots,b_{n}$ 。若该商品由i产地运到j销地的单位运价为 $c_{ij}$ ，问应该如何调运才能使总运费最省？

引入变量 $x_{ij}$ ，其取值为由i产地运往j销地的该商品数量，模型为：
$min\sum_{i=1}^{m}\sum_{j=1}^{n}c_{ij}x_{ij}\\ s.t.\left\{ \begin{aligned} \sum_{j=1}^{n}x_{ij}=a_{j}\\ \sum_{i=1}^{m}x_{ij}=b_{j}\\ x_{ij}\geq0 \end{aligned} \right.$

案例3——收益最大化

在这里插入图片描述

import pulp 
import numpy as np 
from pprint import pprint 
def transportation_problem(costs, x_max, y_max): 
	row = len(costs) 
	col = len(costs[0]) 
	prob = pulp.LpProblem(‘Transportation Problem’, sense = pulp.LpMaximize) 
	var = [[pulp.LpVariable(f’x{i}{j}’, lowBound=0, cat=pulp.LpInteger) for j in range(col)] for i in range(row)] 
	#j在列循环上，i在行循环上
	flatten = lambda x: [y for l in x for y in flatten(l)] if type(x) is list else [x] 
	#np.flatten把多维数组转化成1维的
	#从i行j列转化成了i×j一维的
	prob += pulp.lpDot(flatten(var), costs.flatten()) 
	for i in range(row): 
		prob += (pulp.lpSum(var[i]) <= x_max[i]) 
	for j in range(col): 
		prob += (pulp.lpSum([var[i][j] for i in range(row)]) <= y_max[j]) 
	prob.solve() 
	return {‘objective’:pulp.value(prob.objective), ‘var’:[[pulp.value(var[i][j])for j in range(col)]for i in range(row)]} 
if __name__ == ‘__main__’: 
	costs = np.array([[500, 550, 630, 1000, 800, 700], [800, 700, 600, 950, 900, 930], [1000, 960, 840, 650, 600, 700], [1200, 1040, 980, 860, 880, 780]]) 
	max_plant = [76, 88, 96, 40] 
	max_cultivation = [42, 56, 44, 39, 60, 59] 
	res = transportation_problem(costs, max_plant, max_cultivation) 
	print(f’最大值为{res[“objective”]}’) 	
	print(’各变量的取值为：’) 
	print(res[‘var’])

在这里插入图片描述

整数规划
若有某个变量不是整数，在松弛模型上分别添加约束： $x\leq floor(A)和x\geq ceil(A)$ 然后再分别求解，这个过程叫做分支
所谓定界，叶子节点产生后，相当于给问题定了一个下界。每次新产生叶子节点，则更新下界。
分支定界代码

import math 
from scipy.optimize import linprog 
import sys 
def integerPro(c, A, b, Aeq, beq, t=1.0E-12): 
	res = linprog(c, A_ub=A, b_ub=b, A_eq=Aeq, b_eq=beq) 
	#判断满足约束吗，不是就形成很大的数值
	if (type(res.x) is float): 
		bestX = [sys.maxsize]*len(c) 
	#是整数就读取最优解
	else:
		bestX = res.x 
	bestVal = sum([x*y for x,y in zip(c, bestX)]) 
	if all(((x-math.floor(x))<t or (math.ceil(x)- x)<t) for x in bestX): 
		return (bestVal, bestX) 
	else:
		ind = [i for i, x in enumerate(bestX) if (x-math.floor(x))>t and (math.ceil(x)-x)>t][0] 
		newCon1 = [0]*len(A[0]) 
		newCon2 = [0]*len(A[0]) 
		newCon1[ind] = -1 
		newCon2[ind] = 1 
		newA1 = A.copy() 
		newA2 = A.copy() 
		newA1.append(newCon1) 
		newA2.append(newCon2) 
		newB1 = b.copy() 
		newB2 = b.copy() 
		newB1.append(-math.ceil(bestX[ind]))
		newB2.append(math.floor(bestX[ind])) 
		r1 = integerPro(c, newA1, newB1, Aeq, beq) 
		r2 = integerPro(c, newA2, newB2, Aeq, beq) 
		if r1[0] < r2[0]：
			return r1 
		else:
			return r2 
c = [3, 4, 1] 
A = [[-1, -6, -2], [-2, 0, 0]] 
b = [-5, -3] 
Aeq = [[0,0,0]] 
beq = [0] 
print(integerPro(c, A, b, Aeq, beq))

在这里插入图片描述
非线性规划

例3：计算 $\frac{1}{x}+x$ 的最小值

from scipy.optimize import minimize 
import numpy as np 
#计算 1/x+x 的最小值 
def fun(args): 
	a=args 
	v=lambda x:a/x[0] +x[0] 
	return v
if __name__ == "__main__": 
	args = (1) #a 
	x0 = np.asarray((2)) # 初始猜测值 
	res = minimize(fun(args), x0, method='SLSQP’) 
	print(res.fun) 
	print(res.success) 
	print(res.x)

在这里插入图片描述

函数的最小是为2点多，可以看出minimize求的局部最优

例4：计算 $\frac{(2+x_{1})}{1+x_{2}}-3x_{1}+4x_{3}$ 的最小值，其中 $x_{1}、x_{2}、x_{3}$ 的范围在0.1到0.9之间

from scipy.optimize import minmize
import numpy as np
def fun(args):
	a,b,c,d=args
	v=lambda x: (a+x[0])/(b+[1])-c*x[0]+d*x[2]
	return v
	
def con(args): 
	# 约束条件 分为eq 和ineq 
	#eq表示 函数结果等于0 ； ineq 表示 表达式大于等于0 					
	x1min,x1max,x2min,x2max,x3min,x3max = args 
	cons = ({'type':'ineq','fun':lambda x:x[0]-x1min},\ {'type':'ineq','fun': lambda x:-x[0]+x1max},\ {'type':'ineq','fun': lambda x:x[1]-x2min},\ {'type':'ineq','fun': lambda x:-x[1]+x2max},\ {'type':'ineq','fun': lambda x:x[2]-x3min},\ {'type':'ineq','fun': lambda x:-x[2]+x3max}) 
	return cons
	
if __name__ == "__main__": 
	#定义常量值 
	args = (2,1,3,4) #a,b,c,d 
	#设置参数范围/约束条件 
	args1 = (0.1,0.9,0.1,0.9,0.1,0.9) #x1min, x1max, x2min, x2max
	cons = con(args1) 
	#设置初始猜测值 
	x0 = np.asarray((0.5,0.5,0.5)) 
	res = minimize(fun(args),x0, method='SLSQP',constraints=cons) 
	print(res.fun) 
	print(res.success) 
	print(res.x)