import numpy as np
def erlangB(A, k):

    invB = 1.0
    for m in xrange(1,k+1):
        invB = 1.0 + (m/A)*invB
    return 1.0/invB

def erlangC(A, k):
    """ Erlang C formula
        Also known as the Erlang delay formula

        Parameters
        ----------
        A : float
            offered load. Should be less than `k`
        k : integer
            number of servers or "trunks"

        Returns
        -------
        float
    """
    return erlangB(A, k)/(1+(A/k)*(erlangB(A, k)-1))

import simpy                                                       #librería funcional para colas
from scipy.stats import expon
import numpy as np
from erlang import erlangC

global request_delayed                                              #solicitud retrasada
global request_count                                                #contador de solicitudes 

def source_packets(env, lambda_p, server, capacity, mu_p):          #Define la función paquetes destino
    """Source generates customers randomly"""
    # Generates packets during simulation
    i=1
    global request_delayed
    global request_count
    while True:                                                     #Bucle infinito
        # Packet arrival time
        delta_arrival = expon.rvs(loc=0, scale= 1.0 / lambda_p )    #Tiempo de llegada expon.rvs -> variable aleatoria exponencial
        #print('Arrival Time: %6.3f' % (env.now + delta_arrival) )
        yield env.timeout(delta_arrival)                            #Cuando el tiempo se agota -> Retorna
        # Packet arrives: Create packet.
        c = packet(env, 'Packet%02d' % i, server, mu_p)
        #print('%7.4f Queue Size %02d' % (env.now, queue_length))
        request_count += 1                                          #Contador aumenta 1 paquete más en los servidores
        if server.count == capacity:                                #Si los servidores están llenos
                request_delayed += 1                                #Contador aumenta 1 paquete retrasado en cola
        env.process(c)
        i+=1 # Next packet counter                                  #Nuevo paquete creado

def ENQ(lambda_p, server, mu_p):
    """funcion para hallar E(NQ)"""
    Ro=lambda_p/(server*mu_p)
    Pq=erlangC(lambda_p/mu_p, server)
    yield Pq*Ro/(1-Ro)
    
def ETQ(lambda_p, incog):
    """función para hallar E(TQ)"""
    yield incog/lambda_p
    
def ET(mu_p, incog):
    """función para hallar E(T)"""
    yield incog+(1/mu_p)
    
def EN(lambda_p, incog):
    """función para hallar E(N)"""
    yield incog*lambda_p
        
        
def packet(env, name, server, mu_p):                                #Metodo de manejo de paquetes en servidor con tasa de atención U
    """Customer arrives, is served and leaves."""
    arrive = env.now
    #print('%7.4f %s: Arrival event' % (arrive, name))
    with server.request() as req:                                   #Cuando la función arroje req
        yield req
        delta_service = expon.rvs(loc=0, scale= 1.0 / mu_p )        #Función exponencial
        yield env.timeout(delta_service)                            #Tiempo agotado
        # Resource released automatically
        #print('%7.4f %s: Departure event' % (env.now, name))
        
def queue(env, NUMBER_SERVERS, MU_P, LAMBDA_P ):                    #Función de cola
    server = simpy.Resource(env, capacity=NUMBER_SERVERS)           #Capacidad de servidores -> Número de servidores
    # Start processes
    env.process(source_packets(env, LAMBDA_P, server, NUMBER_SERVERS, MU_P  )) #Proceso de cola con parametros lambda, Mu, K

if __name__=='__main__':

    global request_delayed
    global request_count
    
    LAMBDA_P = 180.0                                              #Lambda = 180
    MU_P = 15.0                                                   #Mu = 15
    k = 15                                                        #K=15 -> Cantidad de servidores
    A=LAMBDA_P/MU_P                                               #A -> Eficiencia del 1 servidor
    
    Enq=ENQ(LAMBDA_p,k,MU_P)
    Etq=ETQ(LAMBDA_P,Enq)
    Et=ET(MU_P,Etq)
    En=EN(LAMBDA_P,Et)
    
    conteo=0
    for i in range(10):
        TERMINATE_TIME = 100.0                                        #Tiempo de simulación
        print('MMK simulation.')
        # Setup and start the simulation
        request_delayed=0                                             #inicialización de Delays y Contadores
        request_count=0
        env = simpy.Environment()
        queue(env, k, MU_P, LAMBDA_P)                                 #Llamada a la función cola con parametros Mu, Lambda, k
        # Run
        env.run(until=TERMINATE_TIME)                                 #Inicio del programa hasta 100
        # Analytical results
        rho = LAMBDA_P/(k*MU_P)                                       #Eficiencia del sistema
        print 'Valor Teórico para E(NQ) = %g' %Enq
        print 'Valor Teórico para E(TQ) = %g' %Etq
        print 'Valor Teórico para E(T) = %g' %Et
        print 'Valor Teórico para E(N) = %g' %En
        print 'Total calls = %g, Delayed calls = %g, Delay Probability = %g' %(request_count, request_delayed, float(request_delayed)/request_count)
        print 'Theoretical request delay = {}'.format(erlangC(A, k))
        conteo=conteo+1
    print 'Probabilidad de Delay después del décimo experimento = %g' %float(request_delayed/request_count)
    
    
""" 
    CONCLUSIONES 
    Punto 4:
    Para un porcentaje de Delay de 10% se necesitan minimo 2774 servidores
    Para un porcentaje de Delay de 1% se necesitan mínimo 2825 servidores
    
    Al realizar una cantidad de iterciones sobre el mismo código continuo
    podemos concluir que la probabilidad de que se genere retraso tiende a cero.
    
    
"""
        

expon.rvs(loc=0, scale= 1.0 / 180.0)