Little's Law Simulation (L = λW)

The picture behind L = λW

Count arrivals up to time t as A(t) and departures as D(t). The number in the system at any moment is the vertical gap A(t) − D(t). In the lab each customer is drawn as a horizontal strip one unit high that starts when the customer arrives and ends when the customer leaves, so the strips exactly fill the region between the two curves. The area can be added up two ways:

• Column by column: integrating the gap over time gives ∫(A − D) dt = T · L, where L is the time-average number in the system.
• Strip by strip: each strip’s length is one customer’s time in the system, so the same area is N · W for N customers with mean time W.

That is the whole argument when the system starts and ends empty. Press Play to slide the window through the run: the strips show waiting (amber) and service (green) separately, and the lower chart shows the running L(t) with λ̂(t)·W(t) drawn on top of it, the two agreeing at every point.

What the law needs, and what it does not

John Little proved the result in 1961 for stationary systems, and the sample-path argument above extends it much further. It makes no assumption about the arrival process, the service-time distribution, the number of servers or the order of service. Switch the service distribution in the lab: the ratio L/(λ̂W) stays at 1 for all of them, while the waits themselves change a lot.

It does need three things:

• Averages that settle. L, λ and W must have finite long-run values. With λ ≥ μ the queue grows without bound; the identity still holds over any finished window, but there is no steady-state L or W for it to relate.
• One consistent boundary. Count L, λ and W for the same box. Around the waiting line only, it gives L_q = λW_q; around the server only, the mean number in service equals λ × mean service time, which is the utilization ρ.
• Care at the edges. A finite window cuts through the strips of customers still inside at T. The lab averages W over completed customers only, so the ratio differs from 1 by that end effect; counting the partial time of those still inside makes it exact. The error shrinks as T grows.

In operations the same law reads work in process = throughput × cycle time, and it is used in reverse: measure two of the three and the third follows. The M/M/1 queue and M/M/c queue pages use it to connect their formulas for L and W.

A second check: the Pollaczek–Khinchine formula

Little’s law relates averages but does not predict them. For a single server with Poisson arrivals and any service distribution (the M/G/1 queue), the Pollaczek–Khinchine formula does, from the first two moments of the service time S:

All four service options have mean 1/μ; they differ in variability:

Same arrival rate, same mean service time, same utilization, and the mean queue differs by a factor of five: variability, not just load, creates waiting. The hyperexponential option splits the mean equally between a frequent short phase and a rare long one, which gives a squared coefficient of variation of 4. Its runs converge slowly, because rare long services dominate the averages.

Example result

With λ = 0.8, μ = 1, exponential service and T = 5,000, seed 42 produced 4,007 arrivals and 4,005 departures. The time-average number in system is L = 3.513, completed customers spent W = 4.385 on average, and with λ̂ = 0.8014 the ratio L/(λ̂W) is 0.9997. The Pollaczek–Khinchine steady-state value is 4.000: a single run of this length can miss the steady state by a wide margin, yet Little’s law holds on the run itself.

Implementation

The page runs simulation.js from the package below; the Python file is a line-by-line port with the same random-number generator, so a given seed produces the same result in both languages. The excerpts show the core of each.

// A single FIFO server needs no event list: customer k starts service when it arrives
// or when customer k − 1 leaves, whichever is later (Lindley's recursion).
// Prefix sums of the waits and sojourns make every statistic on [0, t] an O(log n) lookup.
export function simulate(p, rng) {
  const a = [], start = [], d = [];
  const pa = [0], pq = [0], pw = [0]; // prefix sums of arrival times, waits, sojourns
  let t = 0, prev = 0;
  for (;;) {
    t += rng.exponential(p.lambda);
    if (t > p.horizon) break;
    const s = Math.max(t, prev);
    const dep = s + drawService(p.service, p.mu, rng);
    a.push(t); start.push(s); d.push(dep);
    pa.push(pa[pa.length - 1] + t);
    pq.push(pq[pq.length - 1] + (s - t));
    pw.push(pw[pw.length - 1] + (dep - t));
    prev = dep;
  }
  return { p, a, start, d, pa, pq, pw };
}

