Trabajo final

1 Ejercicio 1: Soldadura de barras (500 simulaciones)

Enunciado: Una barra de longitud \(x_1\) (\(\mathcal{N}(30, 0.81)\)) se une a \(x_2\) (\(\mathcal{N}(18, 0.3^2)\)). La longitud soldada es \(L = x_1 + x_2\). Las especificaciones son 50 ± 5 cm (LSL = 45, USL = 55).

Objetivos: simular 500 barras, estimar \(\hat p = P(L \notin [45,55])\), y calcular Cp y Cpk. Además, comparar contra los valores analíticos.

# Parámetros
n   <- 500
mu1 <- 30;  sd1 <- sqrt(0.81)    # 0.9
mu2 <- 18;  sd2 <- 0.3
LSL <- 45;  USL <- 55

# Simulación
x1 <- rnorm(n, mean = mu1, sd = sd1)
x2 <- rnorm(n, mean = mu2, sd = sd2)
L  <- x1 + x2

# Estimador de la probabilidad fuera de especificación
p_hat <- mean(L < LSL | L > USL)

# Índices de capacidad (muestrales)
mu_hat  <- mean(L)
s_hat   <- sd(L)
Cp_hat  <- (USL - LSL) / (6 * s_hat)
Cpk_hat <- min((USL - mu_hat) / (3 * s_hat),
               (mu_hat - LSL) / (3 * s_hat))

list(
  media_muestral_L = round(mu_hat, 3),
  sd_muestral_L    = round(s_hat, 3),
  p_hat_fuera_esp  = round(p_hat, 5),
  Cp_hat           = round(Cp_hat, 3),
  Cpk_hat          = round(Cpk_hat, 3)
)

## $media_muestral_L
## [1] 47.994
## 
## $sd_muestral_L
## [1] 1.029
## 
## $p_hat_fuera_esp
## [1] 0
## 
## $Cp_hat
## [1] 1.619
## 
## $Cpk_hat
## [1] 0.97

# Valores analíticos para contraste:
mu_L <- mu1 + mu2
sd_L <- sqrt(sd1^2 + sd2^2)

p_exact <- pnorm(LSL, mean = mu_L, sd = sd_L) +
  (1 - pnorm(USL, mean = mu_L, sd = sd_L))

Cp_th  <- (USL - LSL) / (6 * sd_L)
Cpk_th <- min((USL - mu_L) / (3 * sd_L), (mu_L - LSL) / (3 * sd_L))

list(
  media_teorica_L  = round(mu_L, 3),
  sd_teorica_L     = round(sd_L, 3),
  p_exact_fuera_esp= signif(p_exact, 6),
  Cp_teorico       = round(Cp_th, 3),
  Cpk_teorico      = round(Cpk_th, 3)
)

## $media_teorica_L
## [1] 48
## 
## $sd_teorica_L
## [1] 0.949
## 
## $p_exact_fuera_esp
## [1] 0.000782701
## 
## $Cp_teorico
## [1] 1.757
## 
## $Cpk_teorico
## [1] 1.054

hist(L, breaks = 25,
     main = "Longitud soldada L (500 simulaciones)",
     xlab = "cm",
     col = "#90caf9", border = "#1e88e5")
abline(v = c(LSL, 50, USL), lty = c(2, 1, 2), col = c("#e53935", "#43a047", "#e53935"), lwd = c(2, 2, 2))
legend("topright",
       legend = c("LSL = 45", "Centro = 50", "USL = 55"),
       lty = c(2, 1, 2),
       col = c("#e53935", "#43a047", "#e53935"),
       bty = "n", lwd = 2)
rug(L, col = "#0d47a1")

Notas
- Cp mide capacidad potencial (dispersión vs. ancho de tolerancia).
- Cpk incorpora el descentramiento (qué tan centrado está el proceso).
- Para afirmar control estadístico se requieren cartas de control con datos en el tiempo (X̄–R/X̄–S, etc.). Cp/Cpk por sí mismos no prueban control.

2 Ejercicio 2: Reproceso de piezas (200 piezas, 10 réplicas)

