"""Shared line vs. separate lines: a paired discrete-event simulation.

Two configurations run over IDENTICAL customer streams (same arrival times,
same per-customer service times in both setups):

  A) SHARED  : one FIFO line; the first free counter takes the head customer.
  B) SEPARATE: three dedicated lines; a new arrival joins the line with the
     fewest visible customers (the in-service customer counts as occupying
     its line), ties broken by earliest head completion, then uniformly at
     random; customers never switch lines (no jockeying).

Model assumptions
  - Arrivals: Poisson process (exponential interarrival gaps).
  - Service: exponential, mean 6 min (CV = 1). A customer's service time is
    drawn once and used in both setups.
  - Counters identical; FIFO service; no balking, reneging, or jockeying;
    unlimited queue space; service time does not depend on the counter.
"""

import random, math, json, statistics as st
from collections import deque

SERVICE_MEAN = 6.0        # minutes per customer
RHO          = 0.75       # per-counter offered load lam * SERVICE_MEAN / 3
DAY          = 480.0      # 8-hour operating day (minutes)
N_DAYS       = 30         # continuous simulated trading time = 30 shop-days
SEED_PRIMARY = 20260923
TIE_SEED     = 7
K_REPS       = 200        # paired replication count
REPS_DAYS    = 8          # days per replication


def gen_stream(T, lam, rng):
    """Poisson arrivals in [0, T) plus one exponential service time each."""
    ts = []; t = 0.0
    while True:
        t += rng.expovariate(lam)
        if t >= T: break
        ts.append(t)
    return ts, [rng.expovariate(1.0 / SERVICE_MEAN) for _ in ts]


def sim_shared(arrivals, services, c=3):
    """One shared FIFO line, c counters. O(n) with a min-scan."""
    free = [0.0] * c
    waits = []
    for a, s in zip(arrivals, services):
        j = min(range(c), key=free.__getitem__)
        start = max(a, free[j])
        waits.append(start - a)
        free[j] = start + s
    return waits


def sim_separate(arrivals, services, c=3, tie_rng=None):
    """c dedicated FIFO lines; arrivals join the shortest line and stay.

    Finished customers leave the line before the choice is made, so 'shortest'
    counts only customers still physically queued. Lindley recursion per line
    makes each line exactly FCFS.
    """
    lines = [deque() for _ in range(c)]
    waits = []
    for a, s in zip(arrivals, services):
        counts = []; heads = []
        for d in lines:
            while d and d[0] <= a:          # finished customers left the line
                d.popleft()
            counts.append(len(d))
            heads.append(d[0] if d else 0.0)
        m = min(counts)
        cand = [l for l in range(c) if counts[l] == m]
        if len(cand) > 1:
            h = min(heads[l] for l in cand)
            cand = [l for l in cand if heads[l] == h]
            l = tie_rng.choice(cand) if (tie_rng is not None and len(cand) > 1) else cand[0]
        else:
            l = cand[0]
        d = lines[l]
        start = max(a, d[-1] if d else 0.0)  # Lindley recursion, FCFS per line
        waits.append(start - a)
        d.append(start + s)
    return waits


mean = lambda w: sum(w) / len(w)

def p90(w):
    """The 90th percentile (nearest-rank): the wait exceeded by the slowest 10%."""
    s = sorted(w)
    return s[min(len(s) - 1, math.ceil(0.9 * len(s)) - 1)]

def top_decile_mean(w):
    """Mean wait experienced by the slowest 10% of customers."""
    s = sorted(w)
    k = max(1, math.ceil(0.1 * len(s)))
    return mean(s[-k:])


def erlang_c(lam, mu, c):
    """Erlang-C (probability of waiting) for the M/M/c shared-line check."""
    a = lam / mu
    rho = lam / (c * mu)
    s_ = sum(a**n / math.factorial(n) for n in range(c))
    tail = a**c / (math.factorial(c) * (1 - rho))
    return tail / (s_ + tail)


# --- internal validation: at c=1 both disciplines must be numerically identical
a0, s0 = gen_stream(500.0, 0.25, random.Random(99))
assert sim_shared(a0, s0, 1) == sim_separate(a0, s0, 1)
print("validation: single-counter equivalence of the two simulators PASSED")

# ---------------- primary run: one long identical stream ----------------
mu = 1.0 / SERVICE_MEAN
lam0 = RHO * 3 / SERVICE_MEAN          # arrivals/min so that rho = 0.75
arr, serv = gen_stream(DAY * N_DAYS, lam0, random.Random(SEED_PRIMARY))
w_shared = sim_shared(arr, serv)
w_sep    = sim_separate(arr, serv, 3, random.Random(TIE_SEED))
assert min(w_shared) >= 0 and min(w_sep) >= 0
n = len(arr)
print(f"primary stream: {n} customers over {DAY*N_DAYS:.0f} min "
      f"(lam={lam0:.4f}/min, mean service {SERVICE_MEAN} min, rho={RHO})")

print(f"\n{'setup':20s}{'mean wait':>11s}{'p90 wait':>10s}"
      f"{'mean of slowest 10%':>22s}{'zero-wait %':>13s}")
