finite_queue.m

基于matlab的m/m/n排队仿真源码

% Title      : Finite Source Queue

% Authors    : Dr. Antonio A. Trani (231-4418, [e-mail])
%              Dr. Hojong Baik      (231-7021, [e-mail]) 
% Date       : Nov. 24, 1999
% Revisions  : 
% 				1) changed wording in K parameter (Trani: October 2003)
% 				2) Changed function call to FCT.m to avoid use of FCT. Instead use of the Matlab
% 					factorial function (FACTORIAL) - October 23, 2003 (Trani)

% Course No. : CEE3604

%
%     Given    :
%       
%     1) Earthwork operation using s power dozers loading n trucks
%          ie, # of servers (S)       : 2
%              # of calling units (n) : 6
%     2) Mean time for a single calling unit in the queueing system 
%        to complete a cycle (T1) : 30 (min)
%     3) Mean time it takes a calling unit to be served (T2) : 6 (min)
%     4) Random variable T1, T2 are conforming 
%        to a negative exponential distribution   
%
% Typical Questions asked about these systems:    
%
%     1) the fraction of time that a dozer is idle (FTI)
%     2) the expected # of trucks that are not idle (K)
%     3) the operating efficiency for trucks (K/N)
%     4) the operating efficiency for dozers (R/S)
%

%  
%    Equation of steady-state probability that n calling units are in the queueing system
%
%       Pn= factorial(N)/(factorial(n)*factorial(N-n))*(Lambda/Myu)^n*Po,         if n<=S 
%         = factorial(N)/(factorial(S)*factorial(N-n)*S^(n-S))*(Lambda/Myu)^n*Po, if S<=n<=N 
%
%    Steady-state MOEs(measures of effectiveness)
%
%       L  : average # of C.U. in the Q-S. 
%       Lq : average # of C.U. in the Queue
%       N  : total # of calling units
%       K  : average # of C.U. that are working (=N-L)
%       R  : average # of servers that are working (=L-Lq=S-Sum[(s-n)*Pn, n=0..s-1])
%
%       K/N : operating efficiency for C.U.
%       R/S : operating efficiency for servers
%       FTI = Fraction of idle time 
%             (=1-R/S=P0+((s-1)/s)*p1+....+(1/s)*P(s-1))
%
%  Definition OF Variables
%
%    Mu    : Service rate (=1/T1*60, trucks/hr) 
%    Lambda : arrival rate (=1/T2*60, trucks/hr)
%    Pn     : Probability n trucks are in the Q-S(Queueing System) 
%    OE_CU  : Operating efficiency for C.U. (=K/N)
%    OE_Ser : Operating efficiency for servers (=R/S)
%    < See the above problem description >

% Source Code

% Initial conditions
  S=2;
  N=10;			% total number of trucks
  T1=30; 		% cycle time ( minutes)
  T2=5;			% service time
  Lambda=1/T1*60;
  Mu=1/T2*60;

% Prob. that there is no calling units in the queueing system (Po)
    Po_inverse=0;
%    sum_den=0;

  for n=0:S % for the first part of denominator (den_1) 
    den_1=factorial(N)/(factorial(n)*factorial(N-n))*(Lambda/Mu)^n;
    Po_inverse=Po_inverse+den_1;
    disp(den_1);
  end

  for n=S+1:N  % for the second part of denominator (den_2) 
    den_2=factorial(N)/(factorial(S)*factorial(N-n)*S^(n-S))*(Lambda/Mu)^n;
    Po_inverse=Po_inverse+den_2;
    disp(den_2);
  end


    Po=1/Po_inverse;
%
% Steady-state MOEs(measures of effectiveness)
%

% Fraction of idle time for C.U. (FTI)
  FTI=0;
  for n=0:S-1
      Pn=factorial(N)/(factorial(n)*factorial(N-n))*(Lambda/Mu)^n*Po;
      FTI=FTI+1/S*(S-n)*Pn
  end

% L : average # of C.U. in the Q-S. 
  L=0;
  P_sum=0
  for n=0:S     
      P(n+1,1)=n;
      P(n+1,2)=factorial(N)/(factorial(n)*factorial(N-n))*(Lambda/Mu)^n*Po;
      P(n+1,3)=n*P(n+1,2);
      P_sum=P_sum+P(n+1,2);
      L=L+P(n+1,3); 
  end

  for n=S+1:N     
      P(n+1,1)=n;
      P(n+1,2)=factorial(N)/(factorial(S)*factorial(N-n)*S^(n-S))*(Lambda/Mu)^n*Po;
      P(n+1,3)=n*P(n+1,2);
      P_sum=P_sum+P(n+1,2);
      L=L+P(n+1,3);       
  end

% K : average # of C.U. that are working (=N-L)
  K=N-L;

% R : average # of servers that are working (=L-Lq)
  R=S;
  for n=0:S-1
      Pn=factorial(N)/(factorial(n)*factorial(N-n))*(Lambda/Mu)^n*Po;
      R=R-(S-n)*Pn
  end

% K/N : operating efficiency for C.U.
  OE_CU=K/N;

% R/S : operating efficiency for servers
  OE_Ser=R/S;
  
% Average arrival rate of entities is

	lambda_bar = Lambda * K;

  clc
  disp([blanks(5),'Limited Source Queueing Model '])
  disp(' ')
  disp(' ')
  disp([blanks(5),'   arrival rate (trucks/hr) = ', num2str(Lambda)])
  disp([blanks(5),'   service rate (trucks/hr) = ', num2str(Mu)])
  disp(' ')

  disp([blanks(5),'   Po                           = ', num2str(Po)])
  disp([blanks(5),'1) Fraction of idle time (FTI)  = ', num2str(FTI)])
  disp(' ')  
  disp([blanks(5),' ----------    -------------    ---------'])
  disp([blanks(5),'  # of C.U.   Prob. of n C.U.    ']) 
  disp([blanks(5),' in the Q-S     in the Q-S    ']) 
  disp([blanks(5),'    (n)           (Pn)            (n*Pn)  ']) 
  disp([blanks(5),' ----------    -------------    ---------']) 
  for i=1:N+1
      disp([blanks(10),num2str(P(i,1)),'           ',num2str(P(i,2)),'          ',num2str(P(i,3))]) 
  end
  disp([blanks(5),' ----------    -------------    ---------']) 
  disp([blanks(5),'  SUM               ',num2str(P_sum),'              ',num2str(L),' (=L)']) 
  disp([blanks(5),' ----------    -------------    ---------'])
  disp(' ')
  disp([blanks(5),'2) Expected # of working trucks(K=N-L) = ', num2str(K)])
  disp([blanks(5),'3) Operating efficiency of truck (K/N) = ', num2str(OE_CU)])
  disp(' ')
  disp([blanks(5),'4) Expected # of dozers that are working (R) = ', num2str(R)])
  disp([blanks(5),'   Operating efficiency of a dozer (R/S) = ', num2str(OE_Ser)])
  disp(' ')
    disp([blanks(5),' Average arrival rate = ', num2str(lambda_bar)])
  disp(' ')
  
  pause
  
  stem(P(:,1),P(:,2),'o')
  xlabel('Number of Trucks in the System')
  ylabel('Probability')
  grid