M/M/S-混合制排队Matlab实现

%=======================================================================
% M/M/S-混合制排队系统
% 组成：
%    顾客：服从指数分布
%    服务台：多服务台，服务时间服从指数分布
%    服务机制：混合制（损失制、等待制）
%    服务规则：先到先服务(FCFS)
%
% 算法：计算机仿真(事件步长法)
%
% 输入参数：
%    totalNum_client：顾客总人数
%    lambda：单位时间顾客到达数
%    mu：单个服务台单位时间接待顾客数
%    num_serviceTable：服务台数量
%    maxNum_waitedQueue：单行等待顾客的最大排队容量
%
% 输出参数：
%    Pbusy：排队系统繁忙的概率
%    meaWait_peo：平均等待队长  
%    meaWait_time：平均等待时间
%    mean_peo：系统中平均总乘客数
%    lost_peo：损失的乘客数
%=======================================================================

% 清空工作区变量、清屏
clear
clc

% 初始化排队系统参数
totalNum_client = 100;    % 顾客总人数
lambda = 9;               % 单位时间顾客到达数
mu = 2;                   % 单个服务台单位时间接待顾客数
num_serviceTable = 4;     % 服务台数量
maxNum_waitedQueue = 10; % 单行等待顾客的最大排队容量

% 存储每位顾客的信息 (7*total_num_client 矩阵)
list_client = zeros(7, totalNum_client);  
% 产生指数分布的顾客到达时间间隔
list_client(1,:) = exprnd(1/lambda, 1, totalNum_client);  
list_client(1,:) = cumsum(list_client(1,:));            
list_client(2,:) = exprnd(1/mu, 1, totalNum_client); 
% 1.顾客到达时间
% 2.服务台顾客的服务时间
% 3.顾客等待时间
% 4.顾客离开时间
% 5.顾客到达时排队系统当前的人数（包括正在被服务的）
% 6.系统中正在等待中的人数
% 7.该顾客是否直接离开 (1;直接离开，0:排队等待)
pro_client = [];  % 储存每个顾客到达时在系统中总人数的序号

% 计算每个顾客到达时系统的队长
for i = 1:totalNum_client
    % 统计顾客到来时还在系统中的顾客的序号
    % 若顾客的离开时间>到达时间,则认为顾客还在排队系统中
    curPro_client = find(list_client(4,pro_client) > list_client(1,i));  
    % 更新当前系统中顾客的序号
    pro_client = pro_client(curPro_client);  
    % 统计排队系统中总顾客数
    num_queue = length(curPro_client);  
    % 损失制：若排队顾客数量超过单行总顾客数,则该顾客损失
    % 单行总顾客数 = 单行排队容量max_num_waitedQueue + 正在被服务的顾客数(即服务台数量)
    if( num_queue >= maxNum_waitedQueue + num_serviceTable )  
        list_client(2:7,i) = [0; 0; list_client(1,i); num_queue; num_queue-num_serviceTable; 1];
        continue
    else if( num_queue < num_serviceTable )  % 无人排队等待,则等待时间为0
            list_client(3:6,i) = [0, sum(list_client(1:3,i)), num_queue, 0];
            pro_client = [pro_client,i];  % 将当前顾客加入到系统
        else  % 有人排队
            leaveTime_client = sort(list_client(4, pro_client));  % 排序：正在排队系统中顾客的离开时间
            % 排队等待的顾客 = 顾客到来时还在系统中的人的序号cur_pro_client - 服务台数量num_serviceTable
            list_client(3,i) = leaveTime_client(end - num_serviceTable + 1) - list_client(1,i);  
            % 此时,顾客的等待时间 = 当前系统中倒数第num_serviceTable(服务台数量)个离开的顾客的时间 - 该顾客的到达时间
            list_client(4:6,i) = [sum(list_client(1:3,i)); num_queue; num_queue-num_serviceTable];
            pro_client = [pro_client,i];  % 将当前的顾客加入到系统
        end
    end
end

% 计算时间范围内排队系统的队长变化
list_time = zeros(4, 2*totalNum_client);
% 1.时间
% 2.对应时间的总人数
% 3.对应时间正在排队的人数
% 4.是否繁忙(1:繁忙 0:不繁忙)