Sistema: Llegadas determinísticas cada 20 min. Por pieza, número de defectos \(d\sim \text{Binomial}(3, 0.8)\). El tiempo de reparación es la suma de d exponenciales con tasa 0.2 (equivalente a \(\text{Gamma}(d, 0.2)\)). Un solo servidor (FIFO).

Objetivo: simular 10 réplicas y reportar el tiempo total para completar 200 piezas, con IC 95%.

tiempo_total_200 <- function(n_piezas = 200, ia = 20, rate = 0.2) {
  # Llegadas determinísticas
  llegada <- seq(0, by = ia, length.out = n_piezas)
  
  # Número de defectos por pieza
  d <- rbinom(n_piezas, size = 3, prob = 0.8)
  
  # Tiempos de servicio (suma de d exponentiales; 0 si d=0)
  serv <- numeric(n_piezas)
  for (i in seq_len(n_piezas)) {
    if (d[i] > 0) serv[i] <- sum(rexp(d[i], rate = rate))
  }
  
  # Dinámica de cola de un servidor (FIFO)
  ini <- fin <- numeric(n_piezas)
  for (i in seq_len(n_piezas)) {
    if (i == 1) {
      ini[i] <- llegada[i]
    } else {
      ini[i] <- max(llegada[i], fin[i - 1])
    }
    fin[i] <- ini[i] + serv[i]
  }
  fin[n_piezas]
}

R <- 10
tiempos <- replicate(R, tiempo_total_200())
tiempos

##  [1] 3990.693 3986.006 4007.530 3993.484 3986.172 3997.947 3988.592 3985.322
##  [9] 4015.174 3989.471

media <- mean(tiempos)
sdres <- sd(tiempos)
tcrit <- qt(0.975, df = R - 1)
IC_95 <- c(media - tcrit * sdres / sqrt(R),
           media + tcrit * sdres / sqrt(R))

resumen <- data.frame(
  Replica = 1:R,
  Tiempo_min = round(tiempos, 2)
)
resumen

list(
  media_min = round(media, 2),
  sd_min    = round(sdres, 2),
  IC95_inf  = round(IC_95[1], 2),
  IC95_sup  = round(IC_95[2], 2)
)

## $media_min
## [1] 3994.04
## 
## $sd_min
## [1] 10.05
## 
## $IC95_inf
## [1] 3986.85
## 
## $IC95_sup
## [1] 4001.23

Notas
- El último arribo es en \((200-1)\times 20 = 3980\) min.
- El total suele estar cerca de ~4000 min (depende de los tiempos de servicio y de espera de las últimas piezas).

3 Ejercicio 3: Camión de reparto (1 andén de carga)

Un camión tarda 30 ± 10 min en ser cargado (Uniforme(20, 40)), 20 ± 5 min en descargar (Uniforme(15, 25)) y 40 min con distribución exponencial en trasladarse (tanto ida como regreso), es decir, tiempos ~ Exp(rate = 1/40).

3.1 (a) Simular 10 horas con 5 réplicas

# Simulación de un camión por 10 horas (600 min)
viajes_dia_1camion <- function(dia = 600) {
  t <- 0; viajes <- 0
  repeat {
    # Carga
    t <- t + runif(1, 20, 40)
    # Traslado ida
    t <- t + rexp(1, rate = 1/40)
    # Descarga
    t <- t + runif(1, 15, 25)
    # Traslado regreso
    t <- t + rexp(1, rate = 1/40)
    if (t <= dia) viajes <- viajes + 1 else break
  }
  viajes
}

set.seed(525)
R3 <- 5
viajes <- replicate(R3, viajes_dia_1camion())
viajes

## [1] 4 3 4 5 4

barplot(viajes,
        names.arg = paste("Rep", 1:R3),
        ylim = c(0, max(viajes) + 2),
        main = "Número de viajes por réplica (10 horas, 1 camión)",
        xlab = "Réplica", ylab = "Viajes",
        col = "#ffcc80", border = "#ef6c00")
