置换流水车间调度问题的MATLAB求解

物流运筹实务课程设计题目：置换流水车间调度问题的MATLAB求解置换流水车间调度问题的MATLAB求解目录一、前言 (5)二、问题描述 (6)三、算法设计 (7)四、实验结果 (15)摘要自从Johnson 1954年发表第一篇关于流水车间调度问题的文章以来．流水车间调度问题引起了许多学者的关注。

安排合理有效的生产调度是生产活动能井然有序开展，生产资源得到最佳配置，运作过程简明流畅的有力保证。

流水车间调度问题是许多实际流水线生产调度问题的简化模型。

它无论是在离散制造工业还是在流程工业中都具有广泛的应用。

因此，对进行研究具有重要的理论意义和工程价值。

流水线调度问题中一个非常典型的问题，而置换流水线调度问题作为FSP 问题的子问题，是一个著名的组合优化问题。

该问题是一个典型的NP难问题，也是生产管理的核心内容。

随着生产规模的扩大，流水线调度问题的优化对提高资源利用率的作用越来越大，因此对其研究具有重要的理论和现实意义。

关键字：流水车间，单件小批量生产，jsp模型，Matlab前言企业资源的合理配置和优化利用很大程度上体现在车间一层的生产活动中,所以加强车间层的生产计划与控制一直在企业生产经营活动中占有十分重要的地位。

车间生产计划与控制的核心理论是调度理论。

车间调度问题是一类重要的组合优化问题。

为适应订货式、多品种、小批量生产的需要,引进了置换流水车间调度概念。

在置换流水车间调度优化后,可以避免或大大减少流程工作时间、提高生产效率。

因此,研究成组技术下车间调度问题是很有必要的。

生产调度，即对生产过程进行作业计划，是整个个先进生产制造系统实现管理技术、优化技术、白动化与计算机技术发展的核心。

置换流水车间调度问题是许多实际生产调度问题的简化模型。

生产计划与调度直接关系着企业的产出效率和生产成本，有效的计划与调度算法能最大限度地提高企业的效益。

调度问题是组合优化问题，属于NP问题，难以用常规力一法求解。

随着制造业的快速发展，大规模定制生产、全球化制造等思想的提出，使车间调度问题呈现出以下的新特点:约束条件多，时间复杂度高，空问复杂度高。

这将导致在许多情况下，求解所建立的数学模型的快速性无法满足，如果采用适度线形化处理之后求解，将会因简化太多而使结果严承失真。

所以需选择功能强大的数值计算工具来实现这一问题的求解。

MATLAB恰好提供了这样的平台。

MATLAB是一个高度集成的系统，集科学计算、图像处理、声音处理于一体，具有极高的编程效率。

典型JSP模型分析与Matlab的应用结合使流水车间调度问题迎刃而解。

最大完工时间是生产调度中最常用的性能度量指标之一，最大完工时间越短，则说明产品总的生产周期越短，生产能力越大；此类调度问题的优化研究有助于提高企业的生产效率与资源利用率。

一、问题描述流水车间调度问题通常可以描述为n个工件要在m台机器上加工，每个工件有m道工序，每道工序都要在不同的机器上加工，所有工件的加工顺序都相同，问题的目标是确定每台机器上工件的加工顺序及开工时间，使得特定的性能指标最优。

置换流水车间调度问题PFSP是对流水车间调度问题的进一步约束，即约定每台机器上所有工件的加工顺序相同，其解空间的规模为!n，远远小于流水车间调n。

度问题的规模(!)m本次课程实验主要研究PFSP中的最小化最大完工时间问题，利用三元组表示F prmu C）求解Carlier (1978)提出的8个算例、以及Reeves (1995)提法（max出的21个算。

由于三台机器以上的调度问题被证明是NP难问题，对于大规模的调度，至今仍未出现求解最优的方法，常常采用启发式算法来求解近优解。

本案例主要采用instance car2进行求解。

案例：某产品，需要经过4道工序对13个工件进行加工，这13个工件的生产流程是一样的。

