ece610.cpp

一个模拟泊松随机过程的小程序。可以计算系统平均等待时间

// ECE610.cpp: implementation of the ECE610 class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "ECE610.h"
#include "stdio.h"
#include "math.h"
#include "time.h"
#include <deque>
#include <queue>

using namespace std ;
typedef queue<int>  INTQUEUE;


//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

ECE610::ECE610()
{
	printf("happy start!\n");



	erlangRV=0;
	exponRV=0;
	uniRV=0;
	expectedPn=0;
	totalArrival=0; 
	PTime=0;

	for(int j=0;j<1000;j++)
	{
		Ni[j]=0;
		timeI[j]=0;

	}
	currentDepartureTime=0;
	currentArrivalTime=GetExponRv(lamda);

	//NextDepartureTime=currentArrivalTime+1.0/6.0;  //deterministic
	NextDepartureTime=currentArrivalTime+GetErlangRv(mu);

	NextArrivalTime=currentArrivalTime+GetExponRv(lamda);
	//NextArrivalTime=currentArrivalTime+GetErlangRv(lamda);

	queue.push(1);
	totalArrival++;
	Print(false, queue.size()-1);


	fp = fopen("\\log.txt","w");

	Simulation();
	ComputationResult();

}


ECE610::~ECE610()
{

}


float ECE610::GetExponRv(float rate)
{
	uniRV=rand()*1.0/RAND_MAX;
	//generating exponential random variables according to the uniform r.v.
	exponRV=-1*(1.0/rate)*log(1-uniRV);
	
	return exponRV;
}


float ECE610::GetErlangRv(float rate)
{

	//generating uniform random variables
    int i;
	erlangRV = 0;
	for(i=0;i<factorK;i++)
	{
	    erlangRV=erlangRV + GetExponRv(factorK*rate);
	}
	
	return erlangRV;

}

void ECE610::Simulation()
{
	int i=0;
	bool flag=true;
	while(i<10001)
	{
		i++;
		if(NextDepartureTime > NextArrivalTime  ) //a new arrival
		{
plea:
		     if(queue.size()==0)
			 {
   			     timeI[0] = timeI[0] + NextArrivalTime - NextDepartureTime;
			 }
			 queue.push(2);
             timeI[queue.size()] = timeI[queue.size()]+ NextArrivalTime - currentArrivalTime;

			 totalArrival++;
			 currentArrivalTime=NextArrivalTime;
			 Print(false, queue.size()-1);
			 NextArrivalTime=currentArrivalTime+GetExponRv(lamda);
			 //NextArrivalTime=currentArrivalTime+GetErlangRv(lamda);

			 if(queue.size()==1)
			 {
				 //NextDepartureTime=currentArrivalTime+GetExponRv(mu);
                 NextDepartureTime=currentArrivalTime+GetErlangRv(mu);
                 //NextDepartureTime=currentArrivalTime+0.2;
			 }

             continue;
		}
		else if(NextDepartureTime < NextArrivalTime ) // a new departure
		{

			currentDepartureTime=NextDepartureTime;
			queue.pop();
            
			timeI[queue.size()] = timeI[queue.size()]+NextDepartureTime - currentDepartureTime;

			Print(true, queue.size()-1);
			if(queue.size()==0)    
			{
                //flag = false;
				goto plea;
			}
			else if(queue.size()>0)
			{
				//NextDepartureTime=currentDepartureTime+GetExponRv(mu);
                NextDepartureTime=currentArrivalTime+GetErlangRv(mu);
                //NextDepartureTime=currentArrivalTime+0.2;


			}
            continue;
		}
		else if ( NextDepartureTime == NextArrivalTime )
		{

			NextDepartureTime=NextDepartureTime+GetErlangRv(mu);
			//NextDepartureTime=NextDepartureTime+GetExponRv(mu);
            //NextDepartureTime=currentArrivalTime+0.2;

			NextArrivalTime=NextArrivalTime+GetExponRv(lamda);
        }

	}


}

void ECE610::Print(bool departure, int num)
{

	if(departure)
        printf("\n %c, %f", 'd', currentDepartureTime);
    else
		printf("\n %c, %f", 'a', currentArrivalTime);

	for(int i=0;i< num+1;i++)
	{
		printf(" %c",'+');
	}

	if(!departure)
	{
		Ni[num]++;
	}

	
}

void ECE610::ComputationResult()
{
	for(int j=0; j<50;j++)
	{
	    expectedPn=expectedPn+ j * Ni[j]*1.0 / totalArrival;
		PTime=PTime+j*timeI[j];
		//printf("\n prob of state %d is %.4f", j, Ni[j] *1.0/ totalArrival);

		fprintf(fp," %f,\n",Ni[j] *1.0/ totalArrival);
	}
	fclose(fp);
	
	printf("\n Time average is %.4f", PTime);///currentArrivalTime);
    printf("\n%s %.4f\n", "Expected Number:", expectedPn );

}