本文提出了一种结合粒子群优化的遗传算法，用于解决多式联运路径规划问题，旨在降低物流成本、提高运输效率。通过对军事运输中多式联运的应用分析，建立了以时间最短、路线最短、成本最低为目标的模型，并通过Matlab代码实现优化过程。仿真结果显示了算法的有效性。
1 简介
物流运输方式由公路、水路、空运及管道等 3 种方式组成，3 种运输方式在技术上、经济上各有长短，都有适宜的 使用范围，每种运输方式单独运用很难实现节约资源、降本增效。随着我国经济不断发展以及布局网络技术的不断深化，多式 联运通过把传统的、单一的运输方式进行择优组合，充分利用了各个运输方式现有的设施设备，实现了运输过程中的资源整 合，有利于运输过程中的可持续发展及达成规模经济中降本增效的目的，同时提高了物流行业竞争力。特别是通过公铁水多式 联运路径优化，构建以运输时间最少、运输线路距离最短、运输成本最低的公铁水多式联运模式，对于物流企业节约资源、降本增效意义重大。本文提出了一个以粒子群结合遗传算法为主框架的解决方案,用来求解多式联运的路径规划问题。本文从运输需求内容、运输过程、应用场景等角度对多式联运在军事运输中的应用进行分析,定义多式联运路径规划问题,建立分别以时问最短、路线最短、成本最低为目标的多式联运路径规划模型。
2 部分代码
clear
clc
close all
tic
%% 用importdata这个函数来读取文件
shuju= xlsread('shuju.xlsx', 'Sheet1');
bl=0;
cap=100; %车辆最大装载量
%% 提取数据信息
zuobiao=shuju(:,2:3); %所有点的坐标x和y
customer=zuobiao(2:end,:); %顾客坐标
cusnum=size(customer,1); %顾客数
v_num=8; %车辆最多使用数目
demands=shuju(2:end,4); %需求量
a=shuju(2:end,5); %顾客时间窗开始时间[a[i],b[i]]
b=shuju(2:end,6); %顾客时间窗结束时间[a[i],b[i]]
s=shuju(2:end,7); %客户点的服务时间
%% 距离矩阵
h=pdist(zuobiao);
lldist=squareform(h).*1.5; %路路距离矩阵
htdist=squareform(h); %飞机距离矩阵
hydist=squareform(h).*1.2; %%火车距离矩阵
%% 遗传算法参数设置
alpha=100000; %违反的容量约束的惩罚函数系数
belta=0.5;%违反时间窗约束的惩罚函数系数
belta2=0.5;
chesu=[1,5,2];
NIND=300; %种群大小
MAXGEN=500; %迭代次数
Pc=0.9; %交叉概率
Pm=0.05; %变异概率
GGAP=0.9; %代沟(Generation gap)
N=cusnum+v_num-1; %染色体长度=顾客数目+车辆最多使用数目-1
%% 粒子群参数
lx=3;
w=1; %惯性因子
wdamp=0.99; %惯性因子衰减率
c1=1.5; %个体学习因子
c2=2.0; %全局学习因子
XvMin=1; %Xv下限
XvMax=lx; %Xv上限
VvMin=-(lx-1); %Vv下限
VvMax=lx-1; %Vv上限
%% 初始化种群
Chrom=InitPopCW(NIND,N,cusnum,a,demands,cap,XvMin,XvMax,VvMin,VvMax); %构造初始解
%% 输出随机解的路线和总距离
disp('初始种群中的一个随机值:')
[VC,NV,PD]=decode(Chrom{1},cusnum);
% disp(['总距离：',num2str(TD)]);
disp(['车辆使用数目：',num2str(NV)]);
disp('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~')
%% 优化
gen=1;
figure;
hold on;box on
xlim([0,MAXGEN])
title('优化过程')
xlabel('代数')
ylabel('最优值')
ObjV=calObj(Chrom,cusnum,cap,demands,a,b,s,lldist,htdist,hydist,alpha,belta,belta2,chesu,bl); %计算种群目标函数值
preObjV=min(ObjV);
[~,minInd]=min(ObjV);
gbest_pos=Chrom{minInd,1}(2,:); %假设第一个粒子位置为全局最优位置
gbest_obj=preObjV; %第一个粒子位置的目标函数值
pbest_pos=cell(NIND,1); %初始化各个粒子的个体最优位置
pbest_obj=ObjV; %初始化各个粒子的个体最优的目标函数值
for i=1:NIND
particle=Chrom{i,1}; %第i个粒子
position=particle(2,:); %第i个粒子的位置
pbest_pos{i,1}=position; %初始化这个粒子的个体最优
if pbest_obj(i,1)<gbest_obj
%更新初始种群中的全局最优粒子
gbest_obj=pbest_obj(i,1);
gbest_pos=position;
end
end
%%
while gen<=MAXGEN
%% 计算适应度
ObjV=calObj(Chrom,cusnum,cap,demands,a,b,s,lldist,htdist,hydist,alpha,belta,belta2,chesu,bl); %计算种群目标函数值
line([gen-1,gen],[preObjV,min(ObjV)]);pause(0.0001)%画图 最优函数
preObjV=min(ObjV);
FitnV=Fitness(ObjV);
%% 选择
SelCh=Select(Chrom,FitnV,GGAP);
%% OX交叉操作
SelCh=Recombin(SelCh,Pc);
%% 变异
SelCh=Mutate(SelCh,Pm);
%% 重插入子代的新种群
Chrom=Reins(Chrom,SelCh,ObjV);
ObjV=calObj(Chrom,cusnum,cap,demands,a,b,s,lldist,htdist,hydist,alpha,belta,belta2,chesu,bl);
%% 更新各个粒子的位置和速度
Chrom1=Chrom;
for i=1:NIND
particle=Chrom{i,1}; %第i个粒子
position=particle(2,:); %第i个粒子的位置
Xv=position(1,:);
velocity=particle(3,:); %第i个粒子的速度
Vv=velocity(1,:);
%% 更新速度
velocity=w*velocity+ +c1*rand([1,N]).*(pbest_pos{i,1}-position)...
+c2*rand([1,N]).*(gbest_pos-position);
%% 速度越界处理
velocity(1,:)=max(velocity(1,:),VvMin);
velocity(1,:)=min(velocity(1,:),VvMax);
%% 更新位置
position=position+velocity;
position(1,:)=ceil(position(1,:)); %对Xv向上取整
%% 速度镜像影响
IsOutside=(position(1,:)<XvMin | position(1,:)>XvMax );
velocity(IsOutside)=-velocity(IsOutside);
%% 位置越界处理
position(1,:)=max(position(1,:),XvMin);
position(1,:)=min(position(1,:),XvMax);
Chrom1{i,1}(2,:)=position;
end
ObjV1=calObj(Chrom1,cusnum,cap,demands,a,b,s,lldist,htdist,hydist,alpha,belta,belta2,chesu,bl);
%% 重插入子代的新种群
Chrom=Copy_of_Reins(Chrom,Chrom1,ObjV1,NIND,ObjV);
%% 打印当前最优解
ObjV=calObj(Chrom,cusnum,cap,demands,a,b,s,lldist,htdist,hydist,alpha,belta,belta2,chesu,bl); %计算种群目标函数值
[~,minInd]=min(ObjV);
for i=1:NIND
if min(ObjV)<ObjV(i)
particle=Chrom{minInd,1};
position=particle(2,:);
pbest_pos{i,1}=position;
pbest_obj(i,1)=ObjV(i);
%% 更新全局最优
if preObjV<gbest_obj
gbest_pos=pbest_pos{i,1};
gbest_obj=pbest_obj(i,1);
end
end
end
disp(['第',num2str(gen),'代最优解:'])
fprintf('\n')
%% 更新迭代次数
gen=gen+1 ;
end
%% 画出最优解的路线图
ObjV=calObj(Chrom,cusnum,cap,demands,a,b,s,lldist,htdist,hydist,alpha,belta,belta2,chesu,bl); %计算种群目标函数值
[minObjV,minInd]=min(ObjV);
%% 输出最优解的路线和总距离
disp('最优解:')
bestChrom=Chrom{minInd(1)};
[VC,NV,PD]=decode(bestChrom,cusnum);
disp('-------------------------------------------------------------')
[cost]=costFuction(VC,NV,PD,a,b,s,lldist,htdist,hydist,demands,cap,alpha,belta,belta2,chesu,bl);
disp(['总成本：',num2str(cost)]);
disp(['路路：',num2str(1)]);
disp(['航天：',num2str(2)]);
disp(['海运：',num2str(3)]);
%% 画出最终路线图
draw_Best(VC,zuobiao,PD);
% toc
function [A,c]=cdzwz(dist,route,dis)
dis=0.9*dis;
n=length(route);
p_l=0;
DL=1;
c=[];
A=[];
for i=1:n
if i==1
p_l=p_l+dist(1,route(i)+1);
DL=DL-(dist(1,route(i)+1)/dis);
if DL<0
c=[c,i];
DL=1;
end
else
p_l=p_l+dist(route(i-1)+1,route(i)+1);
DL=DL-dist(route(i-1)+1,route(i)+1)/dis;
if DL<0
[a,b]=min(dist(route(i-1)+1,end-9:end));
DL=DL+dist(route(i-1)+1,route(i)+1)/dis-a/dis;
if DL<0
c=[c,i-1];
DL=1;
[a,b]=min(dist(route(i-1)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i-1)+1,b+41)/dis;
else
c=[c,i];
DL=1;
[a,b]=min(dist(route(i)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i)+1,b+41)/dis;
end
end
end
end
p_l=p_l+dist(route(end)+1,1);
DL=DL-dist(route(end)+1,1)/dis;
if DL<0
c=[c,n];
DL=1;
[a,b]=min(dist(route(i)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i)+1,b+41)/dis;
end
end
function [A,c]=cdzwz(dist,route,dis)
dis=0.9*dis;
n=length(route);
p_l=0;
DL=1;
c=[];
A=[];
for i=1:n
if i==1
p_l=p_l+dist(1,route(i)+1);
DL=DL-(dist(1,route(i)+1)/dis);
if DL<0
c=[c,i];
DL=1;
end
else
p_l=p_l+dist(route(i-1)+1,route(i)+1);
DL=DL-dist(route(i-1)+1,route(i)+1)/dis;
if DL<0
[a,b]=min(dist(route(i-1)+1,end-9:end));
DL=DL+dist(route(i-1)+1,route(i)+1)/dis-a/dis;
if DL<0
c=[c,i-1];
DL=1;
[a,b]=min(dist(route(i-1)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i-1)+1,b+41)/dis;
else
c=[c,i];
DL=1;
[a,b]=min(dist(route(i)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i)+1,b+41)/dis;
end
end
end
end
p_l=p_l+dist(route(end)+1,1);
DL=DL-dist(route(end)+1,1)/dis;
if DL<0
c=[c,n];
DL=1;
[a,b]=min(dist(route(i)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i)+1,b+41)/dis;
end
end
function [A,c]=cdzwz(dist,route,dis)
dis=0.9*dis;
n=length(route);
p_l=0;
DL=1;
c=[];
A=[];
for i=1:n
if i==1
p_l=p_l+dist(1,route(i)+1);
DL=DL-(dist(1,route(i)+1)/dis);
if DL<0
c=[c,i];
DL=1;
end
else
p_l=p_l+dist(route(i-1)+1,route(i)+1);
DL=DL-dist(route(i-1)+1,route(i)+1)/dis;
if DL<0
[a,b]=min(dist(route(i-1)+1,end-9:end));
DL=DL+dist(route(i-1)+1,route(i)+1)/dis-a/dis;
if DL<0
c=[c,i-1];
DL=1;
[a,b]=min(dist(route(i-1)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i-1)+1,b+41)/dis;
else
c=[c,i];
DL=1;
[a,b]=min(dist(route(i)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i)+1,b+41)/dis;
end
end
end
end
p_l=p_l+dist(route(end)+1,1);
DL=DL-dist(route(end)+1,1)/dis;
if DL<0
c=[c,n];
DL=1;
[a,b]=min(dist(route(i)+1,end-9:end));
A=[A,b+40];
DL=DL-dist(route(i)+1,b+41)/dis;
end
end
%% 根据vehicles_customer整理出final_vehicles_customer，将vehicles_customer中空的数组移除
function [ final_VC,VU ] = deal_vehicles_customer( VC )
vecnum=size(VC,1); %车辆数
final_VC={}; %整理后的vehicles_customer
count=1; %计数器
for i=1:vecnum
par_seq=VC{i}; %每辆车所经过的顾客
%如果该辆车所经过顾客的数量不为0，则将其所经过的顾客数组添加到final_vehicles_customer中
if ~isempty(par_seq)
final_VC{count}=par_seq;
count=count+1;
end
end
%% 为了容易看，将上述生成的1行多列的final_vehicles_customer转置了，变成多行1列的了
final_VC=final_VC';
VU=size(final_VC,1); %所使用的车辆数
end
%% 根据vehicles_customer整理出final_vehicles_customer，将vehicles_customer中空的数组移除
function [ final_VC,VU ] = deal_vehicles_customer( VC )
vecnum=size(VC,1); %车辆数
final_VC={}; %整理后的vehicles_customer
count=1; %计数器
for i=1:vecnum
par_seq=VC{i}; %每辆车所经过的顾客
%如果该辆车所经过顾客的数量不为0，则将其所经过的顾客数组添加到final_vehicles_customer中
if ~isempty(par_seq)
final_VC{count}=par_seq;
count=count+1;
end
end
%% 为了容易看，将上述生成的1行多列的final_vehicles_customer转置了，变成多行1列的了
final_VC=final_VC';
VU=size(final_VC,1); %所使用的车辆数
end
3 仿真结果
4 参考文献
[1]丁建伟. 基于遗传算法的多式联运应急管理研究[J]. 电脑开发与应用, 2009, 22(1):3.
博主简介：擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真，相关matlab代码问题可私信交流。
部分理论引用网络文献，若有侵权联系博主删除。
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
