mmm排队系统性能仿真代码matlab.txt

M/M/m排队系统性能仿真代码matlab function r=randexp(lambda) r = -lambda*log(rand)

M/M/m排队系统性能仿真代码matlab 
function r=randexp(lambda)        r = -lambda*log(rand);

说明同上

%mmm simulation in matlab
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=[]
% 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_next_event=zeros(1,EV_LEN);
%仿真参数
num_events=EV_LEN-1;
num_server=3;                     %M/M/m  m=2
mean_interarrival=1;
mean_service=.5;
num_delays_required=2000;           %

outfile=fopen('MM m.txt','w');
fprintf(outfile, 'Multiple-server queueing system\n\n');
fprintf(outfile, 'Mean interarrival time%11.3f minutes\n\n',mean_interarrival);
fprintf(outfile, 'Mean service time%16.3f minutes\n\n', mean_service);
fprintf(outfile, 'Number of servers%20d\n\n', num_server);
fprintf(outfile, 'Number of customers%14d\n\n', num_delays_required);

 

 

%part1 initialize
sim_time=0.0;
    %/* Initialize the state variables. */

    server_status   =zeros(1,num_server);   %idle
    num_in_q        = 0;
    time_last_event = 0.0;
    
    %/* Initialize the statistical counters. */

    num_custs_delayed  = 0;
    total_of_delays    = 0.0;
    total_of_time    = 0.0;
   
    area_num_in_q      = 0.0;
    area_server_status = 0.0;
    
    %/* Initialize event list.  Since no customers are present, the departure
    %(service completion) event is eliminated from consideration. */
    time_next_event(EV_Arrive) = sim_time + randexp(mean_interarrival);
    time_next_event(EV_Depart) = 1.0e+230;
    
    time_depart=zeros(1,num_server);
    next_depart=0;

 
    
    
    
%part2 
while (num_custs_delayed < num_delays_required)
%Run the simulation while more delays are still needed.
    %/* 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(2)=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) 
        
    %/* The event list is empty, so stop the simulation. */

        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=area_num_in_q +  num_in_q * time_since_last_event;

    %/* Update area under server-busy indicator function. */
    
    for i=1:num_server
        area_server_status =area_server_status + server_status(i) * time_since_last_event;
    end
    
    %/* Invoke the appropriate event function. */
%arrival

    if(next_event_type==EV_Arrive)
        double delay;
        
    %/* Schedule next arrival. */

        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(i) == ST_Idle)
                s_idle=i;
                break;
            end
        end
            
        %/* all Server is busy, so increment number of customers in queue. */
        if(s_idle== -1 )
            num_in_q=1+num_in_q;
            
        %/* Check to see whether an overflow condition exists. */

            if (num_in_q > Q_LIMIT)                 
                %/* The queue has overflowed, so stop the simulation. */

                fprintf(outfile, '\nOverflow of the array time_arrival at');
                fprintf(outfile, ' time %f', sim_time);
                exit(2);
            end
                %/* There is still room in the queue, so store the time of arrival of the arriving customer at the (new) end of time_arrival. */

            time_arrival(length(time_arrival)+1)=sim_time;

        else 
        %/* Server is idle, so arriving customer has a delay of zero.  (The following two statements are for program clarity 
        %and do not affect the results of the simulation.) */

            delay = 0.0;
            total_of_delays =total_of_delays + delay;
            
        %/* Increment the number of customers delayed, and make server busy. */

            num_custs_delayed = 1 + num_custs_delayed;
            server_status(s_idle) = ST_Busy;
            
        %/* Schedule a departure (service completion). */

            time_depart(s_idle) = sim_time + randexp(mean_service);
        end %    if (server_status == ST_Busy) 
%depart
    else
         double delay;

        %/* Check to see whether the queue is empty. */

        if (num_in_q == 0) 

        % /* The queue is empty so make the server idle and eliminate the departure (service completion) event from consideration. */

            server_status(next_depart)      = ST_Idle;
            time_depart(next_depart) = 1.0e+230;
            %check_depart()
            min_time_dapart=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(2)=min_time_depart;
        else

        %/* The queue is nonempty, so decrement the number of customers in queue. */

            num_in_q=num_in_q-1;

        %/* Compute the delay of the customer who is beginning service and update the total delay accumulator. */

            delay = sim_time - time_arrival(1);
            total_of_delays =total_of_delays + delay;

        %/* Increment the number of customers delayed, and schedule departure. */

            num_custs_delayed = 1 + num_custs_delayed;
            serv_time=randexp(mean_service);
            time_depart(next_depart) = sim_time + serv_time;
            total_of_time = total_of_time + delay + serv_time;
            %check_depart()
            min_time_dapart=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(2)=min_time_depart;
            
        %/* Move each customer in queue (if any) up one place. */

            tempForPop=time_arrival(2:length(time_arrival));
            time_arrival=tempForPop;
        end %if (num_in_q == 0)
        
    end %if(next_event_type==EV_Arrive)
    
end %while
               

 

 


%%%%%%%%%% part 3
%/* Invoke the report generator and end the simulation. */

    fprintf(outfile, '\n\nAverage delay in queue%11.3f minutes\n\n',total_of_delays / num_custs_delayed);
    fprintf(outfile, '\n\nAverage delay in system%11.3f minutes\n\n',total_of_time / num_custs_delayed);
    
    fprintf(outfile, 'Average number in queue%10.3f\n\n',area_num_in_q / sim_time);
    fprintf(outfile, 'Server utilization%15.3f\n\n',area_server_status / sim_time);
    fprintf(outfile, 'Time simulation ended%12.3f minutes', sim_time);
    fclose(outfile);