FlowShop流水线作业
在企业在面临大量多样化的生产任务时，如何合理地安排流水线作业以提高生产效率及确保交货期成为了一个重要的问题。
一个典型的问题就是FlowShop流水线作业安排问题,也有称为生产下料问题。它涉及到多台机器、多个工序以及多个作业的调度安排。在这个问题中，我们需要对多个作业在一组流水线上的处理顺序进行安排，以使得完成所有作业的总时间最短。
与FlowShop相似的还有JobShop调度问题，它们都涉及到多个作业在多台机器上的处理顺序安排，但FlowShop是顺序安排，处理顺序一致，JobShop是作业处理顺序可以不同，安排更灵活，调度约束相对复杂。通常JobShop调度问题通常比FlowShop调度问题更难求解。
本案例将着重介绍FlowShop问题，将向大家介绍如何运用数学规划方法来解决这一问题，帮助企业优化生产过程，实现生产目标。
1. 问题描述
有一个公司生产几套产品，这些产品各不一样，生产各个环节的耗时会不一致，即在每个工序上的耗时不一致。工厂的机器资源有限，如何合理安排产品的生产，让总生产耗时最短？
- 假设所有产品都设计成按同样顺序来经过这些加工工序。
- 各个产品在各个工序的耗时如下表，当有些工序不需要的时候，我们用耗时0来表示
| 产品ID | 产品任务1 | 产品任务2 | 产品任务3 | 产品任务4 | 产品任务5 | 产品任务6 | 产品任务7 | 产品任务8 | 产品任务9 | 产品任务10 | 产品任务11 | 产品任务12 | 产品任务13 | 产品任务14 | 产品任务15 | 产品任务16 | 产品任务17 | 产品任务18 | 产品任务19 | 产品任务20 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 工序机器1 | 54 | 83 | 15 | 71 | 77 | 36 | 53 | 38 | 27 | 87 | 76 | 91 | 14 | 29 | 12 | 77 | 32 | 87 | 68 | 94 |
| 工序机器2 | 79 | 3 | 11 | 99 | 56 | 70 | 99 | 60 | 5 | 56 | 3 | 61 | 73 | 75 | 47 | 14 | 21 | 86 | 5 | 77 |
| 工序机器3 | 16 | 89 | 49 | 15 | 89 | 45 | 60 | 23 | 57 | 64 | 7 | 1 | 63 | 41 | 63 | 47 | 26 | 75 | 77 | 40 |
| 工序机器4 | 66 | 58 | 31 | 68 | 78 | 91 | 13 | 59 | 49 | 85 | 85 | 9 | 39 | 41 | 56 | 40 | 54 | 77 | 51 | 31 |
| 工序机器5 | 58 | 56 | 20 | 85 | 53 | 35 | 53 | 41 | 69 | 13 | 86 | 72 | 8 | 49 | 47 | 87 | 58 | 18 | 68 | 28 |
- 加工完后，工序资源能立即释放出来去生产下一个产品。（此处时简化问题，把工序准备时间加一起了，实际中工序准备可以分开，同一个产品可以在上一个工序还没有生产完的时候就开始准备，即时间可以有重叠，本案例不讨论此情况）
- 每个工序同时处理的产品只有一种。
请问，如何安排生产，可以让最后完工的时间最小?
2. 数据准备
根据问题的描述，我们整理相应的数据为方便程序读的数据，如CSV的表格。
另外，这里我们做两份数据，方便大家理解怎么将自己业务数据整理后输入。
数据1：Raw数据：根据原问题的表格整理，prductionTime_Raw.csv文件
| 产品ID | 产品任务1 | 产品任务2 | 产品任务3 | 产品任务4 | 产品任务5 | 产品任务6 | 产品任务7 | 产品任务8 | 产品任务9 | 产品任务10 | 产品任务11 | 产品任务12 | 产品任务13 | 产品任务14 | 产品任务15 | 产品任务16 | 产品任务17 | 产品任务18 | 产品任务19 | 产品任务20 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 工序机器1 | 54 | 83 | 15 | 71 | 77 | 36 | 53 | 38 | 27 | 87 | 76 | 91 | 14 | 29 | 12 | 77 | 32 | 87 | 68 | 94 |
| 工序机器2 | 79 | 3 | 11 | 99 | 56 | 70 | 99 | 60 | 5 | 56 | 3 | 61 | 73 | 75 | 47 | 14 | 21 | 86 | 5 | 77 |
| 工序机器3 | 16 | 89 | 49 | 15 | 89 | 45 | 60 | 23 | 57 | 64 | 7 | 1 | 63 | 41 | 63 | 47 | 26 | 75 | 77 | 40 |
| 工序机器4 | 66 | 58 | 31 | 68 | 78 | 91 | 13 | 59 | 49 | 85 | 85 | 9 | 39 | 41 | 56 | 40 | 54 | 77 | 51 | 31 |
| 工序机器5 | 58 | 56 | 20 | 85 | 53 | 35 | 53 | 41 | 69 | 13 | 86 | 72 | 8 | 49 | 47 | 87 | 58 | 18 | 68 | 28 |
数据2：修改后数据，修改内容主要是：
- 前一天已经生产到一半的任务，已完成的工序耗时设为0；
- 不同产品生产工序有区别，将有些不需要执行的工序，耗时设置为0。
修改后的数据如下，同 prductionTime.csv文件
| 产品ID | 产品任务1 | 产品任务2 | 产品任务3 | 产品任务4 | 产品任务5 | 产品任务6 | 产品任务7 | 产品任务8 | 产品任务9 | 产品任务10 | 产品任务11 | 产品任务12 | 产品任务13 | 产品任务14 | 产品任务15 | 产品任务16 | 产品任务17 | 产品任务18 | 产品任务19 | 产品任务20 |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 工序机器1 | 0 | 0 | 0 | 0 | 77 | 36 | 53 | 38 | 27 | 87 | 76 | 91 | 14 | 29 | 12 | 77 | 32 | 87 | 68 | 94 |
| 工序机器2 | 0 | 0 | 0 | 99 | 56 | 70 | 99 | 60 | 5 | 56 | 0 | 61 | 73 | 75 | 47 | 14 | 21 | 86 | 5 | 77 |
| 工序机器3 | 0 | 0 | 49 | 15 | 89 | 0 | 60 | 23 | 57 | 64 | 7 | 0 | 63 | 41 | 63 | 47 | 26 | 75 | 0 | 40 |
| 工序机器4 | 0 | 58 | 31 | 68 | 78 | 0 | 13 | 59 | 49 | 85 | 85 | 9 | 39 | 41 | 56 | 40 | 54 | 77 | 51 | 31 |
| 工序机器5 | 58 | 56 | 20 | 85 | 53 | 0 | 53 | 41 | 69 | 13 | 86 | 72 | 8 | 49 | 47 | 0 | 58 | 18 | 68 | 28 |
3. 数学建模
这是个经典的FlowShop问题，最小化生产时间。决策的变量是不同产品/工件的生产顺序。
数学模型如下：
min
e
M
m
a
x
,
J
m
a
x
s.t.
∑
k
∈
J
x
i
,
k
=
1
∀
i
∈
J
∑
i
∈
J
x
i
,
k
=
1
∀
k
∈
J
e
m
,
k
≥
s
m
,
k
+
∑
i
∈
J
P
i
,
m
⋅
x
i
,
k
∀
m
∈
M
,
k
∈
J
s
m
,
k
+
1
≥
e
m
,
k
∀
m
∈
M
,
k
∈
{
1
,
2
,
.
.
,
J
m
a
x
−
1
}
s
m
+
1
,
k
≥
e
m
,
k
∀
k
∈
J
,
M
∈
{
1
,
2
,
.
.
,
M
m
a
x
−
1
}
s
1
,
1
=
0
s
m
,
k
,
e
m
,
k
≥
0
∀
m
∈
M
,
k
∈
J
x
i
,
k
∈
{
0
,
1
}
∀
k
,
i
∈
I
\begin{array}{rll} \min & e_{{M_{max},J_{max}}} & \\ \text{s.t.} & \sum_{k \in J} x_{i,k} = 1 & \forall i \in J \\ & \sum_{i \in J} x_{i,k} = 1 & \forall k \in J \\ & e_{m,k} \geq s_{m,k} + \sum_{i \in J} P_{i,m} \cdot x_{i,k} & \forall m \in M, k \in J\\ & s_{m,k+1} \geq e_{m,k} & \forall m \in M , k \in \{1,2,..,J_{max}-1\} \\ & s_{m+1,k} \geq e_{m,k} & \forall k \in J , M \in \{1,2,..,M_{max}-1\}\\ & s_{1,1}=0 &\\ & s_{m,k},e_{m,k}\geq 0 & \forall m \in M, k\in J \\ & x_{i,k}\in\{0,1\} & \forall k,i \in I \end{array}
mins.t.eMmax,Jmax∑k∈Jxi,k=1∑i∈Jxi,k=1em,k≥sm,k+∑i∈JPi,m⋅xi,ksm,k+1≥em,ksm+1,k≥em,ks1,1=0sm,k,em,k≥0xi,k∈{0,1}∀i∈J∀k∈J∀m∈M,k∈J∀m∈M,k∈{1,2,..,Jmax−1}∀k∈J,M∈{1,2,..,Mmax−1}∀m∈M,k∈J∀k,i∈I
更细节的数学公式请查看案例
4. MindOpt APL 建模和求解
MindOpt APL是一款代数建模语言，它可以方便地将数学语言描述成程序，然后调用多种求解器求解。MindOpt Solver支持混合整数线性规划（MILP）问题的求解，可选用它。
改写上面的数据图和数学模型，其中数据为了读取和索引加序号方便，我们把上面的表存储为kkv的长表单形式 prductionTime_long.csv和 prductionTime_long_Raw.csv，并把产品清单和工序清单名字以数字的形式分别存入：Products_name.csv、Machines_name.csv。
（请注意，由于我们里面需要用到m+1这样的顺序概念，因此建议把工序的ID用从1开始的顺序排的数字来记录，编程序更方便。）
代码解析
决策变量
var x[<i,k> in Jobs*Jobs] binary;
- 这是一个二进制决策变量用于表示作业i是否排在位置k上。
- 如果作业i排在位置k，则x[i,k]等于1，否则等于0。
- 注意这里使用了JobsJobs，正常情况下应该是JobsMachines或者其他表示作业与位置之间关系的集合。使用Jobs*Jobs可能是为了表示作业顺序而非作业和机器之间的关系。
var s[<m,k> in Machines * Jobs] >= 0;
- 这个变量’s’代表作业k在机器m上的开始时间。
- 约束条件保证’s’的值大于或等于0，表示不能有负的开始时间。
var e[<m,k> in Machines * Jobs] >= 0;
- 这个变量’e’代表作业k在机器m上的结束时间。
- 同样地，约束条件保证’e’的值大于或等于0，表示不能有负的结束时间。
决策目标
minimize TotalTime: e[lastMachines,lastJobs];
其中 e[m,k] 表示作业 k 在机器 m 上的结束时间。这里 lastMachines 和 lastJobs 分别指最后一台机器和最后一个作业。
约束
- 该约束条件表示每个产品必须在一个位置上进行处理。对于每个产品i，通过对其所有位置k上的处理决策进行求和，得到的结果必须等于1，表示每个产品只能在一个位置上进行处理。
subto constraint_1: forall { i in Jobs} sum {k in Jobs} x[i,k] == 1;
- 该约束条件表示每个位置只能处理一个产品。对于每个位置k，通过对其所有产品i的处理决策进行求和，得到的结果必须等于1，表示每个位置只能处理一个产品。
subto constraint_1_1: forall { k in Jobs} sum {i in Jobs} x[i,k] == 1;
- 该约束条件表示产品在每个工序的结束时间必须大于等于产品在该工序的开始时间加上该工序的耗时。对于每个工序m和产品k，通过对产品i在工序m上的处理决策乘以其耗时，求和得到的结果必须小于等于产品在该工序的结束时间。
subto constraint_2: forall {<m,k> in Machines * Jobs} e[m,k] >= s[m,k] + sum {i in Jobs} (ProductionTime[i,m] * x[i,k]);
- 该约束条件表示产品在每个工序的下一个位置上的开始时间必须大于等于产品在当前位置的结束时间。对于每个工序m和产品k，下一个位置k+1上的开始时间必须大于等于当前位置k上的结束时间。
subto constraint_3: forall {<m,k> in Machines * (Jobs without {lastJobs})} s[m,k+1] >= e[m,k];
- 该约束条件表示产品在下一个工序的同一位置上的开始时间必须大于等于产品在当前工序的结束时间。对于每个工序m和产品k，下一个工序m+1上的同一位置k的开始时间必须大于等于当前工序m上的结束时间。
subto constraint_3_1: forall {<m,k> in (Machines without {lastMachines}) * Jobs } s[m+1,k] >= e[m,k];
- 该约束条件表示第一个工序中第一个产品的开始时间必须为0，即从0时刻开始生产第一个产品。
subto constraint_4: s[1,1] == 0;
完整代码
代码第11和12行是切换不同的问题数据。
然后如下代码，在Notebook的cell中运行它，进行建模、求解，并整理存储结果：
clear model;
option modelname flowshop_demo; #定义存储文件名
# ----------建模--------Start----
# flowshop_production.mapl
# 读取数据
param fileDir = "data/";
set Jobs = {read fileDir+"Products_name.csv" as "<1n>" skip 1}; # 读取产品任务列表
set Machines = {read fileDir+"Machines_name.csv" as "<1n>" skip 1}; # 读取工序机器列表
#两种任务数据，可改注释#符号切换
param ProductionTime[Jobs*Machines] = read fileDir+"prductionTime_long_Raw.csv" as "<1n,2n> 3n" skip 1; #读取每个产品在每个工序的耗时
#param ProductionTime[Jobs*Machines] = read fileDir+"prductionTime_long.csv" as "<1n,2n> 3n" skip 1; #读取每个产品在每个工序的耗时
param lastJobs = card(Jobs);
param lastMachines = card(Machines);
print "总共有{}个产品任务，{}个工序机器" % lastJobs,lastMachines;
# 声明变量
var x[<i,k> in Jobs*Jobs] binary;
var s[<m,k> in Machines * Jobs] >= 0;
var e[<m,k> in Machines * Jobs] >= 0;
# 申明目标
minimize TotalTime: e[lastMachines,lastJobs];
# 约束
subto constraint_1:
forall { i in Jobs}
sum {k in Jobs} x[i,k] == 1;
subto constraint_1_1:
forall { k in Jobs}
sum {i in Jobs} x[i,k] == 1;
subto constraint_2:
forall {<m,k> in Machines * Jobs}
e[m,k] >= s[m,k] + sum {i in Jobs} (ProductionTime[i,m] * x[i,k]);
subto constraint_3:
forall {<m,k> in Machines * (Jobs without {lastJobs})}
s[m,k+1] >= e[m,k];
subto constraint_3_1:
forall {<m,k> in (Machines without {lastMachines}) * Jobs }
s[m+1,k] >= e[m,k];
subto constraint_4:
s[1,1] == 0;
#求解
option solver mindopt; # （可选）指定求解用的求解器，默认是MindOpt
#option mindopt_options 'max_time=600'; #设置求解器参数
#option mindopt_options 'print=0'; #设置求解器输出级别，减少过程打印
solve; # 求解
# 结果打印和检查结果
print "----------------结果打印和存文件--------------";
#display; #打印太多，注释了
print "最小生产耗时 = ",e[lastMachines,lastJobs];
运行上述代码结果如下：
总共有20个产品任务，5个工序机器
Running mindoptampl
wantsol=1
MindOpt Version 1.0.1 (Build date: 20231114)
Copyright (c) 2020-2023 Alibaba Cloud.
Start license validation (current time : 28-DEC-2023 20:40:53).
License validation terminated. Time : 0.006s
Model summary.
- Num. variables : 600
- Num. constraints : 315
- Num. nonzeros : 3350
- Num. integer vars. : 400
- Bound range : [1.0e+00,1.0e+00]
- Objective range : [1.0e+00,1.0e+00]
Branch-and-cut method started.
Original model: nrow = 315 ncol = 600 nnz = 3350
Tolerance: primal = 1e-06 int = 1e-06 mipgap = 0.0001 mipgapAbs = 1e-06
Limit: time = 1.79769313486232e+308 node = -1 stalling = -1 solution = -1
presolver terminated; took 0 ms
presolver terminated; took 163 ms
Parallelism: root=4, tree=4
accept new sol: obj 1448 bnd vio 1.11022302462516e-16 int vio 1.11022302462516e-16 mipgap 1 time 3
Model summary.
- Num. variables : 476
- Num. constraints : 192
- Num. nonzeros : 7883
- Bound range : [1.0e+00,5.4e+01]
- Objective range : [1.0e+00,8.7e+01]
- Matrix range : [1.0e+00,3.5e+02]
Presolver started.
Presolver terminated. Time : 0.002s
Simplex method started.
Model fingerprint: =Y2dmN2bgdnYgN2dm5mZ
Iteration Objective Dual Inf. Primal Inf. Time
0 1.45592e+02 0.0000e+00 4.0472e+01 0.00s
595 1.24863e+03 0.0000e+00 0.0000e+00 0.02s
Postsolver started.
Simplex method terminated. Time : 0.018s
Root relaxation: 1248.62782173302 iteration = 595 time = 0.02
#node(P:0 Q:0) #(dual:1248.63 best:1448 gap:13.77%) #time = 3
accept new sol: obj 1410 bnd vio 2.8421709430404e-14 int vio 2.22044604925031e-16 mipgap 0.114448353380835 time 3
accept new sol: obj 1385 bnd vio 1.95761774137987e-14 int vio 1.45008721583694e-16 mipgap 0.098463666618756 time 3
accept new sol: obj 1376 bnd vio 2.22044604925031e-16 int vio 2.22044604925031e-16 mipgap 0.0925669900196054 time 3
accept new sol: obj 1349 bnd vio 4.2632564145606e-14 int vio 2.22044604925031e-16 mipgap 0.0744048764025033 time 3
accept new sol: obj 1344 bnd vio 2.8421709430404e-13 int vio 1.11022302462516e-16 mipgap 0.0709614421629293 time 4
accept new sol: obj 1304 bnd vio 2.8421709430404e-13 int vio 9.99200722162641e-16 mipgap 0.0424633268918535 time 5
accept new sol: obj 1291 bnd vio 3.33066907387547e-16 int vio 3.33066907387547e-16 mipgap 0.03282120702322 time 7
accept new sol: obj 1290 bnd vio 2.22044604925031e-16 int vio 2.22044604925031e-16 mipgap 0.0320714560209124 time 13
accept new sol: obj 1278 bnd vio 2.79866832656529e-14 int vio 3.33066907387547e-16 mipgap 0.0223161708476798 time 24
Branch-and-cut method terminated. Time : 31.571s
OPTIMAL; objective 1278.00
Completed.
----------------结果打印和存文件--------------
最小生产耗时 = 1278
产品,序号,……,|第1个工序s,耗时,|第2个工序s,耗时,|第3个工序s,耗时,|第4个工序s,耗时,|第5个工序s,耗时,|最后结束时间
1,9,……,|343,54,|397,79,|488,16,|562,66,|645,58,|703
2,8,……,|260,83,|394,3,|399,89,|504,58,|589,56,|645
3,11,……,|474,15,|532,11,|636,49,|714,31,|759,20,|779
4,7,……,|189,71,|285,99,|384,15,|436,68,|504,85,|589
5,10,……,|397,77,|476,56,|532,89,|628,78,|706,53,|759
6,4,……,|71,36,|129,70,|199,45,|255,91,|362,35,|397
7,12,……,|489,53,|543,99,|685,60,|745,13,|779,53,|832
8,15,……,|647,38,|722,60,|799,23,|884,59,|1001,41,|1042
9,2,……,|32,27,|74,5,|79,57,|150,49,|199,69,|268
10,18,……,|849,87,|942,56,|998,64,|1062,85,|1147,13,|1160
11,13,……,|542,76,|642,3,|751,7,|758,85,|866,86,|952
12,20,……,|1030,91,|1135,61,|1196,1,|1197,9,|1206,72,|1278
13,6,……,|175,14,|212,73,|321,63,|397,39,|496,8,|504
14,14,……,|618,29,|647,75,|758,41,|843,41,|952,49,|1001
15,3,……,|59,12,|82,47,|136,63,|199,56,|268,47,|315
16,17,……,|772,77,|928,14,|951,47,|1020,40,|1060,87,|1147
17,1,……,|0,32,|32,21,|53,26,|79,54,|133,58,|191
18,16,……,|685,87,|782,86,|868,75,|943,77,|1042,18,|1060
19,5,……,|107,68,|207,5,|244,77,|346,51,|397,68,|465
20,19,……,|936,94,|1030,77,|1107,40,|1147,31,|1178,28,|1206