加工时间表见下：表4-3 某产品加工时间表计算步骤如下：首先确定n/m/F/C max 的最大完工时间为： 1,11)1,(c j t j =k j t k j c i j c 1)1,(),(11+-= k=2,...,m 111)1,(c )1,(c i j i t j j +=- i=2,...,n k j i i i i t k j c k j k j +-=-)}1,();,(c m ax {),(c 1则 C max =),(c n m j二、算法设计 (一)假设工件在机器上的加工顺序是相同的，同时假定各工件准备就绪，机器一开动就投入生产，开工时间为0，则最大完工时间等于最大流程时间。

同时3台机器以上的流水车间调度是NP 难问题，所以本文只考虑了2台、3台机器的情况，解决3台机器以上的问题方法也可运用人工智能算法，解的质量更高，但因该类算法需良好的软件编程能力，故本文不加探究。

n 个工件在m 台机器上的加工顺序相同。

工件在机器上的加工时间是给定的。

问题的目标是求n 个工件在每合机器上的最大完工时间等于最大流程时间。

这种流水线调度问题要在满足以下两个约束条件的前提下，使得加工完所有的工件所花的时间尽可能地少： 1、工件约束每个工件在每台机器上恰好加工一次，每个工件在各机器上加工顺序相同。

不失一般性，假设各工件按机器1至m 的顺序进行加工。

各工件在各机器上的加工时间已知。

2、机器约束每台机器在任何时刻至多加工一个工件，每台机器加工的各工件的顺序相同。

置换流水线调度问题实质是如何调整加工工件的序列，提高机器的利用率的问题，即在同一时刻正在加工的机攫数越多，机器利用率越大口根据该原则，我们根据下面规则安排工件的加工顺序:(l)在前面机器加工时间较短、后面机器加工时间较长的工件，安排在序列前。

这样可以使得后面的机器尽快参加工作，并且后面的机器不需要作空等待，(2)机器加工时间较为平均且加工时间较长的工件，安排在序列的中部。

这样可以使得各个机器在中期的时候都能得到运作。

(3〕前面加工时间较长，后面加一〔时间较短的上件女排在序列尾部。

这样使得前面的机器能“延迟”完工，后面的机器尽快完工。

（二）利用Matlab软件对上面的案例进行求解，编程如下:软件输出相应的结果，如下：（三）绘制甘特图编程如下：三、实验结果1、根据上面matlab的求解得到以下实验结果：最优排序为13 4 12 6 7 10 1 14 8 9 3 11 25 min(Cmax)= 84232、甘特图如下所示：四、流水线型车间作业调度问题遗传算法MATLAB源码流水线型车间作业调度问题可以描述如下：n个任务在流水线上进行m个阶段的加工，每一阶段至少有一台机器且至少有一个阶段存在多台机器，并且同一阶段上各机器的处理性能相同，在每一阶段各任务均要完成一道工序，各任务的每道工序可以在相应阶段上的任意一台机器上加工，已知任务各道工序的处理时间，要求确定所有任务的排序以及每一阶段上机器的分配情况，使得调度指标(一般求Makespan)最小。

下面的源码是求解流水线型车间作业调度问题的遗传算法通用MATLAB源码function [Zp,Y1p,Y2p,Y3p,Xp,LC1,LC2]=JSPGA(M,N,Pm,T,P)% 输入参数列表% M 遗传进化迭代次数% N 种群规模(取偶数)% Pm 变异概率% T m×n的矩阵，存储m个工件n个工序的加工时间% P 1×n的向量，n个工序中，每一个工序所具有的机床数目% 输出参数列表% Zp 最优的Makespan值% Y1p 最优方案中，各工件各工序的开始时刻，可根据它绘出甘特图% Y2p 最优方案中，各工件各工序的结束时刻，可根据它绘出甘特图% Y3p 最优方案中，各工件各工序使用的机器编号% Xp 最优决策变量的值，决策变量是一个实数编码的m×n矩阵% LC1 收敛曲线1，各代最优个体适应值的记录% LC2 收敛曲线2，各代群体平均适应值的记录% 最后，程序还将绘出三副图片：两条收敛曲线图和甘特图（各工件的调度时序图）%第一步：变量初始化[m,n]=size(T);%m是总工件数，n是总工序数Xp=zeros(m,n);%最优决策变量LC1=zeros(1,M);%收敛曲线1LC2=zeros(1,N);%收敛曲线2%第二步：随机产生初始种群farm=cell(1,N);%采用细胞结构存储种群for k=1:NX=zeros(m,n);for j=1:nfor i=1:mX(i,j)=1+(P(j)-eps)*rand;endendfarm{k}=X;endcounter=0;%设置迭代计数器while counter<M%停止条件为达到最大迭代次数%第三步：交叉newfarm=cell(1,N);%交叉产生的新种群存在其中Ser=randperm(N);for i=1:2:(N-1)A=farm{Ser(i)};%父代个体Manner=unidrnd(2);%随机选择交叉方式if Manner==1cp=unidrnd(m-1);%随机选择交叉点%双亲双子单点交叉a=[A(1:cp,:);B((cp+1):m,:)];%子代个体b=[B(1:cp,:);A((cp+1):m,:)];elsecp=unidrnd(n-1);%随机选择交叉点b=[B(:,1:cp),A(:,(cp+1):n)];endnewfarm{i}=a;%交叉后的子代存入newfarmnewfarm{i+1}=b;end%新旧种群合并FARM=[farm,newfarm];%第四步：选择复制FITNESS=zeros(1,2*N);fitness=zeros(1,N);plotif=0;for i=1:(2*N)X=FARM{i};Z=COST(X,T,P,plotif);%调用计算费用的子函数FITNESS(i)=Z;end%选择复制采取两两随机配对竞争的方式，具有保留最优个体的能力Ser=randperm(2*N);for i=1:Nf2=FITNESS(Ser(2*i));if f1<=f2farm{i}=FARM{Ser(2*i-1)};fitness(i)=FITNESS(Ser(2*i-1));elsefarm{i}=FARM{Ser(2*i)};endend%记录最佳个体和收敛曲线minfitness=min(fitness)meanfitness=mean(fitness)LC1(counter+1)=minfitness;%收敛曲线1，各代最优个体适应值的记录 LC2(counter+1)=meanfitness;%收敛曲线2，各代群体平均适应值的记录 pos=find(fitness==minfitness);Xp=farm{pos(1)};%第五步：变异for i=1:Nif Pm>rand;%变异概率为PmX=farm{i};I=unidrnd(m);J=unidrnd(n);X(I,J)=1+(P(J)-eps)*rand;farm{i}=X;endendfarm{pos(1)}=Xp;counter=counter+1end%输出结果并绘图figure(1);plotif=1;X=Xp;[Zp,Y1p,Y2p,Y3p]=COST(X,T,P,plotif);figure(2);plot(LC1);figure(3);plot(LC2);function [Zp,Y1p,Y2p,Y3p]=COST(X,T,P,plotif)% JSPGA的内联子函数，用于求调度方案的Makespan值% 输入参数列表% X 调度方案的编码矩阵，是一个实数编码的m×n矩阵% T m×n的矩阵，存储m个工件n个工序的加工时间% P 1×n的向量，n个工序中，每一个工序所具有的机床数目% plotif 是否绘甘特图的控制参数% 输出参数列表% Zp 最优的Makespan值% Y1p 最优方案中，各工件各工序的开始时刻% Y2p 最优方案中，各工件各工序的结束时刻% Y3p 最优方案中，各工件各工序使用的机器编号%第一步：变量初始化[m,n]=size(X);Y1p=zeros(m,n);Y2p=zeros(m,n);Y3p=zeros(m,n);%第二步：计算第一道工序的安排Q1=zeros(m,1);Q2=zeros(m,1);R=X(:,1);%取出第一道工序Q3=floor(R);%向下取整即得到各工件在第一道工序使用的机器的编号%下面计算各工件第一道工序的开始时刻和结束时刻for i=1:P(1)%取出机器编号pos=find(Q3==i);%取出使用编号为i的机器为其加工的工件的编号lenpos=length(pos);if lenpos>=1Q1(pos(1))=0;if lenpos>=2for j=2:lenposQ1(pos(j))=Q2(pos(j-1));Q2(pos(j))=Q2(pos(j-1))+T(pos(j),1);endendendendY1p(:,1)=Q1;Y3p(:,1)=Q3;%第三步：计算剩余工序的安排for k=2:nR=X(:,k);%取出第k道工序Q3=floor(R);%向下取整即得到各工件在第k道工序使用的机器的编号%下面计算各工件第k道工序的开始时刻和结束时刻for i=1:P(k)%取出机器编号pos=find(Q3==i);%取出使用编号为i的机器为其加工的工件的编号lenpos=length(pos);if lenpos>=1EndTime=Y2p(pos,k-1);%取出这些机器在上一个工序中的结束时刻POS=zeros(1,lenpos);%上一个工序完成时间由早到晚的排序for jj=1:lenposPOS(jj)=ppp(1);EndTime(ppp(1))=Inf;end%根据上一个工序完成时刻的早晚，计算各工件第k道工序的开始时刻和结束时刻Q1(pos(POS(1)))=Y2p(pos(POS(1)),k-1);Q2(pos(POS(1)))=Q1(pos(POS(1)))+T(pos(POS(1)),k);%前一个工件的结束时刻if lenpos>=2for j=2:lenposQ1(pos(POS(j)))=Y2p(pos(POS(j)),k-1);%预定的开始时刻为上一个工序的结束时刻if Q1(pos(POS(j)))<Q2(pos(POS(j-1)))%如果比前面的工件的结束时刻还早Q1(pos(POS(j)))=Q2(pos(POS(j-1)));endendendendendY1p(:,k)=Q1;Y2p(:,k)=Q2;Y3p(:,k)=Q3;end%第四步：计算最优的Makespan值Y2m=Y2p(:,n);Zp=max(Y2m);%第五步：绘甘特图if plotiffor i=1:mfor j=1:nmPoint1=Y1p(i,j);mPoint2=Y2p(i,j);mText=m+1-i;PlotRec(mPoint1,mPoint2,mText);Word=num2str(Y3p(i,j));%text(0.5*mPoint1+0.5*mPoint2,mText-0.5,Word);hold onx1=mPoint1;y1=mText-1;x2=mPoint2;y2=mText-1;x4=mPoint1;y4=mText;%fill([x1,x2,x3,x4],[y1,y2,y3,y4],'r');fill([x1,x2,x3,x4],[y1,y2,y3,y4],[1,0.5,1]);text(0.5*mPoint1+0.5*mPoint2,mText-0.5,Word);endendendfunction PlotRec(mPoint1,mPoint2,mText)% 此函数画出小矩形% 输入:% mPoint1 输入点1,较小,横坐标% mPoint2 输入点2,较大,横坐标% mText 输入的文本,序号,纵坐标vPoint = zeros(4,2) ;vPoint(1,:) = [mPoint1,mText-1];vPoint(2,:) = [mPoint2,mText-1];vPoint(3,:) = [mPoint1,mText];vPoint(4,:) = [mPoint2,mText];plot([vPoint(1,1),vPoint(2,1)],[vPoint(1,2),vPoint(2,2)]);hold on ;plot([vPoint(1,1),vPoint(3,1)],[vPoint(1,2),vPoint(3,2)]);plot([vPoint(2,1),vPoint(4,1)],[vPoint(2,2),vPoint(4,2)]);plot([vPoint(3,1),vPoint(4,1)],[vPoint(3,2),vPoint(4,2)]);参考文献(1)《MATLAB基础教程》清华大学出版社 2011.3 第一版(2)《计算机工程与应用》 2009年第36期。