% 系统队长的变化必定发生在顾客到来或离去的时刻
% 排序：顾客的到来时间与离开时间
[list_time(1,:),idxs] = sort([list_client(1,:), list_client(4,:)]); 
% 第一位顾客到达,系统中人数开始变为 1,正在排队等待的为 0,系统不繁忙
list_time(2,1) = 1; 
for i = 2:2*totalNum_client    
    if idxs(i) <= totalNum_client  % 此时刻是有顾客到达
        list_time(2,i) = list_time(2,i-1)+1;
    else  % 此时刻有顾客离去
        list_time(2,i) = list_time(2,i-1) - 1;
    end
    % 判断是否繁忙（繁忙:所有的服务台全满）
    if list_time(2,i) >= num_serviceTable
        % 若繁忙则得出正在排队的人数，繁忙
        list_time(3:4,i) = [list_time(2,i)-num_serviceTable, 1];  
    end
end

% 计算各种统计指标(忽略首尾的非稳态过程)
all_peo = 0;
all_wait = 0;
busyTime = 0;
for i = 1:(2*totalNum_client-1)
    block = list_time(1,i+1) - list_time(1,i);     % 相邻两个“转变时刻”的时间间隔
    all_peo = all_peo + block * list_time(2,i);    % 每个时间段的总人数之和
    all_wait = all_wait + block * list_time(3,i);  % 每个时间段的等待人数之和
    busyTime = busyTime + block * list_time(4,i);  % 总繁忙时刻
end
Pbusy = busyTime / list_time(1,end);        % 系统繁忙的概率
mean_peo = all_peo / list_time(1,end);      % 系统中的平均总人数=每个时间段的总人数之和/总时间
meanWait_peo = all_wait / list_time(1,end); % 平均等待人数=每个时间段的等待人数之和/总时间
meanWait_time = sum(list_client(3,:)) / totalNum_client;  % 平均等待时间=每个人的等待时间之和/总人数
lost_peo = sum(list_client(7,:));           % 失去的顾客数
num_queue = find(list_client(7,:)==1);      % 失去的顾客的序号

% 每个顾客的到达时间和离开时间
figure
plot(1:totalNum_client, list_client(1,:), 1:totalNum_client, list_client(4,:), num_queue, list_client(1,num_queue), '*r');
title('顾客到达与离开时间')
% legend('到达时间 ','离开时间 ','损失的顾客');
% 各顾客等待时间与停留总时间
figure
plot(1:totalNum_client, list_client(3,:), 1:totalNum_client, sum(list_client(2:3,:)), num_queue,list_client(3,num_queue), '* r')
% legend('等待时间','停留总时间','损失的顾客',0)
title('各顾客等待时间与停留总时间')
xlabel('顾客');  ylabel('时间')
% 每个顾客到来时的系统总人数
figure
plot(1:totalNum_client, list_client(5,:), ':.', num_queue, list_client(5,num_queue), '*r')
% legend('每个顾客到来时的系统总人数','损失的顾客',0)
xlabel('顾客');  ylabel('系统总人数')
% 系统中总人数随时间的变化
figure
stairs(list_time(1,:), list_time(2,:));
hold on
plot(list_client(1,totalNum_client), list_time(2,idxs==totalNum_client), '.r', 'MarkerSize', 20);
title('系统中总人数随时间的变化')
if sum(list_client(7,:))
    plot([0,list_time(1,end)], [maxNum_waitedQueue + num_serviceTable, maxNum_waitedQueue + num_serviceTable],'--r')
    text(list_time(1,3), maxNum_waitedQueue + num_serviceTable, '系统最大容量')
end
% legend('人数','最后一个顾客到达',0)
xlabel('时间');  ylabel('人数')
fprintf('排队系统繁忙概率    =%.3f\n', Pbusy);
fprintf('平均等待队长        =%.3f\n', meanWait_peo);
fprintf('平均等待时间        =%.3f\n', meanWait_time);
fprintf('排队系统中平均总人数=%.3f\n', mean_peo);
fprintf('损失的顾客数        =%d\n',   lost_peo);
fprintf('损失率              =%.3f\n', lost_peo/totalNum_client);