线性规划-pulp-复杂矩阵

最新推荐文章于 2024-09-14 11:14:07 发布

cij

最新推荐文章于 2024-09-14 11:14:07 发布

阅读量5.1k

点赞数 7

分类专栏： le 文章标签：线性规划 pulp python 复杂操作

本文链接：https://blog.csdn.net/u010277233/article/details/89071564

版权

le 专栏收录该内容

1 篇文章 0 订阅

订阅专栏

线性规划-pulp-复杂矩阵

1. 简介-线性规划

在数学中，线性规划（Linear Programming，简称LP）特指目标函数和约束条件皆为线性的最优化问题。线性规划是最优化问题中的一个重要领域。在作业研究中所面临的许多实际问题都可以用线性规划来处理，特别是某些特殊情况。最简单情况就，解线性方程组。

举个最简单的例子:

工厂生产A和B两种物品，需要原料配比分别是1：2：1和2：1：3.

现有原料是100，40，8。A和B的价格分别是5和10。

请问如何分配原材料才能获益最大？

求解一个线性方程组即可

2. pulp 包

Python中有很多的工具可以解决线性规划问题。由于pulp因其常用和方便性，故这里介绍pulp包的使用。pulp能够解决大部分的线性规划问题，如""整数规划（Integer Program）、01规划（Binary Program）还是混合整数线性规划（MILP），都能很好的解决。
pulp的github网址和官方文档，说实话官方文档做的并不好，但是里面介绍的例子还是很好的，如果只做最简单的线性规划可以看看学习一下。
安装：按照一般套路pip install pulp 即可。

3. pulp使用流程

获得问题的描述：求解最大值还是最小值，这个要明确了。输入是的数据是什么。
确定需要的决策变量：需要求解的变量，定义好，大小，数据类型。
制定约束条件：决策变量和输入数据之间存在什么样的约束条件。
输入数据，求解得到结果。

那个实际例子来说明：

获得问题的描述：我们有矩阵X和矩阵Y，大小是1008*6。如图最上面的例子：是1008个产品，在6+6个指标上的需求。拆分成6+6两部分可以说明更加复杂的关系，以便更具有扩展性。

目标，我们要求最大值： $\max \sum_{j=1}^{n} \varphi_{j}=\sum_{j=1}^{n} \sum_{r=1}^{6} \mu_{j r} y_{j r}$ 其中 $n = 1008$ ，$ \mu_{j r} $是要求解的决策变量
制定约束条件如下：
1. $\sum_{i=1}^{6} v_{j i} x_{j i}=1, \forall j$ ，$v_{j i} $是要求解的决策变量
2. $\sum_{r=1}^{6} \mu_{j r} y_{k r}-\sum_{i=1}^{6} v_{j i} x_{k i} \leq 0, \forall j, k$ ，$ \mu_{j r} $和$ v_{j i} $是要求解的决策变量
3. $\sum_{j=1}^{n} v_{j i}=\rho_{i}, i=1, \ldots, 6$ ， $\rho_{i}$ 是要求解的决策变量
4. $\sum_{j=1}^{n} \mu_{j r}=5*\rho_{r}, r=1, \ldots, 6$ ， $\rho_{r}$ 是要求解的决策变量
5. $0.15<=\rho_{r}<=0.18$ ……
6. $ \mu_{j r} >=0 $和$ v_{j i}>=0 $
输入数据矩阵X和Y进行求解，即可。

4. 转换为代码求解

废话不多说，直接上代码。

# -*- coding: UTF-8 -*-
import pulp
import numpy as np

