【原创】车间作业调度问题遗传算法通用Matlab程序（附图）

aiwa

GreenSim团队的博客开通了，更多原创程序请访问我们的博客
http://blog.sina.com.cn/greensim

1. function [Zp,Y1p,Y2p,Y3p,Xp,LC1,LC2]=JSPGA(M,N,Pm,T,P)
   
2. %--------------------------------------------------------------------------
   
3. % JSPGA.m
   
4. % 车间作业调度问题遗传算法
   
5. %--------------------------------------------------------------------------
   
6. % 输入参数列表
   
7. % M 遗传进化迭代次数
   
8. % N 种群规模(取偶数)
   
9. % Pm 变异概率
   
10. % T m×n的矩阵，存储m个工件n个工序的加工时间
    
11. % P 1×n的向量，n个工序中，每一个工序所具有的机床数目
    
12. % 输出参数列表
    
13. % Zp 最优的Makespan值
    
14. % Y1p 最优方案中，各工件各工序的开始时刻，可根据它绘出甘特图
    
15. % Y2p 最优方案中，各工件各工序的结束时刻，可根据它绘出甘特图
    
16. % Y3p 最优方案中，各工件各工序使用的机器编号
    
17. % Xp 最优决策变量的值，决策变量是一个实数编码的m×n矩阵
    
18. % LC1 收敛曲线1，各代最优个体适应值的记录
    
19. % LC2 收敛曲线2，各代群体平均适应值的记录
    
20. % 最后，程序还将绘出三副图片：两条收敛曲线图和甘特图（各工件的调度时序图）
    
21. 
    
22. %第一步：变量初始化
    
23. [m,n]=size(T);%m是总工件数，n是总工序数
    
24. Xp=zeros(m,n);%最优决策变量
    
25. LC1=zeros(1,M);%收敛曲线1
    
26. LC2=zeros(1,N);%收敛曲线2
    
27. 
    
28. %第二步：随机产生初始种群
    
29. farm=cell(1,N);%采用细胞结构存储种群
    
30. for k=1:N
    
31. X=zeros(m,n);
    
32. for j=1:n
    
33. for i=1:m
    
34. X(i,j)=1+(P(j)-eps)*rand;
    
35. end
    
36. end
    
37. farm{k}=X;
    
38. end
    
39. 
    
40. counter=0;%设置迭代计数器
    
41. while counter
    
42. 
    
43. %第三步：交叉
    
44. newfarm=cell(1,N);%交叉产生的新种群存在其中
    
45. Ser=randperm(N);
    
46. for i=1:2:(N-1)
    
47. A=farm{Ser(i)};%父代个体
    
48. B=farm{Ser(i+1)};
    
49. Manner=unidrnd(2);%随机选择交叉方式
    
50. if Manner==1
    
51. cp=unidrnd(m-1);%随机选择交叉点
    
52. %双亲双子单点交叉
    
53. a=[A(1:cp,:);B((cp+1):m,:)];%子代个体
    
54. b=[B(1:cp,:);A((cp+1):m,:)];
    
55. else
    
56. cp=unidrnd(n-1);%随机选择交叉点
    
57. a=[A(:,1:cp),B(:,(cp+1):n)];%双亲双子单点交叉
    
58. b=[B(:,1:cp),A(:,(cp+1):n)];
    
59. end
    
60. newfarm{i}=a;%交叉后的子代存入newfarm
    
61. newfarm{i+1}=b;
    
62. end
    
63. %新旧种群合并
    
64. FARM=[farm,newfarm];
    
65. 
    
66. %第四步：选择复制
    
67. FITNESS=zeros(1,2*N);
    
68. fitness=zeros(1,N);
    
69. plotif=0;
    
70. for i=1:(2*N)
    
71. X=FARM{i};
    
72. Z=COST(X,T,P,plotif);%调用计算费用的子函数
    
73. FITNESS(i)=Z;
    
74. end
    
75. %选择复制采取两两随机配对竞争的方式，具有保留最优个体的能力
    
76. Ser=randperm(2*N);
    
77. for i=1:N
    
78. f1=FITNESS(Ser(2*i-1));
    
79. f2=FITNESS(Ser(2*i));
    
80. if f1<=f2
    
81. farm{i}=FARM{Ser(2*i-1)};
    
