综合算法03—FrankWolfe_BPR配流算法
2024-05-18 21:16:10
%% 算法符号及程序说明
%说明：本程序为采用美国联邦公路阻抗函数BPR时的frankwolfe算法，考虑了换乘（已经将等待时
%间考虑在内并在K短路的确定过程中计算）及拥挤附加时间，在路网情况已知时配流并求出费用值。
%计算条件：根据K短路算法计算出考虑换乘的路径-路段关系矩阵和边权。
%缺点:BPR函数或许不能准确的描述轨道交通系统的路段阻抗，列车的能力也不能用程序中那么简单
%代替（但是，这无关紧要，暂时可以认为是正确的，后续可以用已知的路段函元胞和能力矩阵替换）
%下一步研究方向：网络配流--所有OD对之间的配流。
%==============================================================================
%符号体系：
% Q-OD需求量
% Cs-路段额定能力（列车座位数）
%Cmax-路段极限能力（列车座位数+人均面积达标的时的站立人数）
% t-路段阻抗函数（值）
%numf-起点到终点的路径数量（自己按照K短路算法进行定义）
%numx-网络中的路段数量（直接连接两节点的边数，包含所有路段）
% W-边权（0流量时的阻抗）,在该程序中，通过调用K短路换乘算法自动标号索引得出。
%Transfer_time-某路径的路段换乘时间
% G-路段拥挤附加时间
% A-一般拥挤放大系数(介于0到等车时间之间），可以看作一个综合系数，指部分人选择挤车，部分人选择下趟列车。
% B-特别拥挤放大系数（等车时间，与发车频率对应），如果过度拥挤，则无法乘车的人必须等待后车到来，等1趟时间加B，等2趟，时间加2B，等几趟由无法乘车人数与极限能力决定。
% L-路网路线矩阵（0-1矩阵，表示路线与路段的关系，1为包含，0为不包含）
% Mxf-路径与路段的关系矩阵（每一行表示一个路径，矩阵值为1表示经过该路段，否则不经过该路段）
%cont-当前配流次数
%Ckrs-路径配流结果（更新后的阻抗值）
% X1-当前配流结果
% Y1-辅助流量
% S-搜索方向
%lmbda-搜索步长
% X2-新的迭代点
% Z-目标函数（通用配流目标函数，可认为是费用）
%==============================================================================
clear all
format
clc
warning off
disp('========================================================================');
disp(' 《FrankWolfe_BPR配流算法》');
disp('运行环境：MATLAB 8.3.0.532 ');
disp('制 作 人：兰州交通大学 刘志祥');
disp('完成时间：2015-12-30');
disp('Q Q:531548824');
disp('请提出宝贵的意见!');
disp('=========================================================================');disp('----------------------------------------------------------------------')
fprintf('问题描述：已知路网数据（Road_Net）和线路数据（Line_Station）及OD需求（Q），\n求指定OD对的配流结果.\n\n');
disp('----------------------------------------------------------------------')%% 定义变量及常数
syms lambda real
for i=1:100syms x(i) real;
end
cont=0;
e=inf;
A=2;
B=5;
cs=1460;
cmax=2070;
Q=7500; %% 输入路网数据
%例1
% Road_Net=[0 15 5 inf;15 0 inf 20;5 inf 0 25;inf 20 25 0]; %路网关联矩阵，表示点到点的距离（或时间、费用）
%例2
Road_Net =[0 2 Inf 1 Inf Inf Inf Inf Inf2 0 3 Inf 3 Inf Inf Inf InfInf 3 0 Inf Inf 2 Inf Inf Inf1 Inf Inf 0 3 Inf 5 Inf InfInf 3 Inf 3 0 1 Inf 2 InfInf Inf 2 Inf 1 0 Inf Inf 3Inf Inf Inf 5 Inf Inf 0 2 InfInf Inf Inf Inf 2 Inf 2 0 4Inf Inf Inf Inf Inf 3 Inf 4 0];
Line_Station={[1 2 3 6 9 8 7 4],[2 5 8 9],[1 4 5 6]};%% 调用KSP（K shortest path)计算出W,L,Mxf
[W_Section,Line_Section_01Mat,Mxf]=KSP(Road_Net,Line_Station);
W=cell2mat(W_Section)
L=Line_Section_01Mat
Mxf=Mxf;
if isempty(Mxf)Mxf=zeros(1,length(W));disp('警告：无路径可达.');
end
disp('路径：（矩阵Mxf的行表示路径号，列表示路段号）');Mxf%% 计算路网基础数据及路网表达
numf=size(Mxf,1);
numx=length(W);
Cs=cs.*ones(1,numx);
Cmax=cmax.*ones(1,numx);
x=x(1:numx);
disp('路网介绍')
disp(['总需求：',num2str(Q)]);
disp(['路径数：',num2str(numf)]);
disp(['路段数：',num2str(numx)]);
disp(['初始路段阻抗：',num2str(W)]);
disp(['路段额定能力：',num2str(Cs)]);
disp(['路段极限能力：',num2str(Cmax)]);
delete('E:\MATLAB\自己的算法\FW_BPR算法计算过程.txt');
fid=fopen('E:\MATLAB\自己的算法\FW_BPR算法计算过程.txt','a');%% step1_初始化
X0=zeros(1,numx);
t=zeros(1,numx);
t=W.*(1+0.15*(x./Cmax).^4); %路段走行时间函数
tt=t;
disp('路段阻抗函数t=');
disp(vpa(t.',2));
t=W.*(1+0.15*(X0./Cmax).^4); %路段的初始阻抗函数
Ckrs=(Mxf*t')'; %路径的走行时间初值
disp(['初始路径阻抗：',num2str(Ckrs)]);
%% step2_更新流量和阻抗
[Min,index]=min(Ckrs);
X1=Mxf(index,:).*Q; %全有全无法为最短路径上的路段分配流量fprintf(fid,'\r\n==============================================================\r\n');
while e>1e-3cont=cont+1;%% 计算路段附加拥挤时间GG=zeros(1,length(W));a1=(X1<=Cs); %a1-指向流量不大于额定能力的下标a2=(X1>Cs&X1<=Cmax); %a2-指向流量介于额定能力和极限能力的下标a3=(X1>Cmax); %a3-指向流量大于极限能力的下标G(a1)=0;G(a2)=(X1(a2)-Cs(a2))./Cs(a2).*A;G(a3)=(Cmax(a3)-Cs(a3))./Cs(a3).*A+(X1(a3)-Cmax(a3))./Cs(a3).*B;%% 计算路段阻抗函数并求出路径阻抗输出至指定文件t=(W+G).*(1+0.15*(X1./Cmax).^4); %路段时间fprintf(fid,' <第%d次的计算结果>\r\n',cont);Ckrs=(Mxf*t')'; %路径时间fprintf(fid,'阻 抗 值：');fprintf(fid,'%9.2f ',Ckrs);fprintf(fid,'\r\n');%% 输出当前流量至指定文件fprintf(fid,'当前流量：');fprintf(fid,'%9.2f ',X1);fprintf(fid,'\r\n');%% 输出输出拥挤附加时间至指定文件fprintf(fid,'拥挤附加：');fprintf(fid,'%9.2f ',G);fprintf(fid,'\r\n');%% step3_计算辅助流量并输出至指定文件[Min,index]=min(Ckrs);Y1=Mxf(index,:).*Q; %全有全无法求辅助流量fprintf(fid,'辅助流量：');fprintf(fid,'%9.2f ',Y1);fprintf(fid,'\r\n');%% step4_求搜索方向和步长并输出至指定文件S=Y1-X1; %搜索方向fprintf(fid,'搜索方向：');fprintf(fid,'%9.2f ',S);fprintf(fid,'\r\n'); X2=X1+lambda*S; %先将X2用X1和lambda进行表达t=(W+G).*(1+0.15*(X2./Cmax).^4); %含lambda的阻抗表达f=sum(S.*t,2); %2表示按行求和 if S==0lambda2=0;else lambda1=double(solve(f)); %求解方程，确定步长。k=length(lambda1); %如步长lambda1的解不唯一，取实数，且大于0小于1；for m=1:kif lambda1(m,1)>=0&&lambda1(m,1)<=1lambda2=lambda1(m,1);endendendfprintf(fid,'迭代步长：');fprintf(fid,' lambda=%9.2f ',lambda2);fprintf(fid,'\r\n');%% step5_确定新的迭代点并输出至指定文件X2=X1+lambda2*S; %得到下一步的流量值，且进行下一次迭代fprintf(fid,'新迭代点：');fprintf(fid,'%9.2f ',X2);fprintf(fid,'\r\n'); %% step6_收敛性判断e=sqrt(sum((X2-X1).^2))/sum(X1);X1=X2;fprintf(fid,'\r\n==============================================================\r\n');
end%% 求目标函数值
Xx=zeros(numx,1); %积分下界
Xn=X1; %积分上界
Z=zeros(numx,1);
for i=1:numxZ(i)=int(tt(i),Xx(i),Xn(i)); %对每一个路径积分
end
Z=sum(Z);%% step7_输出结果
disp([' 迭代次数：',num2str(cont)]);
disp([' 误 差 值：',num2str(e)]);
disp([' 配流结果：',num2str(Xn)]);
disp([' 路径阻抗：',num2str(Ckrs)]);
disp([' 目 标 值：',num2str(Z)]);
toc
disp('**************************************************************************************')
disp('提示：详细计算过程请打开 “E:\MATLAB\自己的算法\FW_BPR算法计算过程.txt” 查看。');
open('E:\MATLAB\自己的算法\FW_BPR算法计算过程.txt')
status=fclose('all');
算例分析：
路网如图，节点1到节点9的OD需求为7500，其他数据及说明见程序。求配流结果。
运行结果：
========================================================================《FrankWolfe_BPR配流算法》
运行环境：MATLAB 8.3.0.532
制 作 人：兰州交通大学 刘志祥
完成时间：2015-12-30
Q Q:531548824
请提出宝贵的意见!
=========================================================================
----------------------------------------------------------------------
问题描述：已知路网数据（Road_Net）和线路数据（Line_Station）及OD需求（Q），
求指定OD对的配流结果.----------------------------------------------------------------------
Original=1
Destination=9
起点和终点间有直达线路:1
提示：若选择换乘,可能的换乘站有:2 4 5 6 8
提示：满足距离要求的路径最多只有10条，其中直达线路2条，1次换乘可达的3条，2次换乘可达的2条，3次及以上换乘可达的3条(已舍去）。KSP = Road_Net: {[9x9 double]}Line_Station: {[1 2 3 6 9 8 7 4] [2 5 8 9] [1 4 5 6]}Section_Station: {1x24 cell}W_Section: {[2] [1] [2] [3] [3] [3] [2] [1] [3] [5] [3] [3] [1] [2] [2] [1] [3] [5] [2] [2] [2] [4] [3] [4]}Line_Section: {[1 4 7 17 24 21 18] [5 14 22] [2 9 13]}Line_Section_01Mat: [3x24 double]W_Line_Section: {{1x7 cell} {1x3 cell} {1x3 cell}}Original: {[1]}Destination: {[9]}Kmax: {[100]}Cost_of_Transfer: {[2.5000]}H: {[1.5000]}Transfer_Station: {[2 4 5 6 8]}Irrespective_of_the_Transfer_KPaths: {1x10 cell}Irrespective_of_the_Transfer_KCosts: {[8] [9] [10] [10] [11] [12] [14] [14] [15] [19]}Transfer_Station_of_Path: {[6] [2 5 6] [5] [] [2] [] [6 5] [5 6 8] [5 2] [2 4 5]}Times_of_Transfer: {[1] [3] [1] [0] [1] [0] [2] [3] [2] [3]}KCosts_of_Transfer: {[2.5000] [Inf] [2.5000] [0] [2.5000] [0] [5] [Inf] [5] [Inf]}Consider_the_Transfer_KCosts: {[10.5000] [Inf] [12.5000] [10] [13.5000] [12] [19] [Inf] [20] [Inf]}Consider_not_consider_KPaths_number: {[4 1 6 3 5]}Consider_not_consider_KPaths_Costs: {[10 10.5000 12 12.5000 13.5000]}Consider_the_Transfer_KPaths: {[1 2 3 6 9] [1 4 5 6 9] [1 4 7 8 9] [1 4 5 8 9] [1 2 5 8 9]}Mxf: {[5x24 double]}计算耗时:
时间已过 0.379715 秒。W =2 1 2 3 3 3 2 1 3 5 3 3 1 2 2 1 3 5 2 2 2 4 3 4L =1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 1 1 0 0 1 0 0 10 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 00 1 0 0 0 0 0 0 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0路径：（矩阵Mxf的行表示路径号，列表示路段号）Mxf =1 0 0 1 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 00 1 0 0 0 0 0 0 1 0 0 0 1 0 0 0 1 0 0 0 0 0 0 00 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 1 0 00 1 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 1 0 01 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0路网介绍
总需求：7500
路径数：5
路段数：24
初始路段阻抗：2 1 2 3 3 3 2 1 3 5 3 3 1 2 2 1 3 5 2 2 2 4 3 4
路段额定能力：1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460 1460
路段极限能力：2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070 2070
路段阻抗函数t=1.6e-14*x(1)^4 + 2.08.2e-15*x(2)^4 + 1.01.6e-14*x(3)^4 + 2.02.5e-14*x(4)^4 + 3.02.5e-14*x(5)^4 + 3.02.5e-14*x(6)^4 + 3.01.6e-14*x(7)^4 + 2.08.2e-15*x(8)^4 + 1.02.5e-14*x(9)^4 + 3.04.1e-14*x(10)^4 + 5.02.5e-14*x(11)^4 + 3.02.5e-14*x(12)^4 + 3.08.2e-15*x(13)^4 + 1.01.6e-14*x(14)^4 + 2.01.6e-14*x(15)^4 + 2.08.2e-15*x(16)^4 + 1.02.5e-14*x(17)^4 + 3.04.1e-14*x(18)^4 + 5.01.6e-14*x(19)^4 + 2.01.6e-14*x(20)^4 + 2.01.6e-14*x(21)^4 + 2.03.3e-14*x(22)^4 + 4.02.5e-14*x(23)^4 + 3.03.3e-14*x(24)^4 + 4.0初始路径阻抗：10 8 12 10 11迭代次数：117误 差 值：0.00092332配流结果：3699.3642 3800.6358 0 1901.3204 1798.0438 0 1901.3204 0 2226.3406 1574.2952 0 0 1914.7689 2109.6155 0 0 3816.0893 0 1574.2952 0 0 3683.9107 0 0路径阻抗：54.494 54.8772 55.0724 56.1326 54.5126目 标 值：88930.865
时间已过 111.671728 秒。
**************************************************************************************
提示：详细计算过程请打开 “E:\MATLAB\自己的算法\FW_BPR算法计算过程.txt” 查看。
>> 详细计算过程可去http://download.csdn.net/detail/lzx19901012/9383749查看。
综合算法03—FrankWolfe_BPR配流算法相关推荐
- 配流05—增量配流算法
说明:指定两点间的客流需求总量,建立费用函数(阻抗函数),一般情况下费用是流量的函数,就可以运用增量配流法配流,核心还是全有全无算法,只是该方法把流量等分为N份,每次全有全无配流1份,直至流量全部被加 ...
- 微服务限流及熔断一：四种限流算法（计数器算法、滑动窗口算法、令牌限流算法、漏桶限流算法）
引言 本篇内容根据<spring cloud alibaba 微服务原理与实战>中内容摘取,希望和大家分享限流的思想,本篇不涉及代码层面的实现. 限流的目的 目的:通过限制并发访问数或者限 ...
- 配流03—全有全无配流算法（1）
说明:指定两点间的客流需求总量,建立费用函数(阻抗函数),一般情况下费用是流量的函数,就可以运用全有全无算法进行配流. step1:建立费用函数m文件:feiyonghanshu.m function ...
- 算法9-5：最大流算法的Java代码
残留网络 在介绍最大流算法之前先介绍一下什么是残留网络.残余网络的概念有点类似于集合中的补集概念. 下图是残余网络的样例. 上面的网络是原始网络.以下的网络是计算出的残留网络.残留网络的作用就是用来描 ...
- 接口限流算法：漏桶算法令牌桶算法
工作中对外提供的API 接口设计都要考虑限流,如果不考虑限流,会成系统的连锁反应,轻者响应缓慢,重者系统宕机,整个业务线崩溃,如何应对这种情况呢,我们可以对请求进行引流或者直接拒绝等操作,保持系统的可 ...
- 面试必备：4种经典限流算法讲解
最近,我们的业务系统引入了Guava的RateLimiter限流组件,它是基于令牌桶算法实现的,而令牌桶是非常经典的限流算法.本文将跟大家一起学习几种经典的限流算法. 公众号:捡田螺的小男孩 限流是什 ...
- 接口限流算法：漏桶算法amp;令牌桶算法
转载自 接口限流算法:漏桶算法&令牌桶算法 背景 每一个对外提供的API接口都是需要做流量控制的,不然会导致系统直接崩溃.很简单的例子,和保险丝的原理一样,如果用电符合超载就会烧断保险丝断掉电 ...
- 详解4种经典的限流算法
最近,我们的业务系统引入了Guava的RateLimiter限流组件,它是基于令牌桶算法实现的,而令牌桶是非常经典的限流算法.本文将跟大家一起学习几种经典的限流算法. 限流是什么? 维基百科的概念如下 ...
- 面试必备：四种经典限流算法讲解
大家好,我是田螺. 最近一位朋友去拼夕夕面试,被问了这么一道题:限流算法有哪些?用代码实现令牌桶算法.跟星球好友讨论了一波,发现大家都忘记得差不多了.所以田螺哥再整理一波,常见的四种限流算法,以及简单 ...
最新文章
- Spring Cloud Alibaba - 24 Gateway-路由、断言(Predicate)、过滤器(Filter)初体验
- 框架：Spring Aop、拦截器、过滤器的区别
- 软件项目管理0724：见供应商的体会
- python如何判断字典中是否存在某个键_总结：11个Python3字典内置方法大全及示例...
- 用树莓派做蜘蛛机器人，还是3D打印的！
- 计算机考研备考指南,计算机专业考研复习指南篇
- spring 事务持久性_项目学生：Spring数据的持久性
- office高级应用与python综合案例教程_office高级应用与python综合案例实验指导--详细介绍...
- 增大mysql修改表空间_innodb系统表空间维护方法
- qt4 mysql_qt4连接mysql_MySQL
- CDays-3 习题二 （字典及文件读取练习）及相关内容解析。Python 基础教程
- Java中的properties文件中的key不能使用项目中的接口名和Java文件名
- android 辐射动画_Android 四种动画效果的调用实现代码
- JAVAFX 第三方库 布局 小工具 美化 测试 UI 框架 推荐
- 叶三《我们唱》-野孩子（白银饭店）
- DMS专线联通外网测试
- 立创开源 WCHLink/DapLink下载器 沁恒
- 深度学习核心词汇-英文
- 开什么玩笑？股票价格如何经得起AI的推敲？| 技术头条
- 我的2021年度书单（主要教你面试怎么装B）
热门文章
- ACPR'11 Accepted
- MySQL 大作业实训考试题_2020系统综合实践 期末大作业 15组
- 解决“手机能胜场使用校园网 笔记本电脑连接不上校园网或者连接上不可用”的问题
- Attempted to read or write protected memory. This is often an indication that other memory is corrup
- 查询注册表的命令行工具reg
- html视频怎么编辑倍速,浏览器flash/html5视频播放如何倍速(Enounce MySpeed)
- mb63.nte.ios.html,2009 Diagnosis, assessment, and treatment of non-pulmonary arterial hypertension
- Flink CDC 将MySQL的数据写入Hudi实践
- 【Linux】【操作】Linux环境运行Windows程序方式一览（全网最全）
- php 同步微信大量粉丝在数据表,微粉丝—— 微信加粉统计系统/复制统计准确率90%以上...