数学建模——线性规划模型

前言：当学习完线性规划模型，我感觉到了数学建模的“细腻”之处，也可以从中感觉到他“细腻”的美感，为此想记录一下我学习数学建模的一些笔记跟心得。

线性规划模型一般是求解最大值最小值问题，如果目标函数f(x)和约束条件均是决策变量的线性表达式，（即没有平方项和乘积项），那么此时的数学规划问题就属于线性规划问题，下面我将通过对软件得到的结果进行分析，进而改进我们的模型。

首先线性回归具有三个重要因素：目标函数、约束条件、决策变量

1、决策变量为自变量（一般含有多个自变量），即函数中的x。

2、目标函数即f（x)。

3、约束条件是符号，即大于号，小于号，等号。它分为三种，即不等式约束，等式约束，整数约束。

线性规划问题的解的概念

线性规划的目标函数可以求最大值，也可以求最小值，约束条件的不等号可以是小于号也可以是大于号。下面我将分别对目标函数的最大值与最小值进行分析。（matlab版)

matlab中的代码：

经典案例：

求解过程：因为上面的案例是求最大值的，要满足目标函数的最小值，需要变换符号。

备注：1、不等号必须为小于号，若为大于号的需变成小于号，即加个负号。

           2、aeq为等于号左边x前面的系数，beq为等号右边的数，a为左边不等号式子X前的系数，b为不等号右边的数，注意要把不等号变成小于号，f为min中的x的系数，注意是min,不是max,如果为max,需加负号。lb为x的下限，ub为x的上限。xo一般为迭代初始量，一般不用给。如果上面条件哪个没有，在代码linprog中直接写【】。

           3、上面的y等于-y是因为求解最大值问题，求解最小值问题只需要y等于y.

           4、如果形式为x<0,没有等于号，这时我们可以放松一下条件，即x<=1（结合题意放松）。

           5、本道题中的zeros(3,1)表示下限，因为x1,x2,x3都大于0，即写成三行一列都为0，即

matlab运行的结果：

可以转化为线性规划的问题：

案例：

建模的分析和建议：

模型一解法：

matlab的代码解法：

matlab代码运行的结果：

python解法：

import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import linprog

c = [-0.05,-0.27,-0.19,-0.185,-0.185]
A = np.c_[np.zeros(4),np.diag([0.025,0.015,0.055,0.026])]
# np.r_ : 按列连接两个矩阵，即把两矩阵上下相加，要求列数相等
# np.c_ : 按行连接两个矩阵，即把两矩阵左右相加，要求行数相等
Aeq =[[1,1.01,1.02,1.045,1.065]]; beq = [1]
a=0; aa=[]; rr=[]
while a<0.05:
    b = np.ones(4)*a
    res = linprog(c,A,b,Aeq,beq)
    x = res.x; Q = -res.fun
    aa.append(a); rr.append(Q) #保存最优值
    a = a+0.001

plt.rc('font',size=16); plt.rc('font',family='SimHei'); plt.rc('text',usetex=False)
plt.plot(aa,rr,'k')
plt.plot(aa,rr,'r.')
plt.xlabel("风险 a"); plt.ylabel("收益 R",rotation=0)
plt.show()

案例2：

解法：

       1、因为要经过A和B两道工序加工，即产品在A中完成的件数与在B中完成的件数一样

     备注：将上面的min修改为max

matlab代码解法：

format long g

c = zeros(9,1);

c(1) = 1.25-0.25-300/600*5;

c(2) = 1.25-0.25-321/10000*7;

c(3) = -250/4000*6;

c(4) = -783/7000*4;

c(5) = -200/4000*7;

c(6) = -300/6000*10;

c(7) = -321/10000*7;

c(8) = 2-0.35-250/4000*8;

c(9) = 2.8-0.5-321/10000*12-783/7000*11;

c=-c;

A = zeros(5,9);

A (1,1) = 5; A (1,6) = 10;

A (2,2) = 7; A (2,7) = 9; A (2,9) = 12;

A (3.3) = 6; A (3,8) = 8;

A (4,4) = 4; A (4,9) = 11;

A (5,5) = 7;

b = [6000 10000 4000 7000 4000]';

Aeq = [1 1 -1 -1 -1 0 0 0 0;

              0 0 0 0 0 1 1 -1 0];

beq = [0 0]'；

lb = zeros(9,1);

[x fval] = linpro(c,A,b,Aeq,beq,lb)

fval = -fval

python代码解法：

import numpy as np
from scipy.optimize import linprog

# 初始化目标函数的系数向量
c = np.zeros(9)
c[0] = 1.25 - 0.25 - 300 / 6000 * 5
c[1] = 1.25 - 0.25 - 321 / 10000 * 7
c[2] = -250 / 4000 * 6
c[3] = -783 / 7000 * 4
c[4] = -200 / 4000 * 7
c[5] = -300 / 6000 * 10
c[6] = -321 / 10000 * 9
c[7] = 2 - 0.35 - 250 / 4000 * 8
c[8] = 2.8 - 0.5 - 321 / 10000 * 12 - 783 / 7000 * 11

# 不等式约束矩阵和右侧常数
A = np.array([
    [5, 0, 0, 0, 0, 10, 0, 0, 0],
    [0, 7, 0, 0, 0, 0, 9, 0, 12],
    [0, 0, 6, 0, 0, 0, 0, 8, 0],
    [0, 0, 0, 4, 0, 0, 0, 0, 11],
    [0, 0, 0, 0, 7, 0, 0, 0, 0]
])
b = np.array([6000, 10000, 4000, 7000, 4000])

# 等式约束矩阵和右侧常数
A_eq = np.array([
    [1, 1, -1, -1, -1, 0, 0, 0, 0],  # x1 + x2 - x3 - x4 - x5 = 0
    [0, 0, 0, 0, 0, 1, 1, -1, 0]     # x6 + x7 - x8 = 0
])

#b_eq = np.array([0, 0]).reshape(-1, 1)      #-1 表示这一维的大小自动计算，使总元素数量保持不变。1 表示这一维的大小为1（列向量）。
b_eq = np.array([[0],
                           [0]])

# 变量界限
x_bounds = [(0, None) for _ in range(9)]

# 求解线性规划问题
res_lp = linprog(-c, A_ub=A, b_ub=b, A_eq=A_eq, b_eq=b_eq, bounds=x_bounds, method='highs')

print("线性规划最优值:", -res_lp.fun)
print("线性规划最优解:", res_lp.x)