%mm1 simulation in matlab
%有优先级的MMm排队系统
%强插型
clc;clear;
ST_Idle=0;
ST_Busy=1;
EV_NULL=0;
EV_Arrive=1;
EV_Depart=2;
EV_LEN=3;
Q_LIMIT=1000;
% next_event_type=[];
% next_depart=[] %因为有m个窗口,所以要确定下一个离开的是哪一个窗口的顾客
% num_custs_delayed=[]; %已到达的顾客数目
% num_delays_required=[];
% num_events=[];
% num_in_q=[];
% server_status=[];
% area_num_in_q=[];
% area_server_status=[];
% mean_interarrival=[];
% mean_service=[];
% sim_time=[];
% time_last_event=[];
% total_of_delays=[];
%
time_arrival=[]; %顾客到达时刻,只是记录各级等待队列中顾客的到达时刻 ,该数组可能为空
time_arrival_II=[];
time_arrival_III=[];
time_next_event=zeros(1,EV_LEN);
%仿真参数
num_events=EV_LEN-1;
num_server=2; %M/M/m m=2
% mean_interarrival=[30 300 100 50];%总的平均到达间隔和各类平均到达间隔
mean_interarrival=[30 150 100 60];%总的平均到达间隔和各类平均到达间隔
mean_service=20;
num_delays_required=20000; %仿真要求到达的顾客数目
outfile=fopen('doctor.txt','w');
fprintf(outfile, 'Multiple-doctor queueing system\n\n');
fprintf(outfile, 'Mean interarrival time of all patients%26.3f minutes\n\n',mean_interarrival(1));
fprintf(outfile, 'Mean interarrival time of I patients%26.3f minutes\n\n',mean_interarrival(2));
fprintf(outfile, 'Mean interarrival time of II patients%26.3f minutes\n\n',mean_interarrival(3));
fprintf(outfile, 'Mean interarrival time of III patients%26.3f minutes\n\n',mean_interarrival(4));
fprintf(outfile, 'Mean service time%26.3f minutes\n\n', mean_service);
fprintf(outfile, 'Number of servers%26d\n\n', num_server);
fprintf(outfile, 'Number of customers%26d\n\n\n\n\n', num_delays_required);
%%%%%%%%%%%%%%%%%%%%%%%%% part1 initialize
sim_time=0.0;
%/* Initialize the state variables. */
server_status =zeros(2,num_server); %窗口状态
server_status(2,:)=20; %正在服务的顾客级别,初始化为最低
num_in_q = zeros(1,4); %解释类似于 total_of_delays
time_last_event = 0.0;
%/* Initialize the statistical counters. */
num_custs_delayed = zeros(1,4); %解释类似于 total_of_delays
total_of_delays = zeros(1,4); %第一个是所有级别客人的总等待时间,后三个各个级别客人的总的等待时间
total_of_time = 0.0;
area_num_in_q = zeros(1,4);
area_server_status = 0.0;
%/* Initialize event list. Since no customers are present, the departure
%(service completion) event is eliminated from consideration. */
time_next_arrival=randexp(mean_interarrival);
[time_next_arrival(1),rank_next_arrival]=min(time_next_arrival(2:4));
server_status(2,1)=rank_next_arrival; %第一个窗口顾客级别即为第一个将要到达的顾客等级
time_next_event(EV_Arrive) = sim_time + time_next_arrival(1);
time_next_event(EV_Depart) = 1.0e+230;
time_depart=zeros(1,num_server);
next_depart=0; %有顾客将要离开的窗口
%%%%%%%%%%%%%part2 Run the simulation while more delays are still needed.
while (num_custs_delayed(1) < num_delays_required)
%/* Determine the next event. */
min_time_next_event = 1.0e+290;
next_event_type = 0;
% Determine m depart event
min_time_depart=1e290;
next_depart=-1;
for i=1:num_server
if(server_status(i)==1 & time_depart(i)<min_time_depart)
min_time_depart=time_depart(i);
next_depart=i;
end
end
time_next_event(EV_Depart)=min_time_depart;
%
%/* Determine the event type of the next event to occur. */
for i = 1: num_events
if (time_next_event(i) < min_time_next_event)
min_time_next_event = time_next_event(i);
next_event_type = i;
end
end
%/* Check to see whether the event list is empty. */
if (next_event_type == 0)
fprintf(outfile, '\nEvent list empty at time %f', sim_time);
exit(1);
end
%/* The event list is not empty, so advance the simulation clock. */
sim_time = min_time_next_event;
%/* Update time-average statistical accumulators. */
double time_since_last_event;
%/* Compute time since last event, and update last-event-time marker. */
time_since_last_event = sim_time - time_last_event;
time_last_event = sim_time;
%/* Update area under number-in-queue function. */要考虑优先级
area_num_in_q(1)=area_num_in_q(1) + num_in_q(1) * time_since_last_event;
area_num_in_q(2)=area_num_in_q(2) + num_in_q(2) * time_since_last_event;
area_num_in_q(3)=area_num_in_q(3) + num_in_q(3) * time_since_last_event;
area_num_in_q(4)=area_num_in_q(4) + num_in_q(4) * time_since_last_event;
%/* Update area under server-busy indicator function. */
for i=1:num_server
area_server_status =area_server_status + server_status(1,i) * time_since_last_event;
end
%arrival
if(next_event_type==EV_Arrive)
% %/* Schedule next arrival. */
% time_next_arrival=randexp(mean_interarrival);
% [time_next_arrival(1),rank_next_arrival]=min(time_next_arrival(2:4));
% rank_next_arrival=rank_next_arrival-1;
% time_next_event(EV_Arrive) = sim_time + time_next_arrival(1);
% time_next_event(1) = sim_time + randexp(mean_interarrival);
%/* Check to see whether server is busy. */
s_idle=-1;
for i=1:num_server
if (server_status(1,i) == ST_Idle)
s_idle=i;
break;
end
end
%所有窗口忙,考虑优先级,强插情形
if(s_idle== -1 )
if(rank_next_arrival < max(server_status(2,:))) %可以挤掉正在服务的最低级顾客
[temp_rank,temp_server]=max(server_status(2,:));
num_in_q(1)=1+num_in_q(1); %所有排队队列中人数加一
num_in_q(1+temp_rank)=1+num_in_q(1+temp_rank);
if (num_in_q(1) > Q_LIMIT) %检查等待队列是否满
%/* 等待队列溢出,停止仿真*/
fprintf(outfile, '\nOverflow of the array time_arrival at');
fprintf(outfile, ' time %f', sim_time);
exit(2);
else
if temp_rank==1 %被挤掉顾客放在相应等待队列首位,记录时间
time_arrival=[sim_time,time_arrival];
elseif temp_rank==2;
time_arrival_II=[sim_time,time_arrival_II];
else
time_arrival_III=[sim_time,time_arrival_III];
end
server_status(2,temp_server)=rank_next_arrival;
end
%不能挤掉已有顾客
else
num_in_q(1)=1+num_in_q(1);
num_in_q(1+rank_next_arrival)=1+num_in_q(1+rank_next_arrival); %相应等级的等待队列中加一
%/* Check to see whether an overflow condition exists. */
if (num_in_q(1) > Q_LIMIT)
%/* 等待队列溢出,停止仿真*/
fprintf(outfile, '\nOverflow of the array time_arrival at');
fprintf(outfile, ' time %f', sim_time);
exit(2);
else
%/* There is still room in the queue, so store the time of arrival of the arriving customer at the (new) end of time_arrival. */