Search code examples
rsimulationmarkovpoison-queue

Birth Death Code


My question:

At an institute for experimental mathematics there is a computer that helps solve problems.

Problems arrive to the computer at a poisson process with intensity "Landa" per hour.

The time to solve each problem can be seen as a exponential distribution with parameter "mu".

In our world we have four different states. S = (0,1,2,3)

State 0 = 0 problems have arrived to the computer

State 1 = the computer is solving 1 question

State 2 = the computer is solving 1 question + 1 in queue.

State 3 = the computer is solving 1 question + 2 in queue.

If a question comes when we are in state 3, the sender gets an error message and tries again later. The institution has decided that maximum of 5 % of the senders should get this error message.

To decide who should have access to the computer, we are 3 different proposals.

  1. Only the professors are alowed to send questions ( Landa = 2, Mu = 10)
  2. Professors and students are alowed to send questions( Landa = 6, Mu = 10)
  3. Anyone is alowed to send questions( Landa =10, Mu = 10)

We should investigate what of the 3 proposals don't fill upp the computer more than 5% of the time.

I have two things i need help with

First thing: To Solve the question I've been given this structure of code(the code under). What i need help with is if someone can briefly explain to me the purpose of the following code paragraf where i have written "#?".

So what i really need help with is some1 to explain parts of the code.

Secound thing: In two places i have written "...", there i need help to fill in some code.

bd_process <- function(lambda, mu, initial_state = 0, steps = 100) {
    time_now <- 0
    state_now <- initial_state
    time <- 0
    state <- initial_state

    for (i in 1:steps) {

        if (state_now == 3) {
            lambda_now <- 0
        } else {
            lambda_now <- lambda
        }

        if (state_now == 0) {
            mu_now <- 0
        } else {
            mu_now <- mu
        }

        #?
        time_to_transition <- ...

        #? 
        if (...) {
            state_now <- state_now - 1
        } else {
            state_now <- state_now + 1
        }

        #?
        time_now <- time_now + time_to_transition 
        #?
        time <- c(time, time_now) 
        #?
        state <- c(state, state_now) #WHAT DOE THIS VECTOR CONSIST OF?
    }
    list(time = time, state = state)
}

Solution

  • The code appears to be written with an implicit assumption that the interarrival and service distributions are memoryless, i.e., either exponential or geometric. Without memorylessness it's invalid to turn off processing by setting the rates to zero.

    With the memoryless property, you can figure out the time_to_transition as a superposition of the two Poisson processes, and determine whether it was an arrival or a departure by randomizing proportional to the ratio of one of the rates to the combined rate. It's also valid to zero one of the rates then because when you unzero it the elapsed time where the rate was zero doesn't matter due to the memoryless property.