#!/usr/bin/perl -w

use 5.030;

say "";
say "Modelo de colas (M/M/R):(GD/K/K)";

my $K = 22;         # tamaño de la cola y la población
say "Máxima cola y población K = $K";

# my $R = 5;          # número de servidores o reparadores
my $lambda = 0.5;   # tasa de llegadas
say "Tasa de llegadas o averías, lambda = $lambda";
my $mu = 5.0;       # tasa de servicio
say "Tasa de servico o reparación,   mu = $mu";
my $po = $lambda/$mu; # ratio que debe ser < 1
say "";

my (@lambda, @mu, @p);

@lambda = (0) x 1000;

my $R = 5;         # ejemplo

for my $i (0 .. $K){
    $lambda[$i] = ($K - $i) * $lambda;
    if ($i <= $R){
	$mu[$i] = $i * $mu;
    }
    if ($i > $R && $i <= $K){
	$mu[$i] = $R * $mu;
    }elsif ($i > $K) {
	$mu[$i] = 0;
    }
}

format STDOUT =
@### @###.### @###.### @###.### @###.### @####.#### @#.####### @####.#### @####.#### @####.#### @####.#### @#.#######
$R , $lambda,$lambda[$R],$mu,$mu[$R],l_eff($R), p(0) ,   L($R)  ,  Lq($R)  ,  Ws($R),  Wq($R)  , utilization($R)
.
  

say "R      lambda lambda(R)   mu      mu(R)   lambda_eff(R)   p(0)      Ls(R)       Lq         Ws        Wq     utilization"; 
for my $i (1 .. $K){
    $R = $i;
#    say "$i $lambda $lambda[$i] $mu $mu[$i]  ", l_eff($i),"  ", p(0),"  ", L($i),"  ", Lq($i),"  ", Ws($i),"  ", Wq($i),"  ", utilization($i), "";
    write();
}

exit 1;


sub p {
    my $n = shift;
    
    if ($n == 0){
	my $sum;
	for my $i (0 .. $R){
	    $sum += comb($K, $i) * ($po**$i);
	}
	for my $i ($R+1 .. $K){
	    $sum += comb($K, $i) * fact($i) * ($po**$i) / (fact($R) * $R**($i-$R));
	}
	return $sum**(-1);
    }
    
    if ($n > 0 && $n <= $R){
	return comb($K, $n) * ($po**$n) * p(0);
    }
    
    if ($n >= $R && $n <= $K){
	return comb($K, $n) * fact($n) * ($po**$n) * p(0) / (fact($R) * $R**($n-$R));
    }
}

sub L {
    my $s = shift;
    my $sum;
    for my $i (0 .. $K){
	$sum += $i * p($i);
    }
    return $sum;
}

sub l_eff {
    my $s = shift;
    return $lambda * ($K - L($s));
}

sub Ws {
    my $s = shift;
    return L($s) / l_eff($s);
}

sub Lq {
    my $s = shift;
    return L($s) - l_eff($s) / $mu;
}

sub Wq {
    my $s = shift;
    return Lq($s) / l_eff($s);
}

sub utilization {
    my $s = shift;
    return (L($s) - Lq($s)) / $s;
}


sub fact {
    my $i = shift;
    return 1 if $i <= 1;
    return $i*fact($i-1);
}

sub comb {
    my ($m, $n) = @_;
    return fact($m)/(fact($n)*fact($m-$n));
}

__END__