python数学建模的书 数学建模python教程

文章目录

• （一）简单陈述本文章的内容
• （二）常用导入文件方式
• （三）线性规划

• 3.1 线性规划的一般模型
• 3.2 运用python各种库和模块求解线性方程

• 3.2.1 Scipy线性规划模型
• 3.2.2 pulp线性规划模型
• 3.2.3 cvxopt.solvers 模块求解
• 3.2.4 用cvxpy库求解

（一）简单陈述本文章的内容

python建模会持续更新，用途是只作为个人笔记。我博客中的所有资料都可通过我提供的链接永久获取，最后希望大家一起相互促进，相互努力。
本文章只介绍如何用python中的库和模块求解相对应的建模，不介绍相应的理论知识，理论知识可以查看链接下载：《数学建模算法与应用》 提取码：fx4q

本文章所有（源码、文件）获取：链接：https://pan.baidu.com/s/1SFazXZ5uqJ2A6rIYXiG0cA 提取码：kozt

（二）常用导入文件方式

• 为了方便初学者使用，下面是把几种常用的导入文件方法写入到一个定义函数中。
  books文件下载：books提取码：ouvw --来自百度网盘的分享

# 数据导入例程

import pandas as pd

# 读取数据文件
def readDataFile(readPath):  # readPath: 数据文件的地址和文件名
    # readPath = "../data/youcansxupt.csv"  # 文件路径也可以直接在此输入
    try:
        if (readPath[-4:] == ".csv"):
            dfFile = pd.read_csv(readPath, header=0, sep=",")  # 间隔符为逗号，首行为标题行
            # dfFile = pd.read_csv(filePath, header=None, sep=",")  # sep: 间隔符，无标题行
        elif (readPath[-4:] == ".xls") or (readPath[-5:] == ".xlsx"):  # sheet_name 默认为 0
            dfFile = pd.read_excel(readPath, header=0)  # 首行为标题行
            # dfFile = pd.read_excel(filePath, header=None)  # 无标题行
        elif (readPath[-4:] == ".dat"):  # sep: 间隔符，header：首行是否为标题行
            dfFile = pd.read_table(readPath, sep=" ", header=0)  # 间隔符为空格，首行为标题行
            # dfFile = pd.read_table(filePath,sep=",",header=None) # 间隔符为逗号，无标题行
        else:
            print("不支持的文件格式。")
    except Exception as e:
        print("读取数据文件失败：{}".format(str(e)))
        return
    return dfFile


# 主程序
def main():
    # 读取数据文件
    readPath = "books.csv"  # 数据文件的地址和文件名
    dfFile = readDataFile(readPath)  # 调用读取文件子程序

    print(type(dfFile))  # 查看 dfFile 数据类型
    print(dfFile.shape)  # 查看 dfFile 形状（行数，列数）
    print(dfFile.head())  # 显示 dfFile 前 5 行数据

    return


if __name__ == '__main__':
    main()

OUT:

<class 'pandas.core.frame.DataFrame'>
(10000, 23)
   id  ...                                    small_image_url
0   1  ...  https://images.gr-assets.com/books/1447303603s...
1   2  ...  https://images.gr-assets.com/books/1474154022s...
2   3  ...  https://images.gr-assets.com/books/1361039443s...
3   4  ...  https://images.gr-assets.com/books/1361975680s...
4   5  ...  https://images.gr-assets.com/books/1490528560s...

[5 rows x 23 columns]

（三）线性规划

3.1 线性规划的一般模型

max(min) = $\sum_{j=1}^n$ $c_j$$x_j$,
$s.t. \begin{cases} \sum_{j=1}^n a_{ij}x_{j} \le(\ge) b_{i}, &i = 1,2,3,···,m\\ x_{j}>=0, &j = 1,2,3,···,m \end{cases}$

3.2 运用python各种库和模块求解线性方程

min z = $c^Tx$,
$s.t. \begin{cases} A \bullet x \le b, &\\ Aeq \bullet x = beq, &\\ LB \le x \le UB \end{cases}$
linprog 的基本调用格式为

from scipy.optimize import linprog

