function dialsuanfa(T)
%程序说明
clc
disp(‘========================================================================================‘);
disp(‘ 《基于LOGIT的STOCH配流法——dial算法》‘);
disp(‘运行环境:MATLAB 8.3.0.532 ‘);
disp(‘作者信息:兰州交通大学 刘志祥 QQ:531548824‘);
disp(‘说明:本程序用于进行静态配流,缺点是有效路径的界定太过严,请读者根据需要进行修改。‘);
disp(‘ dial算法分四大步骤:一是求最短路,二是求边权似然数,三是求路权,四是配流‘);
disp(‘========================================================================================‘);
disp(‘按任意键继续...‘);
pause;
%数据获取,人机交互
disp(‘请按照提示输入以下参数‘);
Q=input(‘总需求量 :‘);
thita=input(‘参数thita :‘);
r=input(‘起 点 :‘);
s=input(‘终 点 :‘);
n=size(T,1);
%初始化
L=zeros(n,n);
W=zeros(n,n);
X=zeros(n,n);
%求最短距离矩阵
disp(‘第一步:求最短距离,其中‘);
disp(‘------------------------------------------------------------------------------------------‘);
disp(‘ R—起点r到其他点的最短距离‘);
disp(‘ S—其他点到终点s的最短距离‘);
disp(‘按任意键继续...‘);
pause;
Tmin=zuiduanjulijuzhen(T);
R=Tmin(r,:)
S=Tmin(:,s)‘%注意因为方向性,这里作转置处理
%找上游节点和下游节点(显然up和down互为转置——对称,因为若j是i的下游节点,则i必是j的上游节点)
for i=1:n
for j=1:n
if T(i,j)~=0&&T(i,j)~=inf
down(i,j)=1;
up(j,i)=1;
else
down(i,j)=0;
up(j,i)=0;
end
end
end
down;
up;
%计算边权
disp(‘第二步:计算边权似然值(任意键继续)‘);
disp(‘------------------------------------------------------------------------------------------‘);
pause;
for i=1:n
for j=1:n
if down(i,j)~=0
if R(i)<R(j)&&S(i)>S(j)
P=1;
else
P=0;
end
L(i,j)=P*exp(thita*(R(j)-R(i)-T(i,j)));
end
end
end
L
%计算路权
disp(‘第三步:计算路权(任意键继续)‘);
disp(‘------------------------------------------------------------------------------------------‘);
pause;
for i=1:n
for j=1:n
if down(i,j)~=0
if R(i)<R(j)&&S(i)>S(j)
if i==r
W(i,j)=L(i,j);
else
W(i,j)=L(i,j)*(up(i,:)*W(:,i));
%这是核心句,先找到上游节点,再写出上游节点到i的W值(若还没有,则递归直到能够算出),请仔细查阅路权的算法好好理解。
end
end
end
end
end
W
%配流
disp(‘第四步:配流(任意键继续)‘);
disp(‘------------------------------------------------------------------------------------------‘);
pause;
for i=n:-1:1
for j=n:-1:1
if down(i,j)~=0
if R(i)<R(j)&&S(i)>S(j)
if j==s
X(i,j)=Q*W(i,j)/((up(j,:)*W(:,j)));
else
X(i,j)=X(j,:)*down(j,:)‘*W(i,j)/(up(j,:)*W(:,j));
%注意X(j,:)是1*n行向量,down(j,:)是1*n行向量,因此要将down向量转置为1*n的列向量才能相乘。
end
end
end
end
end
X
disp(‘==========================================================================================‘);
disp(‘<程序运行完毕>‘);
例:某路网如图所示,阻抗已标注于边上(这里阻抗用距离值表示),起点为v1,终点为v9,θ=1。请用经典dial算法进行配流。
解:(1)写权值矩阵
quanzhijuzhen=[
0 2 Inf 2 Inf Inf Inf Inf Inf
Inf 0 2 Inf 2 Inf Inf Inf Inf
Inf Inf 0 Inf Inf 2 Inf Inf Inf
Inf Inf Inf 0 1 Inf 2 Inf Inf
Inf Inf Inf Inf 0 1 Inf 2 Inf
Inf Inf Inf Inf Inf 0 Inf Inf 2
Inf Inf Inf Inf Inf Inf 0 2 Inf
Inf Inf Inf Inf Inf Inf Inf 0 2
Inf Inf Inf Inf Inf Inf Inf Inf 0];
(2)带入程序(格式整理后输出如下)
说明:图形表示如下
版权声明:博主文章可以被非商用转载,但请务必注明出处,因水平有限,难免出错,在此免责。