82. fitness(i)=FITNESS(Ser(2*i-1));
    
83. else
    
84. farm{i}=FARM{Ser(2*i)};
    
85. fitness(i)=FITNESS(Ser(2*i));
    
86. end
    
87. end
    
88. %记录最佳个体和收敛曲线
    
89. minfitness=min(fitness)
    
90. meanfitness=mean(fitness)
    
91. LC1(counter+1)=minfitness;%收敛曲线1，各代最优个体适应值的记录
    
92. LC2(counter+1)=meanfitness;%收敛曲线2，各代群体平均适应值的记录
    
93. pos=find(fitness==minfitness);
    
94. Xp=farm{pos(1)};
    
95. 
    
96. %第五步：变异
    
97. for i=1:N
    
98. if Pm>rand;%变异概率为Pm
    
99. X=farm{i};
    
100. I=unidrnd(m);
     
101. J=unidrnd(n);
     
102. X(I,J)=1+(P(J)-eps)*rand;
     
103. farm{i}=X;
     
104. end
     
105. end
     
106. farm{pos(1)}=Xp;
     
107. 
     
108. counter=counter+1
     
109. end
     
110. 
     
111. %输出结果并绘图
     
112. figure(1);
     
113. plotif=1;
     
114. X=Xp;
     
115. [Zp,Y1p,Y2p,Y3p]=COST(X,T,P,plotif);
     
116. figure(2);
     
117. plot(LC1);
     
118. figure(3);
     
119. plot(LC2);
     
120. 
     
121. function [Zp,Y1p,Y2p,Y3p]=COST(X,T,P,plotif)
     
122. % JSPGA的内联子函数，用于求调度方案的Makespan值
     
123. % 输入参数列表
     
124. % X 调度方案的编码矩阵，是一个实数编码的m×n矩阵
     
125. % T m×n的矩阵，存储m个工件n个工序的加工时间
     
126. % P 1×n的向量，n个工序中，每一个工序所具有的机床数目
     
127. % plotif 是否绘甘特图的控制参数
     
128. % 输出参数列表
     
129. % Zp 最优的Makespan值
     
130. % Y1p 最优方案中，各工件各工序的开始时刻
     
131. % Y2p 最优方案中，各工件各工序的结束时刻
     
132. % Y3p 最优方案中，各工件各工序使用的机器编号
     
133. 
     
134. %第一步：变量初始化
     
135. [m,n]=size(X);
     
136. Y1p=zeros(m,n);
     
137. Y2p=zeros(m,n);
     
138. Y3p=zeros(m,n);
     
139. 
     
140. %第二步：计算第一道工序的安排
     
141. Q1=zeros(m,1);
     
142. Q2=zeros(m,1);
     
143. R=X(:,1);%取出第一道工序
     
144. Q3=floor(R);%向下取整即得到各工件在第一道工序使用的机器的编号
     
145. %下面计算各工件第一道工序的开始时刻和结束时刻
     
146. for i=1:P(1)%取出机器编号
     
147. pos=find(Q3==i);%取出使用编号为i的机器为其加工的工件的编号
     
148. lenpos=length(pos);
     
149. if lenpos>=1
     
150. Q1(pos(1))=0;
     
151. Q2(pos(1))=T(pos(1),1);
     
152. if lenpos>=2
     
153. for j=2:lenpos
     
154. Q1(pos(j))=Q2(pos(j-1));
     
155. Q2(pos(j))=Q2(pos(j-1))+T(pos(j),1);
     
156. end
     
157. end
     
158. end
     
159. end
     
160. Y1p(:,1)=Q1;
     
161. Y2p(:,1)=Q2;
     
162. Y3p(:,1)=Q3;
     
163. 
     
164. %第三步：计算剩余工序的安排
     
165. for k=2:n
     
166. R=X(:,k);%取出第k道工序
     
167. Q3=floor(R);%向下取整即得到各工件在第k道工序使用的机器的编号
     
168. %下面计算各工件第k道工序的开始时刻和结束时刻
     
169. for i=1:P(k)%取出机器编号
     
170. pos=find(Q3==i);%取出使用编号为i的机器为其加工的工件的编号
     
171. lenpos=length(pos);
     
172. if lenpos>=1
     
173. POS=zeros(1,lenpos);%上一个工序完成时间由早到晚的排序
     