# 默认每个决策变量下界为0，上界为正无穷
res = linprog(c,A,b,Aeq,beq) 
res = linprog(c, A = None, b = None, Aeq = None, beq = None, bounds = None, mothod = 'simplex')
print(res.fun)  # 显示目标函数的最小值
print(res.x)    # 显示最优解

3.2.1 Scipy线性规划模型

max z = $2x_1 + 3x_2 - 5x_3$,
$s.t. \begin{cases} x_1 + 3x_2 + x_3 \le 12,&\\ 2x_1 - 5x_2 + x_3 \ge 10,\\ x_1 + x_2 +x_3 = 7,\\ x_1, x_2, x_3 \ge 0 \end{cases}$

# Scipy线性规划模型
from scipy.optimize import linprog

c = [-2,-3,5]; A = [[1,3,1],[-2,5,-1]]   # 或者 c = np.array([-2,-3,5])
b = [[12],[-10]]; Aeq = [[1,1,1]]; beq = [7]
LB = [0,0,0]
UB = [None]*len(c)   # 生成3个None的列表
bound = tuple(zip(LB,UB))   #生成决策向量界限的元组
res = linprog(c,A,b,Aeq,beq,bound)
print("目标函数的最小值:",res.fun)
print("最优解为:",res.x)
# 目标函数的最小值: -14.571428565645057, 所以最大值为:14.571428565645057
# 最优解为: [6.42857143e+00 5.71428571e-01 2.35900788e-10]

OUT：

目标函数的最小值: -14.571428565645057
最优解为: [6.42857143e+00 5.71428571e-01 2.35900788e-10]

Process finished with exit code 0

3.2.2 pulp线性规划模型

• 与2.2.1一样的题目，用不同模型求解
  max z = $2x_1 + 3x_2 - 5x_3$,
  $s.t. \begin{cases} x_1 + 3x_2 + x_3 \le 12,&\\ 2x_1 - 5x_2 + x_3 \ge 10,\\ x_1 + x_2 +x_3 = 7,\\ x_1, x_2, x_3 \ge 0 \end{cases}$

# pulp线性规划模型
import pulp

# 定义一个规划问题
MyProbLP = pulp.LpProblem("LPProbDemo1", sense=pulp.LpMaximize)  # 求最大值

# 定义决策变量
x1 = pulp.LpVariable('x1', lowBound=0, upBound=7, cat='Continuous')
x2 = pulp.LpVariable('x2', lowBound=0, upBound=7, cat='Continuous')
x3 = pulp.LpVariable('x3', lowBound=0, upBound=7, cat='Continuous')

# 添加目标函数
MyProbLP += 2*x1 + 3*x2 - 5*x3  	# 设置目标函数
# 添加约束条件
MyProbLP += (2*x1 - 5*x2 + x3 >= 10)  # 不等式约束
MyProbLP += (x1 + 3*x2 + x3 <= 12)  # 不等式约束
MyProbLP += (x1 + x2 + x3 == 7)  # 等式约束

MyProbLP.solve()
print("Status:", pulp.LpStatus[MyProbLP.status]) # 输出求解状态
for v in MyProbLP.variables():  # youcans
    print(v.name, "=", v.varValue)  # 输出每个变量的最优值
print("Max F(x) = ", pulp.value(MyProbLP.objective))  #输出最优解的目标函数值

OUT：

Status: Optimal
x1 = 6.4285714
x2 = 0.57142857
x3 = 0.0
Max F(x) =  14.57142851

Process finished with exit code 0

3.2.3 cvxopt.solvers 模块求解

• cvxopt.solvers 模块求解线性规划模型的标准型如下：
  min z = $c^Tx$,
  $s.t. \begin{cases} A \bullet x \le b, &\\ Aeq \bullet x = beq. \end{cases}$
• 求解线性规划
  min z = $-4x_1 - 5x_2$,
  $s.t. \begin{cases} 2x_1 + x_2 \le 3,\\ x_1 + 2x_2 \le 3,&\\ x_1 \ge 0, x_2 \ge 0. \end{cases}$

# 线性规划 使用cvxopy.solvers模块求解
import numpy
from cvxopt import matrix,solvers

