元胞自动机N a S c h模型
及其M A T L A B代码 This manuscript was revised by the office on December 22, 2012
元胞自动机N a S c h模型及其M A T L A B代码
作业要求
根据前面的介绍,对NaSch模型编程并进行数值模拟:
模型参数取值:Lroad=1000,p=0.3,Vmax=5。
边界条件:周期性边界。
数据统计:扔掉前50000个时间步,对后50000个时间步进行统计,需给出的结果。
基本图(流量-密度关系):需整个密度范围内的。
时空图(横坐标为空间,纵坐标为时间,密度和文献中时空图保持一致,画500个时间步即可)。
指出NaSch模型的创新之处,找出NaSch模型的不足,并给出自己的改进思路。
流量计算方法:
密度=车辆数/路长;
流量flux=density×V_ave。
在道路的某处设置虚拟探测计算统计时间T内通过的车辆数N;
流量flux=N/T。
在计算过程中可都使用无量纲的变量。
1、NaSch模型的介绍
作为对184号规则的推广,Nagel和Schreckberg在1992年提出了一个模拟车辆交通的元胞自动机模型,即NaSch模型(也有人称它为NaSch模型)。
时间、空间和车辆速度都被整数离散化。道路被划分为等距离的离散的格子,即元胞。
每个元胞或者是空的,或者被一辆车所占据。
车辆的速度可以在(0~Vmax)之间取值。
2、NaSch模型运行规则
在时刻t到时刻t+1的过程中按照下面的规则进行更新:
(1)加速:vnmin(vn1,vmax)
规则(1)反映了司机倾向于以尽可能大的速度行驶的特点。
(2)减速:vnmin(vn,dn)
规则(2)确保车辆不会与前车发生碰撞。
(3)随机慢化:以随机概率p进行慢化,令:vnmin(vn-1,0)
规则(3)引入随机慢化来体现驾驶员的行为差异,这样既可以反映随机加速行为,又可以反映减速过程中的过度反应行为。这一规则也是堵塞自发产生的至关重要因素。
(4)位置更新:vnxnvn,车辆按照更新后的速度向前运动。其中vn,xn分别表示第n辆车位置和速度;l(l≥1)为车辆长度;
p表示随机慢化概率;dnxn1xn1表示n车和前车n+1之间空的元胞数;
vmax为最大速度。
3、NaSch模型实例
根据题目要求,模型参数取值:L=1000,p=0.3,Vmax=5,用matlab软件进行编程,扔掉前11000个时间步,统计了之后500个时间步数据,得到如下基本图和时空图。
3.1程序简介
初始化:在路段上,随机分配200个车辆,且随机速度为1-5之间。
图
图3.1.1NaSch模型运行图
图3.1.2NaSch模型
3.2流量密度分析
图3.2描述了交通流量与密度的关系,从图中可知,该模型中,当密度为0——0.185时,流量随密度的增加而增加;当密度超过0.185时,流量开始随密度的增加而下降。
图3.2基于NaSch模型的流量密度图
3.3NaSch模型时空图分析
图
图NaSch模型的时空图
图NaSch模型的时空图
4模型评价
优点:该程序基本实现了NaSch模型的基本功能,并且最大速度、元胞数量、车辆数量以及运行间隔时间都可以修改,程序很灵活,并且可以清晰的看出每一次运行过程。
缺点:当时间步超过20000步时,内存占用量大。
附件
%主程序:NaSch_3.m程序代码
%单车道最大速度3个元胞开口边界条件加速减速随机慢化
clf
clearall
%buildtheGUI
%definetheplotbutton
plotbutton=uicontrol('style','pushbutton',...
'string','Run',...
'fontsize',12,...
'position',[100,400,50,20],...
'callback','run=1;');
%definethestopbutton
erasebutton=uicontrol('style','pushbutton',...
'string','Stop',...
'fontsize',12,...
'position',[100,500,50,20],...
'callback','freeze=1;');
%definetheQuitbutton
quitbutton=uicontrol('style','pushbutton',...
'string','Quit',...
'fontsize',12,...
'position',[100,600,50,20],...
'callback','stop=1;close;');
number=uicontrol('style','text',...
'string','1',...
'fontsize',12,...
'position',[20,400,50,20]);
%CAsetup
n=1000;%数据初始化
z=zeros(1,n);%元胞个数
z=roadstart(z,200);%道路状态初始化,路段上随机分布200辆
cells=z;
vmax=5;%最大速度
v=speedstart(cells,vmax);%速度初始化
x=1;%记录速度和车辆位置
memor_cells=zeros(3600,n);
memor_v=zeros(3600,n);
imh=imshow(cells);%初始化图像白色有车,黑色空元胞
set(imh,'erasemode','none')
axisequal
axistight
stop=0;%waitforaquitbuttonpush
run=0;%waitforadraw
freeze=0;%waitforafreeze(冻结)
while(stop==0&x<11502)
if(run==1)
%边界条件处理,搜素首末车,控制进出,使用开口条件
a=searchleadcar(cells);
b=searchlastcar(cells);
[cells,v]=border_control(cells,a,b,v,vmax);
i=searchleadcar(cells);%搜索首车位置
forj=1:i
ifi-j+1==n
[z,v]=leadcarupdate(z,v);
continue;
else
%======================================加速、减速、随机慢化
ifcells(i-j+1)==0;%判断当前位置是否非空
continue;
elsev(i-j+1)=min(v(i-j+1)+1,vmax);%加速
%=================================减速
k=searchfrontcar((i-j+1),cells);%搜素前方首个非空元胞位置ifk==0;%确定于前车之间的元胞数
d=n-(i-j+1);
elsed=k-(i-j+1)-1;
end
v(i-j+1)=min(v(i-j+1),d);
%==============================%减速
%随机慢化
v(i-j+1)=randslow(v(i-j+1));
new_v=v(i-j+1);
%======================================加速、减速、随机慢化%更新车辆位置
z(i-j+1)=0;
z(i-j+1+new_v)=1;
%更新速度
v(i-j+1)=0;
v(i-j+1+new_v)=new_v;
end
end
end
cells=z;
memor_cells(x,:)=cells;%记录速度和车辆位置
memor_v(x,:)=v;
x=x+1;
set(imh,'cdata',cells)%更新图像
%updatethestepnumberdisplay
pause(0.0001);
stepnumber=1+str2num(get(number,'string'));
set(number,'string',num2str(stepnumber))
end
if(freeze==1)
run=0;
freeze=0;
end
drawnow
end
figure(1)
forl=11001:1:11500
fork=1:1:1000
ifmemor_cells(l,k)>0
plot(k,l,'k.');
holdon;
end
end
end
xlabel('空间位置')
ylabel('时间(s)')
title('时空图')
fori=1:1:500
density(i)=sum(memor_cells(i,:)>0)/1000;
flow(i)=sum(memor_v(i,:))/1000;
end
figure(2)
plot(density,flow,'k.');
title('流量密度图')
xlabel('density')
ylabel('flow')
%
%//////////////////////////////////////////////////////////////////// ///
%
%
%函数:searchlastcar.m程序代码
function[location_lastcar]=searchlastcar(matrix_cells)
%搜索尾车位置
fori=1:length(matrix_cells)
ifmatrix_cells(i)~=0
location_lastcar=i;
break;
else%如果路上无车,则空元胞数设定为道路长度
location_lastcar=length(matrix_cells);
end
end
%函数:searchfrontcar.m程序代码
function[location_frontcar]=searchfrontcar(current_location,matrix_ce lls)
i=length(matrix_cells);
ifcurrent_location==i
location_frontcar=0;
else
forj=current_location+1:i
ifmatrix_cells(j)~=0
location_frontcar=j;
break;
else
location_frontcar=0;
end
end
end
%函数:roadstart.m程序代码
function[matrix_cells_start]=roadstart(matrix_cells,n)
%道路上的车辆初始化状态,元胞矩阵随机为0或1,matrix_cells初始矩阵,n初始车辆数k=length(matrix_cells);
z=round(k*rand(1,n));
fori=1:n
j=z(i);
ifj==0
matrix_cells(j)=0;
else
matrix_cells(j)=1;
end
end
matrix_cells_start=matrix_cells;
%函数:randslow.m程序代码
function[new_v]=randslow(v)
p=0.3;%慢化概率
rand('state',sum(100*clock)*rand(1));%¨ú×p_rand=rand;%产生随机概率ifp_rand<=p
v=max(v-1,0);
end
new_v=v;
%函数:leadcarrupdate.m程序代码
function[new_matrix_cells,new_=leadcarupdate(matrix_cells,v)%第一辆车更新规则
n=length(matrix_cells);
ifv(n)~=0
matrix_cells(n)=0;
v(n)=0;
end
new_matrix_cells=matrix_cells;
new_v=v;
%函数:searchleadcar.m程序代码
function[location_leadcar]=searchleadcar(matrix_cells)i=length(matrix _cells);
forj=1:i
ifmatrix_cells(i-j+1)~=0
location_leadcar=i-j+1;
break;
else
location_leadcar=0;
end
end
%函数:speadstart.m程序代码
function[v_matixcells]=speedstart(matrix_cells,vmax)%道路初始状态车辆速度初始化
v_matixcells=zeros(1,length(matrix_cells)); fori=1:length(matrix_cells)
ifmatrix_cells(i)~=0
v_matixcells(i)=round(vmax*rand(1));
end
end