174. for jj=1:lenpos
     
175. MinEndTime=min(EndTime);
     
176. ppp=find(EndTime==MinEndTime);
     
177. POS(jj)=ppp(1);
     
178. EndTime(ppp(1))=Inf;
     
179. end
     
180. %根据上一个工序完成时刻的早晚，计算各工件第k道工序的开始时刻和结束时刻
     
181. if lenpos>=2
     
182. for j=2:lenpos
     
183. Q1(pos(POS(j)))=Y2p(pos(POS(j)),k-1);%预定的开始时刻为上一个工序的结束时刻
     
184. if Q1(pos(POS(j)))
     
185. Q1(pos(POS(j)))=Q2(pos(POS(j-1)));
     
186. end
     
187. end
     
188. end
     
189. end
     
190. end
     
191. Y1p(:,k)=Q1;
     
192. Y2p(:,k)=Q2;
     
193. Y3p(:,k)=Q3;
     
194. end
     
195. 
     
196. %第四步：计算最优的Makespan值
     
197. Y2m=Y2p(:,n);
     
198. Zp=max(Y2m);
     
199. 
     
200. %第五步：绘甘特图
     
201. if plotif
     
202. for i=1:m
     
203. for j=1:n
     
204. mPoint1=Y1p(i,j);
     
205. mPoint2=Y2p(i,j);
     
206. mText=m+1-i;
     
207. PlotRec(mPoint1,mPoint2,mText);
     
208. Word=num2str(Y3p(i,j));
     
209. %text(0.5*mPoint1+0.5*mPoint2,mText-0.5,Word);
     
210. hold on
     
211. x1=mPoint1;y1=mText-1;
     
212. x2=mPoint2;y2=mText-1;
     
213. x3=mPoint2;y3=mText;
     
214. x4=mPoint1;y4=mText;
     
215. %fill([x1,x2,x3,x4],[y1,y2,y3,y4],'r');
     
216. fill([x1,x2,x3,x4],[y1,y2,y3,y4],[1,0.5,1]);
     
217. text(0.5*mPoint1+0.5*mPoint2,mText-0.5,Word);
     
218. end
     
219. end
     
220. end
     
221. 
     
222. function PlotRec(mPoint1,mPoint2,mText)
     
223. % 此函数画出小矩形
     
224. % 输入:
     
225. % mPoint1 输入点1,较小,横坐标
     
226. % mPoint2 输入点2,较大,横坐标
     
227. % mText 输入的文本,序号,纵坐标
     
228. vPoint = zeros(4,2) ;
     
229. vPoint(1,:) = [mPoint1,mText-1];
     
230. vPoint(2,:) = [mPoint2,mText-1];
     
231. vPoint(3,:) = [mPoint1,mText];
     
232. vPoint(4,:) = [mPoint2,mText];
     
233. plot([vPoint(1,1),vPoint(2,1)],[vPoint(1,2),vPoint(2,2)]);
     
234. hold on ;
     
235. plot([vPoint(1,1),vPoint(3,1)],[vPoint(1,2),vPoint(3,2)]);
     
236. plot([vPoint(2,1),vPoint(4,1)],[vPoint(2,2),vPoint(4,2)]);
     
237. plot([vPoint(3,1),vPoint(4,1)],[vPoint(3,2),vPoint(4,2)]);

复制代码

--------------------------------------------------------------------------------

【GreenSim团队简介】

——更多原创程序，请访问GreenSim团队的主页→http://blog.sina.com.cn/greensim

    GreenSim团队长期从事建模仿真、算法设计、代写程序等业务，为您提供实验、课题、论文、项目等方面的仿真服务。
    GreenSim团队擅长的学科领域主要有
【智能算法】遗传算法、模拟退火、神经网络、蚁群算法、免疫优化、支持向量机、拟物拟人
【运筹优化】数学建模、数值优化算法、大规模/非线性问题、组合优化、NP完全问题
【信息与通信】无线通信、信号处理、计算机网络
【其它学科】任何问题，只要能归结为数学问题或者算法仿真，我们都有实力和信心解决
    GreenSim团队的联系方式：
QQ: 330264704，535115369  Email：greensim@163.com（来 信 必 复）

[ 本帖最后由 aiwa 于 2007-3-18 13:38 编辑 ]