Erlang B

ErlangB_simu.jl

using Compat, Random, Distributions

# to set a markov chain, use a dictionary data structure like this.
example_set = Dict("lambda"=>19.5, "mu"=>1, "k"=>10, "max_na"=>500000)

# Uniform distribution
uniform_d = Uniform()

function Markov_Chain(paratmers)
    global uniform_d
    # initialize paratmers
    # queue size
    Q_size = 0
    # N_a is number of customer arrivals counted so far
    # N_b is number of blocked customers counted so far
    N_a, N_b = 0, 0
    lambda = paratmers["lambda"]
    mu = paratmers["mu"]
    k = paratmers["k"]
    max_na = paratmers["max_na"]
    rand_X = rand(uniform_d, max_na)
    for x in rand_X
        if x < lambda/(lambda+Q_size*mu)
            N_a = N_a + 1
            if Q_size == k # when Q_size is equal to k, the next arrival will be blocked
                N_b = N_b + 1
                #println("Blocking happens")
            else
                # goes into queue
                Q_size = Q_size + 1
            end
        else
            # quit queue
            Q_size = Q_size - 1
        end
    end
    return N_b/N_a
end


function erlangB(test_case)
    result_set = []
    for i=1:11
        x = Markov_Chain(test_case)
        push!(result_set, x)
    end 
    E_x = mean(result_set)
    Var_x = var(result_set)
    sigma = sqrt(Var_x)
    Ur = 2.23*(sigma/sqrt(11))
    lower_bound = E_x-Ur
    upper_bound = E_x+Ur
    length_CI = 2*Ur
    if (2*Ur < (0.05*E_x))
        condition_flag = true
    else
        condition_flag = false
    end
    println("Arrival number for each experiment: $(test_case["max_na"])")
    println("Mean blocking probability: $(E_x)")
    println("Confidence Interval 95%: ($(lower_bound), $(upper_bound))")
    println("The Confidence Interval length is less than 5%: $(condition_flag)\n")

    # return condition flag for judge whether needs next loop
    return condition_flag
end

# replace these paratmers for different kind of experiment
test_lambda = 20
test_mu = 1
test_k = 18
max_na = 100

constantC = (test_lambda/test_mu)/test_k
if (constantC==1)
    println("constantC=1")
elseif (constantC>1)
    println("constantC>1")
else
    println("constantC<1")
end

while (true)
    global test_lambda,test_mu,test_k,max_na
    test_case = Dict("lambda"=>test_lambda, "mu"=>test_mu, "k"=>test_k, "max_na"=>max_na)
    condition_flag = erlangB(test_case)
    if condition_flag
        break
    else
        # multiply max_na with 10 for next experiment
        max_na = max_na * 10
    end
end