% 源程序
clear all
 
 
tic %开始计时,一般与toc连用,表示结束
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 调整参数,影响稳定性
a=1.3;%最大加速度%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%,调整对稳定性有重要影响
b=3;%期望减速度
T=1.6;%安全车头时距%
v00=33;%期望速度
s0=4;%拥堵车距
alpha=0.3;%三个0到1范围内的系数
beta=alpha;
gamma=alpha;
tau1=0.2;  tau2=0;%反应时间
s=5; % 车长
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  调整密度,影响稳定性
road=2000;
N=50;
h=road/N;%平均初始车头间距,对应的密度也不相同
tend=40000;
dt=0.1;
 
allv=cell(5,1);
alldx=cell(5,1);
allx=cell(5,1);
iN=4;
Be=1;
 
allv=cell(5,1);
alldx=cell(5,1);
allx=cell(5,1);
iN=4;
%% 随机参数 ,确定网联车的分布
% for
ii=2;
dx=[];v=[];dv=[];x=[];
df=0.2;
Jy=[1 0];%判断是否网联车,1是,0不是
p1=0.3;p2=1-p1; % 网联车的比p1;  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
prob1=[p1 p2];
L=length(Jy);%Jy的长度 1行2列所以为2
cc=zeros(L,1);%2行1列0阵
i=1;
FVm=[]; % 向量,确定每辆车是否是网联车
fprob=prob1;
prob=prob1;
fVm=Jy;   
fnum1=[];
while i<=N
    vm=randsrc(1,1,[fVm; fprob]);   %按照fprob的概率随机生成fVm中指定的数字 也就是0和1  1代表网联车概率为p1
    FVm(i)=vm;
    nn=find(abs(Jy-vm)<10^(-10));%满足条件的所有vm的所在位置
    fnum1(:,i)=nn(1)';%nn的第一个元素赋值给fnum1的第一列
    cc(nn(1))=cc(nn(1))+1;
    c1=find(abs(cc-N*prob')>1/2);
    fVm=Jy(c1);
    fprob=prob(c1)/sum(prob(c1));
    i=i+1;
end
AJc=zeros(1,N);   % 生成1xN的全0阵  向量,确定当前车和前车是否为网联车,以确定网联是否起效。注意:网联起效比例跟理论有一点差距的,理论值属于期望平均值。
%
for i=1:N
    if i==N
        if FVm(i)==1 && FVm(1)==1  %环形链接
            AJc(i)=1;
        else
            AJc(i)=0;
        end
    else
        if FVm(i)==1 && FVm(i+1)==1  %前后链接
            AJc(i)=1;
        else
            AJc(i)=0;
        end
    end
end
 
AJc=FVm;
%% 稳态车头间距
dx1=h;
aa=0; bb=100;
ve=(aa+bb)/2;
mm1=T*ve;
rf=a* (1-(ve/v00).^4-((s0+mm1)./(dx1-s)).^2);
while abs(rf)>10^(-20)
    if rf>0
        aa=ve;
    else
        bb=ve;
    end
    ve=(aa+bb)/2;
    %     v1=[1:30];
    mm1=T*ve;
    rf=a* (1-(ve/v00).^4-((s0+mm1)./(dx1-s)).^2);
end
 
%% 初始条件
v1=ones(1,N)*ve; dv1=zeros(1,N)*ve;   %稳态速度代入作为初始速度
dx1=ones(1,N)*h;  
x1=[];
for i=1:N
    x1(i)=sum(dx1(1:i-1));%每一列求和
end
Bi1=ceil(tau1/dt); Bi2=ceil(tau2/dt); %%%%%%向前取整
Be=Bi1+1;
x=[];v=[];dx=[];dv=[];
for i=1:Be
    x(i,:)=x1;
    v(i,1:N)=v1*0;
    dv(i,:)=dv1;
    dx(i,:)=dx1;
end
v(Be,1:2)=v(Be,1:2)+1;
  • 1.
  • 2.
  • 3.
  • 4.
  • 5.
  • 6.
  • 7.
  • 8.
  • 9.
  • 10.
  • 11.
  • 12.
  • 13.
  • 14.
  • 15.
  • 16.
  • 17.
  • 18.
  • 19.
  • 20.
  • 21.
  • 22.
  • 23.
  • 24.
  • 25.
  • 26.
  • 27.
  • 28.
  • 29.
  • 30.
  • 31.
  • 32.
  • 33.
  • 34.
  • 35.
  • 36.
  • 37.
  • 38.
  • 39.
  • 40.
  • 41.
  • 42.
  • 43.
  • 44.
  • 45.
  • 46.
  • 47.
  • 48.
  • 49.
  • 50.
  • 51.
  • 52.
  • 53.
  • 54.
  • 55.
  • 56.
  • 57.
  • 58.
  • 59.
  • 60.
  • 61.
  • 62.
  • 63.
  • 64.
  • 65.
  • 66.
  • 67.
  • 68.
  • 69.
  • 70.
  • 71.
  • 72.
  • 73.
  • 74.
  • 75.
  • 76.
  • 77.
  • 78.
  • 79.
  • 80.
  • 81.
  • 82.
  • 83.
  • 84.
  • 85.
  • 86.
  • 87.
  • 88.
  • 89.
  • 90.
  • 91.
  • 92.
  • 93.
  • 94.
  • 95.
  • 96.
  • 97.
  • 98.
  • 99.
  • 100.
  • 101.
  • 102.
  • 103.
  • 104.
  • 105.
  • 106.
  • 107.
  • 108.
  • 109.
  • 110.
  • 111.
  • 112.
  • 113.
  • 114.

异构交通流matlab源码_matlab