def assignment_problem(matrix_x,matrix_y):
    row_matrix_x = len(matrix_x) 
    col_matrix_x = len(matrix_x[0]) 

    row_matrix_y = len(matrix_x)
    col_matrix_y = len(matrix_x[0])
    
    #这里是将变量拉成一条直线，方便求解
    flatten = lambda x: [y for l in x for y in flatten(l)] if type(x) is list else [x] 
    
    #这是说明了是求最大值：LpMaximize，如果需要求最小值：LpMinimize。
    m = pulp.LpProblem('task_name', sense=pulp.LpMaximize)    

    #这里是mu和v2个决策变量的定义，是个矩阵形式，所以需要使用的dict这个函数，lowBound说明了下界大于等于0，和约束条件6对应
    #这里的数据类型默认是浮点型，也可以换成Ingredients，cat=’Ingredients‘就是整数， cat=pulp.LpBinary这是离散型数据就是整数了。
    var_mu = pulp.LpVariable.dicts("mu", (range(row_matrix_x), range(col_matrix_x)), lowBound = 0) 
    var_vi = pulp.LpVariable.dicts("vi", (range(row_matrix_x), range(col_matrix_x)), lowBound = 0) 
  
    #这里是6个变量，由于约束条件5中说明了上下届
    rho1 = pulp.LpVariable('rho1', lowBound = 0.15, upBound = 0.18) 
    rho2 = pulp.LpVariable('rho2', lowBound = 0.02, upBound = 0.17) 
    rho3 = pulp.LpVariable('rho3', lowBound = 0.15, upBound = 0.33) 
    rho4 = pulp.LpVariable('rho4', lowBound = 0.25, upBound = 0.33) 
    rho5 = pulp.LpVariable('rho5', lowBound = 0.08, upBound = 0.17) 
    rho6 = pulp.LpVariable('rho6', lowBound = 0.01, upBound = 0.18)  
    

    #对应这求解最大值，即目标。pulp.lpDot相乘操作
    m += pulp.lpDot(matrix_y.flatten(), flatten(var_mu)) 
   
    #对应这约束条件1，pulp.lpSum当成求和用即可。
    for j in range(row_matrix_x):
        m +=(pulp.lpSum(var_vi[j][i]*matrix_x[j][i] for i in range(col_matrix_x))==1)  
    
    #对应这约束条件2，看起来有点长。
    for j in range(row_matrix_x):
        for k in range(row_matrix_x): 
            m+=(pulp.lpSum(var_mu[j][i]*matrix_y[k][i] for i in range(col_matrix_x))-pulp.lpSum(var_vi[j][i]*matrix_x[k][i] for i in  range(col_matrix_x))<=0.0)
    
    #对应这约束条件3和4，为了方便就展开来写了。也可以改成上述的for循环。
    m += (pulp.lpSum(var_mu[j][0] for j in range(row_matrix_x))==rho1*row_matrix_x) 
    m += (pulp.lpSum(var_mu[j][1] for j in range(row_matrix_x))==rho2*row_matrix_x) 
    m += (pulp.lpSum(var_mu[j][2] for j in range(row_matrix_x))==rho3*row_matrix_x) 
    m += (pulp.lpSum(var_mu[j][3] for j in range(row_matrix_x))==rho4*row_matrix_x) 
    m += (pulp.lpSum(var_mu[j][4] for j in range(row_matrix_x))==rho5*row_matrix_x) 
    m += (pulp.lpSum(var_mu[j][5] for j in range(row_matrix_x))==rho6*row_matrix_x)  

    m += (pulp.lpSum(var_vi[j][0] for j in range(row_matrix_x))==row_matrix_x*rho1*5)
    m += (pulp.lpSum(var_vi[j][1] for j in range(row_matrix_x))==row_matrix_x*rho2*5)
    m += (pulp.lpSum(var_vi[j][2] for j in range(row_matrix_x))==row_matrix_x*rho3*5)
    m += (pulp.lpSum(var_vi[j][3] for j in range(row_matrix_x))==row_matrix_x*rho4*5)
    m += (pulp.lpSum(var_vi[j][4] for j in range(row_matrix_x))==row_matrix_x*rho5*5)
    m += (pulp.lpSum(var_vi[j][5] for j in range(row_matrix_x))==row_matrix_x*rho6*5)  
    
    #求解目标，就这么简单。不用怀疑
    m.solve()
    #可以打印出来各个约束条件，拆分开的很多，可以注释掉。
    print(m)
    
    #把求解的值，取出来的操作。 
    result_mu=[[pulp.value(var_mu[i][j]) for j in range(col_matrix_x)] for i in range(row_matrix_x)]
    result_vi=[[pulp.value(var_vi[i][j]) for j in range(col_matrix_x)] for i in range(row_matrix_x)]
    result_rho1=pulp.value(rho1)
    result_rho2=pulp.value(rho2)
    result_rho3=pulp.value(rho3)
    result_rho4=pulp.value(rho4)
    result_rho5=pulp.value(rho5)
    result_rho6=pulp.value(rho6)
 
    return {'objective':pulp.value(m.objective), 'result_mu':result_mu,'result_vi':result_vi}
 
#构造数据，我们这里为了方便直接写0，自己进行替换即可。
matrix_x=np.zeros((1008,6))
matrix_y=np.zeros((1008,6))
#输入数据求解，打印结果.  
res = assignment_problem(matrix_x,matrix_y) 
print('{res["objective"]}')
print(res['result_mu'])

5. 思考

在代码中我们在m上加上不同的约束添加，看起来不可思议。pulp为了能够这么方便的使用，使用了操作符重载的操作。感谢大佬，让线性求解如此简单。
pulp是做线性规划的，所以不支持2个决策变量相乘。例如： $\rho_{1}*\rho_{2}$ 是不支持。查阅了很多发现真的不支持，这样好像就是非线性规划了，所以做不了。
关系式中>=，<=。但是不支持>和<，请注意。
一些参考链接，感谢下列的博客。

例子很多，适合小白开始，这里举了常见的三个例子，如果不需要我上面这么复杂可以先看这个链接，当然其中代码存在一点问题，就是定义决策变量如果是个矩阵会出错，可能博客中是以前的版本吧。如果需要定义决策变量请参考我的代码，所有的都可以。

官网数独求解这个例子很有意思。