// Everything observed on [0, t]. Area between A(u) and D(u) = time in system of the
// customers who left + time so far of those still inside.
export function statsAt(r, t) {
  const A = countUpTo(r.a, t), B = countUpTo(r.start, t), D = countUpTo(r.d, t);
  const area = r.pw[D] + (A - D) * t - (r.pa[A] - r.pa[D]);
  const areaQ = r.pq[B] + (A - B) * t - (r.pa[A] - r.pa[B]);
  const L = area / t, lambdaHat = A / t, W = D > 0 ? r.pw[D] / D : null;
  return { t, A, D, L, Lq: areaQ / t, lambdaHat, W, ratio: W !== null && W > 0 ? L / (lambdaHat * W) : null };
}

def simulate(p, rng):
    """Lindley's recursion: start = max(arrival, previous departure). Keeps prefix sums."""
    a, start, d = [], [], []
    pa, pq, pw = [0.0], [0.0], [0.0]  # prefix sums of arrival times, waits, sojourns
    t = prev = 0.0
    while True:
        t += rng.exponential(p["lambda"])
        if t > p["horizon"]:
            break
        s = max(t, prev)
        dep = s + draw_service(p["service"], p["mu"], rng)
        a.append(t)
        start.append(s)
        d.append(dep)
        pa.append(pa[-1] + t)
        pq.append(pq[-1] + (s - t))
        pw.append(pw[-1] + (dep - t))
        prev = dep
    return {"a": a, "start": start, "d": d, "pa": pa, "pq": pq, "pw": pw}


def stats_at(r, t):
    """L, Lq, lambda-hat and W observed on [0, t]; area between A(u) and D(u) = total time in system."""
    A, B, D = count_up_to(r["a"], t), count_up_to(r["start"], t), count_up_to(r["d"], t)
    area = r["pw"][D] + (A - D) * t - (r["pa"][A] - r["pa"][D])
    area_q = r["pq"][B] + (A - B) * t - (r["pa"][A] - r["pa"][B])
    L, lam_hat = area / t, A / t
    W = r["pw"][D] / D if D > 0 else None
    ratio = L / (lam_hat * W) if W is not None and W > 0 else None
    return {"A": A, "D": D, "L": L, "Lq": area_q / t, "lambda_hat": lam_hat, "W": W, "ratio": ratio}

Run the full model from a terminal:

node simulation.js --lambda 0.8 --mu 1 --seed 42
python simulation.py --lambda 0.8 --mu 1 --seed 42

Download the model

Assumptions & limitations

• One server, first come, first served. Service start times follow Lindley’s recursion; Little’s law itself would hold for any discipline.
• Poisson arrivals, independent service times. Needed for the Pollaczek–Khinchine comparison, not for Little’s law.
• Starts empty at t = 0, observed on [0, T]. No warm-up is discarded, so the comparison with the steady-state prediction includes initialization bias; the Little’s law ratio is unaffected by it.
• Completed customers only in W. Customers still present at T contribute to L but not to W, which is the source of the small departure of the ratio from 1.

Related simulations

• OS-018M/M/1 queueOne server, random arrivals, random service: why waits explode as utilization ρ = λ/μ approaches 1.
• OS-019M/M/c queueSeveral servers, one shared line: Erlang C waiting probability and the pooling effect.
• OS-021Discrete-event simulationClock, future event list, event routines: the method behind queue simulators, one event at a time.
• OS-015Poisson processEvents at random times with a constant rate: Poisson counts, exponential gaps and the waiting-time paradox.

Background: SimPy Examples: Queues in Python

References

1. Little, J. D. C. (1961). A proof for the queuing formula: L = λW. Operations Research, 9(3), 383–387. doi.org/10.1287/opre.9.3.383
2. Kleinrock, L. (1975). Queueing Systems, Volume 1: Theory. Wiley.
3. Law, A. M. (2015). Simulation Modeling and Analysis, 5th ed. McGraw-Hill.