abline(h = 10, col = "#d32f2f", lty = 2, lwd = 2)
legend("topleft", legend = c("Meta: 10 viajes"), lty = 2, col = "#d32f2f", bty = "n", lwd = 2)

3.2 (b) Intervalo de confianza para los viajes por día

media_v  <- mean(viajes)
sd_v     <- sd(viajes)
tcrit3   <- qt(0.975, df = R3 - 1)
IC95_v   <- c(media_v - tcrit3 * sd_v / sqrt(R3),
              media_v + tcrit3 * sd_v / sqrt(R3))

data.frame(
  media_viajes = round(media_v, 2),
  sd_viajes    = round(sd_v, 2),
  IC95_inf     = round(IC95_v[1], 2),
  IC95_sup     = round(IC95_v[2], 2)
)

3.3 (c) Recomendaciones para lograr ≥ 10 entregas/día con un solo andén

El tiempo de un ciclo esperado por camión es: \[ E[\text{ciclo}] = E[\text{carga}] + E[\text{descarga}] + 2\,E[\text{traslado}] = 30 + 20 + 2\cdot 40 = 130 \text{ min}. \] En 10 horas (600 min), un camión hace en promedio \(600/130 \approx 4.6\) viajes → insuficiente para 10.

Una opción viable es operar k camiones en rotación con un solo andén (cola FIFO para cargar). Abajo simulamos este sistema para distintos \(k\).

# Simulador discreto: k camiones, 1 andén, FIFO
viajes_dia_kcam <- function(k = 3, dia = 600) {
  t <- 0
  cola <- 1:k               # al inicio, todos esperan cargar
  anden_ocupado <- FALSE
  camion_en_carga <- NA
  t_fin_carga <- Inf
  t_regreso <- rep(Inf, k)  # tiempos de regreso de cada camión
  viajes <- 0
  
  # iniciar primera carga si procede
  if (length(cola) > 0) {
    anden_ocupado <- TRUE
    camion_en_carga <- cola[1]; cola <- cola[-1]
    t_fin_carga <- t + runif(1, 20, 40)
  }
  
  repeat {
    t_eventos <- c(t_fin_carga, t_regreso)
    idx <- which.min(t_eventos)
    t <- t_eventos[idx]
    if (!is.finite(t)) break
    
    if (idx == 1) {
      # Fin de carga -> camión sale a ruta
      anden_ocupado <- FALSE
      id <- camion_en_carga
      t_vuelta <- rexp(1, rate = 1/40) + runif(1, 15, 25) + rexp(1, rate = 1/40)
      t_regreso[id] <- t + t_vuelta
      t_fin_carga <- Inf
      # Atender siguiente de la cola si hay y no se acabó el día
      if (length(cola) > 0 && t <= dia) {
        anden_ocupado <- TRUE
        camion_en_carga <- cola[1]; cola <- cola[-1]
        t_fin_carga <- t + runif(1, 20, 40)
      }
    } else {
      # Regresa camión j
      j <- idx - 1
      if (t <= dia) viajes <- viajes + 1
      cola <- c(cola, j)
      t_regreso[j] <- Inf
      if (!anden_ocupado && t <= dia) {
        anden_ocupado <- TRUE
        camion_en_carga <- cola[1]; cola <- cola[-1]
        t_fin_carga <- t + runif(1, 20, 40)
      }
    }
  }
  viajes
}

set.seed(525)
k_vals <- 1:5
B <- 500  # número de días simulados por k
prom <- sapply(k_vals, function(k) mean(replicate(B, viajes_dia_kcam(k))))
prob10 <- sapply(k_vals, function(k) mean(replicate(B, viajes_dia_kcam(k)) >= 10))

tabla_k <- data.frame(
  k_camiones = k_vals,
  promedio_viajes = round(prom, 2),
  prob_cumplir_ge_10 = round(prob10, 3)
)
tabla_k

plot(k_vals, prom, type = "b", pch = 19,
     main = "Viajes promedio vs. número de camiones (1 andén)",
     xlab = "k camiones", ylab = "Viajes promedio en 10 horas",
     col = "#42a5f5")
abline(h = 10, col = "#c62828", lty = 2, lwd = 2)