pri = {}
for name, w in (("shared", w_shared), ("separate", w_sep)):
    pri[name] = {"mean": mean(w), "p90": p90(w),
                 "top10_mean": top_decile_mean(w),
                 "zero_share": sum(1 for x in w if x <= 1e-9) / n}
C = erlang_c(lam0, mu, 3)
th_shared_mean = C / (3 * mu - lam0)
th_shared_p90  = math.log(C / 0.1) / (3 * mu - lam0)
lam_line = lam0 / 3
rho_line = lam_line / mu
th_sep_mean = lam_line / (mu * (mu - lam_line))          # random-assignment M/M/1 x3
th_sep_p90  = math.log(rho_line / 0.1) / (mu - lam_line)
print(f"{'Shared (one line)':20s}{mean(w_shared):>9.2f}m{p90(w_shared):>9.2f}m"
      f"{top_decile_mean(w_shared):>21.2f}m{100*pri['shared']['zero_share']:>12.1f}%")
print(f"{'Separate (3 lines)':20s}{mean(w_sep):>9.2f}m{p90(w_sep):>9.2f}m"
      f"{top_decile_mean(w_sep):>21.2f}m{100*pri['separate']['zero_share']:>12.1f}%")
print(f"M/M theory (random assignment, not shortest-line): shared mean "
      f"{th_shared_mean:.2f} / p90 {th_shared_p90:.2f}; separate mean "
      f"{th_sep_mean:.2f} / p90 {th_sep_p90:.2f} (min)")

# ---------------- sensitivity to traffic level ----------------
print("\nsensitivity: rho | n | shared mean/p90 | separate mean/p90 (min)")
sens_rows = []
for i, rho in enumerate((0.5, 0.65, 0.75, 0.9)):
    lam = rho * 3 / SERVICE_MEAN
    rng = random.Random(4242 + i)
    ar, sv = gen_stream(DAY * N_DAYS, lam, rng)
    w1 = sim_shared(ar, sv)
    w2 = sim_separate(ar, sv, 3, random.Random(TIE_SEED))
    row = {"rho": rho, "n": len(ar),
           "mean_shared": mean(w1), "mean_sep": mean(w2),
           "p90_shared": p90(w1), "p90_sep": p90(w2)}
    sens_rows.append(row)
    print(f"   {rho:.2f} | {row['n']:5d} | {row['mean_shared']:7.2f} / {row['p90_shared']:6.2f}"
          f" | {row['mean_sep']:8.2f} / {row['p90_sep']:8.2f}")

# ---------------- paired replications for uncertainty ----------------
dm = []; dp90 = []; wins = 0
for i in range(K_REPS):
    rng = random.Random(1000 + i)
    ar, sv = gen_stream(DAY * REPS_DAYS, lam0, rng)
    w1 = sim_shared(ar, sv)
    w2 = sim_separate(ar, sv, 3, random.Random(TIE_SEED))
    m1, m2 = mean(w1), mean(w2)
    dm.append(m2 - m1)
    dp90.append(p90(w2) - p90(w1))
    wins += (m1 < m2)
rep = {"k": K_REPS, "mean_gap": st.mean(dm), "sd_gap": st.stdev(dm),
       "min_gap": min(dm), "max_gap": max(dm), "p90_gap_mean": st.mean(dp90),
       "shared_better": wins / K_REPS}
print(f"\npaired replications (n={K_REPS} independent 8-day streams, "
      f"same stream for both setups):")
print(f"mean-wait gap (separate minus shared) = {rep['mean_gap']:.2f} +/- "
      f"{rep['sd_gap']:.2f} min (range {rep['min_gap']:.2f} .. {rep['max_gap']:.2f}); "
      f"shared faster in {wins}/{K_REPS} streams")

# ---------------- save results ----------------
out = {
    "meta": {"run_date": "2026-09-23",
             "model": "discrete-event FCFS simulation; identical streams for both setups"},
    "assumptions": [
        "Poisson arrivals (exponential interarrival times, lam=0.375/min for primary run)",
        "Exponential service times, mean 6.0 min (CV=1); a customer's service time is "
        "fixed and identical in both setups",
        "Counters identical; FIFO service; no balking, reneging, or jockeying; "
        "unlimited queue space",
        "In the separate-lines setup a customer counts the in-service person as "
        "occupying that line",
        "Ties in line length broken by earliest head completion, then uniformly at random",
    ],
    "parameters": {"service_mean_min": SERVICE_MEAN, "rho_primary": RHO,
                   "lam_primary_per_min": lam0, "day_minutes": DAY,
                   "days_primary": N_DAYS, "primary_seed": SEED_PRIMARY,
                   "tie_seed": TIE_SEED, "replications": K_REPS,
                   "replication_days": REPS_DAYS},
    "primary": {"n_customers": n,
                "shared": pri["shared"], "separate": pri["separate"],
                "shared_m_mean_ratio": mean(w_sep) / mean(w_shared),
                "shared_p90_ratio": p90(w_sep) / p90(w_shared)},
    "theoretical_mm_check": {"shared_mean": th_shared_mean, "shared_p90": th_shared_p90,
                             "separate_mean": th_sep_mean, "separate_p90": th_sep_p90},
    "sensitivity": sens_rows,
    "replications": rep,
}
with open("results.json", "w") as f:
    json.dump(out, f, indent=1)
print("\nsaved results.json")
