!_______________________________________________________________________________
!
!     Simulación de colas provocadas por interrupciones de tráfico tales que:
!
!       - Duración de cada interrupción:      T = EX(beta),
!                                             con media dmt=1/beta
!
!       - Tiempo entre llegadas de vehículos: W = Ex(dlam),
!                                             con media dmw=1/dlam
!
!     Para realizar la simulación se generan números aleatorios:
!
!       X  = U(0,1) -->  t=-1.d+00/beta*dlog(1.d+00-x) -> T = EX(beta)
!                        w=-1.d+00/dlam*dlog(1.d+00-x) -> W = EX(dlam)
!
!     El resultado de la simulación se compara con la solución teórica.
!
!     Se puede demostrar que el tamaño teórico de la cola que se forma es
!
!       C = G(p) - 1,                    <--- Geométrica trasladada
!       con p=1/(dmc+1), media dmc=dlam/beta y varianza dvc=(dmc+1)*dmc
!
!     Nota: El modelo no es muy realista, ya que
!             - puede haber interrupciones muy cortas,
!             - las llegadas no se ven afectadas por la longitud de la cola, y
!             - las colas desaparecen instantáneamente al abrirse el tramo.
!_______________________________________________________________________________

      implicit real*8 (a-h,o-z)

    1 write(6,101)
  101 format(/,' Duracion media de cada interrupcion (en minutos): ',$)
        read(5,*) dmt
        beta=1.d+00/dmt
        if (beta.le.0.d+00) goto 1

    2 write(6,102)
  102 format(/,' Intensidad media del trafico (vehiculos/minuto) : ',$)
        read(5,*) dlam
        if (dlam.le.0.d+00) goto 2

    3 write(6,103)
  103 format(/,' Numero de interrupciones a simular              : ',$)
        read(5,*) nintr
        if (nintr.lt.1) goto 3

    4 write(6,104)
  104 format(/,' Semilla para la simulacion (entero positivo)    : ',$)
        read(5,*) iseed
        if (iseed.le.0) goto 4

      write(6,200)
  200 format(/,' COMIENZO DE SIMULACION',/)
   
      write(6,201)
  201 format(/,'        Interrupcion',
     .         '    Duracion  mm: ss',
     .         '   Cola de Vehiculos',/)

      s_mc=0.d+00  ! Inicialización de la media muestral
      s_vc=0.d+00  ! Inicialización de la varianza muestral

      call srand(iseed)

      do iintr=1,nintr                          ! Generación de interrupciones

        x=dble(rand())
        t=-1.d+00/beta*dlog(1.d+00-x)           ! Duración de la interrupción

        nc=0
        x=dble(rand())
        w=-1.d+00/dlam*dlog(1.d+00-x)           ! Llegada del primer vehiculo
        wa=w 
        do while (wa.lt.t)
          nc=nc+1                               ! El vehiculo se pone en cola
          x=dble(rand())
          w=-1.d+00/dlam*dlog(1.d+00-x)         ! Llegada del siguiente vehiculo
          wa=wa+w
        enddo          

        ktm=int(t)                              ! Minutos
        kts=int((t-dble(ktm))*60.d+00)          ! Segundos

        write(6,202) iintr,ktm,kts,nc
        
  202   format(5x,i15,
     .         5x,i10,' : ',i2.2,
     .         5x,i15)

        s_mc=s_mc+dble(nc)
        s_vc=s_vc+dble(nc)**2

      enddo

      s_mc=s_mc/dble(nintr)                     ! Media muestral
      s_vc=s_vc/dble(nintr)-s_mc**2             ! Varianza muestral
      s_dc=sqrt(s_vc)                           ! Desviación muestral

      t_mc=dlam/beta                            ! Media teórica
      t_vc=(t_mc+1.d+00)*t_mc                   ! Varianza teórica
      t_dc=sqrt(t_vc)                           ! Desviación teórica

      write(6,203)
  203 format(/,' FIN DE SIMULACION',/)

      write(6,204) s_mc,t_mc,(1.d+00-s_mc/t_mc)*100.d+00
  204 format(/,' Media de la cola de vehiculos',
     .      //,'        Valor obtenido = ',f12.2,
     .       /,'        Valor teorico  = ',f12.2,
     .       /,'        Error          = ',f12.2,' %',/)
 
      write(6,205) s_dc,t_dc,(1.d+00-s_dc/t_dc)*100.d+00
  205 format(/,' Desviacion tipica de la cola de vehiculos',
     .      //,'        Valor obtenido = ',f12.2,
     .       /,'        Valor teorico  = ',f12.2,
     .       /,'        Error          = ',f12.2,' %',/)


      read(5,'()')
      call exit(0)
      end