c = matrix([2.,1]); A = matrix([[-1.,1],[-1,-1],[1,-2],[0,-1]]).T  # matrix([2.,1]):记住matrix中第一个元素不许是float型
b = matrix([1.,-2,4,0]); Aeq = matrix([1.,2],(1,2))  # Aeq为行向量
beq = matrix(3.5); sol = solvers.lp(c,A,b,Aeq,beq)
print("最优解为:\n",sol['x'])
print("最优值为:\n",sol['primal objective'])

OUT:

D:\ProgramData\Anaconda3\python.exe C:/Users/Administrator/PycharmProjects/pythonProject/数学建模/线性规划/线性规划3.cvxopy.solvers.py
     pcost       dcost       gap    pres   dres   k/t
 0:  5.5556e+00  1.2222e+00  1e+01  0e+00  2e+00  1e+00
 1:  4.6038e+00  3.7995e+00  2e+00  1e-16  4e-01  2e-01
 2:  2.5229e+00  2.4639e+00  2e-01  2e-16  4e-02  4e-02
 3:  2.5002e+00  2.4997e+00  2e-03  1e-16  4e-04  4e-04
 4:  2.5000e+00  2.5000e+00  2e-05  4e-16  4e-06  4e-06
 5:  2.5000e+00  2.5000e+00  2e-07  1e-16  4e-08  4e-08
Optimal solution found.
最优解为:
 [ 5.00e-01]
[ 1.50e+00]

最优值为:
 2.500000024611047

3.2.4 用cvxpy库求解

P183 . 例题: 已知某种商品 6 个仓库的存货量，8个客户对该商品的需求量，单位商品运价如下表所示。试确定6个仓库到8个客户的商品调运数量，使总的运输费用最小。

仓库：单位运价：客户	V1	V2	V3	V4	V5	V6	V7	V8	存货量
W1	6	2	6	7	4	2	5	9	60
W2	4	9	5	3	8	5	3	2	55
W3	5	2	1	9	7	4	3	3	51
W4	7	6	7	3	9	2	7	1	43
W5	2	3	9	5	7	2	6	5	41
W6	5	5	2	2	8	1	4	3	52
需求量	35	37	22	32	41	32	43	38

解 ： 设$x_{ij}$(i = 1,2,·······,8)表示第 i 个仓库运到第 j 个客户的商品数量，$c_{ij}$表示第 i个仓库到底 j 个客户的单位运价， $d_j$ 表示第 j 个客户的需求量，$e_i$ 表示第 i 仓库的存货量，建立如下线性规划模型：
min $\sum_{i=1}^6$$\sum_{j=1}^8$$c_{ij}x_{ij}$,
$s.t. \begin{cases} \sum_{j=1}^8 x_{ij} \le e_i,&i = 1,2,···,6\\ \sum_{i=1}^6 x_{ij} = d_j,&j = 1,2,···,8\\ x_{ij} \ge 0, &i = 1,2,···,6;j = 1,2,···,8. \end{cases}$
把上表中7行9列的数据保存到Excel文件，并命名为工作簿3.xlsx。
链接：工作簿3.xlsx 提取码：kzmp

import cvxpy as cp
import numpy as np
import pandas as pd

d1 = pd.read_excel("工作簿3.xlsx",header=None)
d2 = d1.values; c = d2[:-1,:-1]
d = d2[-1,:-1].reshape(1,-1)
e = d2[:-1,-1].reshape(-1,1)
x = cp.Variable((6,8))
obj = cp.Minimize(cp.sum(cp.multiply(c,x)))
con = [cp.sum(x,axis=1,keepdims=True)<=e,cp.sum(x,axis=0,keepdims=True)==d,x>=0]
prob = cp.Problem(obj,con)
prob.solve(solver='GLPK_MI',verbose=True)
print("最优值为:\n",prob.value)
print("最优解为:\n",x.value)

OUT:

最优值为:
 664.0
最优解为:
 [[-0. 19. -0. -0. 41. -0. -0. -0.]
 [-0. -0. -0. 32. -0. -0. -0.  1.]
 [-0. 12. 22. -0. -0. -0. 17. -0.]
 [-0. -0. -0. -0. -0.  6. -0. 37.]
 [35.  6. -0. -0. -0. -0. -0. -0.]
 [-0. -0. -0. -0. -0. 26. 26. -0.]]

Process finished with exit code 0