matlab求解vrptw,遗传蚁群模拟退火粒子群算法求带时间窗的物流配送路径VRPTW问题MATLAB源码...

在编程求解VRPTW问题之前，首先要明白并完成VRP问题的编程，也就是说VRPTW问题是在VRP问题的进一步计算，也就是根据VRP路径进一步计算时间惩罚成本。

那么VRP具体是怎么编程的？这个我在上一篇博客里有介绍，如果你不会，先把上面一篇博客内容搞明白。

那么VRP问题搞明白了之后，VRPTW就十分简单了。也就是像计算汽车载货量一样，对配送时间进行计算，配送时间包括车辆在路径的行驶时间，以及车辆在客户点的服务时间。计算好车辆的时间后，然后判断车辆到达每个客户的时间是否落在该客户的时间窗内。具体大家可以看一下程序的计算效果。

(2019年8月1日更新)对原VRPTW程序的输出效果进行更新，主要是路径图，不同颜色区分不同的车辆，如下图，其他部分由些许调整，不过不是很大，效果图如下

2019年8月1日之前内容：

如下

%%%% 遗传算法求解带时间窗的VRP问题,VRPTW

%%%% by 圆 一个有心的人在用心做事

%%%% 2018-05-14

%%%% 如有问题请联系

QQ530807088 新浪微博:MATLAB圆创工作室

% ======================= 高 端 定 制 · 华

丽 分 割 ======================

clc

clear

close all

tic

data =

xlsread('数据.xlsx','C3:I12');

JW = data(:,1:2);

dem = data(:,3); % 各网点的需求量

demand =

rem(dem,6); %

各网点除去整车运输的量(整车运输是指网点需求量多于汽车装载量时，需要在该网点多次运输)

ntimes = (dem - demand)/6; %

各网点整车运输的次数

%% 经纬度转大地坐标，单位km

e2=0.0066943799;

N=6386942;

H=44.4; % 海拔

W=JW(:,2)/180*pi; % 维度

J=JW(:,1)/180*pi; % 经度

X=-(N+H)*cos(W).*cos(J)/1000;

Y=(N+H)*cos(W).*sin(J)/1000;

position=[X Y];

%%

CityNum=9; % 需求节点数

CarDistance=300000;% 车最大行驶距离

CarLoad=6; % 车容量

speed=30; %

车辆行驶速度

c0=0; % 发车成本

c1=3; % 车辆行驶单位距离成本

CarNum=18; % 配送车辆总数，实际车辆可以少于该值

fitmax=10000000;% 最大惩罚

%%

NP=80; % 种群个体数量

maxgen=200;

Pc=0.9; % 交叉概率

Pm=0.1; % 变异概率

Gap=0.9; %

代沟(Generation gap)

%%

ET=data(:,6); %

最早服务时间窗

LT=data(:,7); %

最晚服务时间窗

ST=zeros(size(ET));% 服务时间

CE=[0

ones(1,CityNum)]*1000; % 早到惩罚系数

CL=[0

ones(1,CityNum)]*1000; % 晚到惩罚系数

%% 计算各城市之间的距离

distance=zeros(CityNum+1);

for i=1:CityNum+1

for j=i+1:CityNum+1

distance(i,j)=((position(i,1)-position(j,1))^2+(position(i,2)-position(j,2))^2)^0.5;

distance(j,i)=distance(i,j);

end

end

%% 整车运输费用

for j = 2:length(dem)

mT = distance(1,j)/speed;

if mT < ET(j)

dpunish =

(ET(j)-mT)*CE(j);

elseif mT > LT(j)

dpunish =

(mT-LT(j))*CL(j);

else

dpunish = 0;

end

mCost(j) = ntimes(j)*(distance(1,j)*c1*2);

end

dCost = sum(mCost);

%% 路径初始化

X=initpop(NP,CityNum,CarNum);

Xa=X(1,:);

%% 迭代

gen=1;

while gen<=maxgen

gen

% 计算适应值矩阵

[allcost,fit]=fitness(distance,demand,X,ET,LT,CE,CL,ST,CarDistance,CarLoad,speed,fitmax,c0,c1);

%找出最优个体适应值

allcost=allcost+dCost;

[leastcost,bestindex]=min(allcost);

bestindex=bestindex(1);

fpbest(gen)=leastcost; % 最小适应值fit的集

pbest(gen,:)=X(bestindex,:);% 最优个体集

% 选择

XSel=Select(X,fit,Gap);

% 交叉操作

XSel=Cross(XSel,Pc);

% 变异

XSel=Mutate(XSel,Pm);

% 逆转操作

%XSel=Reverse(distance,demand,XSel,ET,LT,CE,CL,ST,CarDistance,CarLoad,speed,fitmax,c0,c1);

% 重插入子代的新种群

X=Reins(X,XSel,fit);

gen=gen+1;

end

% 找出最优的适应值、个体

[minfpbest,minindex]=min(fpbest);%

取最优适应值的位置、最优适应值

minindex=minindex(1);

bestRoute=pbest(minindex,:); %

取最优个体

Path=bestRoute;

for i=1:length(Path)-1

if Path(i)-Path(i+1)==0

Path(i)=0;

end

end

ii=find(Path==0);

Path(ii)=[];

for j=1:length(Path) %

编码各减1，与文中的编码一致

Path(j)= Path(j)-1;

end

fpbest; %查看适应值的变化情况

%% 迭代图

figure

plot(fpbest)

title('遗传算法优化过程')

xlabel('迭代次数')

ylabel('最优值')

%% 计算结果数据输出

clc

toc

disp(['最优成本为：',num2str(minfpbest)]);

OutputResult(distance,demand,Path,ET,LT,CE,CL,ST,CarDistance,CarLoad,speed,fitmax,c0,c1)

disp('各网点满载运输次数')

ntimes'

disp('各网点满载运输费用')

mCost

disp('----------------------------------------------------------------------------------------------------------------------------')

%% 绘制实际路线，如果没有各点坐标时，可以不用绘制

str={'冷链仓库','虹口','嘉定','普陀','浦东','闵行','青浦','松江','金山','嘉兴'};

figure

for i=1:length(Path)-1

plot(JW(Path(i:i+1)+1,1),JW(Path(i:i+1)+1,2),'ro-','MarkerFaceColor','b')

hold on

end

plot(JW(1,1),JW(1,2),'bp','MarkerFaceColor','r','MarkerSize',15)

for i=1:CityNum+1

text(JW(i,1)+0.01,JW(i,2),str{i})

end

title('麦德龙冷链配送路线图')

xlabel('经度')

ylabel('纬度')

所有程序文件如下图