Some jams have a cause you can point to: a crash, a closed lane, a merge. Many have none. When each driver reacts a little late, one tap on the brakes can grow into a wave of stopped cars that moves backward while every car in it moves forward. Three labs: a ring road where a jam appears from nothing, the curve that ties density, speed and flow together, and a lane closure where taking turns shortens the line.
Put cars on a circular road with no exits, no lane changes and nothing in the way. Each driver follows the car ahead using the Intelligent Driver Model (Treiber, Hennecke and Helbing, 2000): speed up toward a desired speed, keep a safe time gap, and brake harder the faster you close in. In 2008 Yuki Sugiyama and colleagues tried this with real drivers, 22 cars on a circle around, and a stop-and-go jam formed with no obstacle at all. Tap a car to make its driver brake for a moment.
How to read this: in the lower chart, time runs left to right and height is where a car is on the ring. Each car draws a rising diagonal line as it drives around. A jam shows up as a red band that slopes down: the jam moves backward even though every car moves forward.
Try this: pick Sluggish, press Brake a car, and watch a red clump form and drift backward around the ring. Time gap is how far back drivers stay, Pickup is how briskly they speed up, and Reaction delay is how late they respond.
v is your speed and v₀ the speed you want ( here). s is the gap to the car ahead and Δv how fast you are closing on it. s* is the gap you would like: a minimum s₀ of , plus T seconds of travel, plus extra if you are closing in. a is how briskly you speed up, b how hard you like to brake.
A jam is a place, and cars pass through it. Each driver at the front pulls away a moment after the car ahead does, so the front of the jam slips back one car at a time, while new cars arrive at its back. Every car moves forward, yet the jam drifts against traffic. On real highways these waves travel upstream at about 9 mph.
Whether a brake tap fades or grows depends on the drivers. With brisk pickup each driver brakes a little less than the one ahead, and the tap dies out. With slow pickup or a reaction delay each driver brakes a little more, and the tap grows until someone stops. Engineers call that string instability.
import math, random
L, N = 230.0, 22 # ring length (m), number of cars
LEN, S0 = 5.0, 2.0 # car length (m), minimum gap (m)
v0, T, a, b = 50/3.6, 1.0, 1.0, 1.5 # desired speed m/s, time gap s, accel and braking m/s^2
def idm(v, dv, s): # dv = my speed - leader's speed
s_star = S0 + max(0.0, v*T + v*dv/(2*math.sqrt(a*b)))
return a*(1 - (v/v0)**4 - (s_star/s)**2)
def v_eq(s): # speed at which idm(v, 0, s) = 0
lo, hi = 0.0, v0
for _ in range(60):
m = (lo + hi)/2
lo, hi = (m, hi) if idm(m, 0, s) > 0 else (lo, m)
return lo
x = [i*L/N + random.uniform(-0.15, 0.15) for i in range(N)] # unwrapped positions
v = [v_eq(L/N - LEN)]*N
dt = 0.05
for step in range(int(300/dt)):
t = step*dt
acc = []
for i in range(N):
j = (i + 1) % N # the car ahead
gap = x[j] + (L if j == 0 else 0) - x[i] - LEN
acc.append(idm(v[i], v[i] - v[j], gap))
if 3 <= t < 4.5:
acc[0] = min(acc[0], -4.0) # one driver taps the brakes
for i in range(N):
vn = v[i] + acc[i]*dt
if vn < 0: # stop, don't roll backward
x[i] += -v[i]**2/(2*acc[i]); v[i] = 0.0
else:
x[i] += v[i]*dt + acc[i]*dt*dt/2; v[i] = vn
print('slowest %.1f m/s, fastest %.1f m/s' % (min(v), max(v)))Traffic engineers describe a lane with three numbers. Density k is how many cars fill a stretch of it. Speed v is how fast they go. Flow q is how many pass a fixed point each hour. Flow is density times speed, so an empty road has speed but little flow and a packed road has cars but little speed. In between sits a critical density where flow peaks. Here the same driver model runs on a ring 0.62 mi long. Pick a density, or sweep them all and plot what the ring actually carries.
How to read this: the strip on top is the test road, with cars colored by speed as in Lab 1 (green cruising, amber slowing, red stopped). The chart below shows how many cars per hour the road carries at each level of crowding. The amber dot is the strip right now; the green dots are past test runs.
Try this: press Sweep densities, then switch to the lower speed limit and sweep again. Compare where the two sets of dots peak.
The gold curve assumes every car drives at the steady speed ve for its spacing, the speed at which the model's acceleration is zero. For this limit it peaks at —. Green dots come from running the ring for 15 minutes and averaging the last 5.
Past the critical density, the green dots fall below the gold curve. Steady flow is possible there in theory, but it is unstable: small speed differences grow into stop-and-go waves, and cars stopped in a wave carry no one anywhere.
A lower limit lowers the theoretical peak a little. It also narrows the speed differences that seed waves, so near capacity the stream stays smooth where the faster one breaks down. That is the idea behind variable speed limits on busy freeways. Far below capacity, the lower limit only slows everyone.
import math, random
LEN, S0 = 5.0, 2.0
T, a, b = 1.2, 0.7, 2.0 # the drivers in this lab
def idm(v, dv, s, v0):
s_star = S0 + max(0.0, v*T + v*dv/(2*math.sqrt(a*b)))
return a*(1 - (v/v0)**4 - (s_star/s)**2)
def v_eq(s, v0):
lo, hi = 0.0, v0
for _ in range(60):
m = (lo + hi)/2
lo, hi = (m, hi) if idm(m, 0, s, v0) > 0 else (lo, m)
return lo
def theory(k, v0): # k in cars per km, v0 in m/s
v = v_eq(1000/k - LEN, v0)
return k*v*3.6 # flow in cars per hour
def measured(k, v0, L=1000.0, dt=0.1): # run the ring and time-average
N = round(k*L/1000)
x = [i*L/N + random.uniform(-0.15, 0.15) for i in range(N)]
v = [v_eq(L/N - LEN, v0)]*N
total, n = 0.0, 0
for step in range(9000): # 900 s; measure the last 300 s
acc = [idm(v[i], v[i] - v[(i+1) % N],
x[(i+1) % N] + (L if i == N-1 else 0) - x[i] - LEN, v0) for i in range(N)]
for i in range(N):
vn = v[i] + acc[i]*dt
if vn < 0: x[i] -= v[i]**2/(2*acc[i]); v[i] = 0.0
else: x[i] += v[i]*dt + acc[i]*dt*dt/2; v[i] = vn
if step >= 6000 and step % 10 == 0:
total += sum(v)/N; n += 1
return k*(total/n)*3.6
for v0 in (120/3.6, 80/3.6):
kc = max(range(5, 140), key=lambda k: theory(k, v0))
print('limit %3.0f km/h: peak %4.0f cars/h at %d cars/km' % (v0*3.6, theory(kc, v0), kc))
print(' at 40 cars/km: theory %4.0f, measured %4.0f' % (theory(40, v0), measured(40, v0)))A lane closure squeezes two lanes into one. Drivers can move over early, as soon as they see the first sign, or use both lanes up to the merge point and take turns there like the teeth of a zipper. This idealized model is a Nagel–Schreckenberg cellular automaton (1992). The road is a row of cells, one car per cell, and cars move up to 3 cells per second (about ). Each second every car speeds up by one, slows to fit its gap, hesitates with 15% chance (a 006 · Probability coin flip), then moves. Both roads below receive exactly the same arriving cars.
How to read this: two copies of the same road get exactly the same cars. Top: everyone moves over at the first sign. Bottom: both lanes are used up to the gold merge line, then drivers take turns. The chart below tracks how long each line gets.
Try this: drag Arriving cars low, press Restart, then try a high value. Watch when the amber line in the chart pulls away from the green one.
In the early version the closed lane sits empty, so every waiting car stands in one file and the line is about twice as long for the same number of cars. In the zipper version two lines feed the merge point. When one driver is slow to pull away, the car from the other lane fills the space, so the single lane past the merge stays fuller. In this model that lifts the cars served per hour as well as shortening the line.
Neither strategy beats the one-lane limit. Below it, nothing queues and the two roads behave the same. Set the arrivals under about 1,500 per hour to see that.
Headway is the time between one car's front bumper passing a point and the next car's. If everyone keeps a 2-second gap, a lane can never pass more than 3600/2 = 1,800 cars an hour, and a bit less once car length counts. Busy freeway lanes do carry around 2,000 an hour, so many drivers follow closer than 2 seconds.
import random
VMAX, P, CELL = 3, 0.15, 7.5 # cells per second, hesitation chance, metres per cell
ROAD, MERGE, ZONE = 240, 200, 10 # road length, where lane B ends, merge zone (cells)
def ahead(lane, i): # index of the next car in front, or None
return next((j for j in range(i+1, len(lane)) if lane[j]), None)
def behind(lane, i):
return next((j for j in range(i-1, -1, -1) if lane[j]), None)
def nasch(car, gap): # Nagel-Schreckenberg: speed up, keep the gap, hesitate
v = min(car['v'] + 1, VMAX, gap)
if v > 0 and random.random() < P: v -= 1
car['v'] = v
def run(mode, demand, minutes=20):
A, B = [None]*ROAD, [None]*MERGE
wait = {'A': 0, 'B': 0}; last = 'A'; out = 0
for t in range(minutes*60):
for i in range(MERGE-1, MERGE-ZONE-1, -1): # zipper: move over near the end
j = behind(A, i)
if B[i] and not A[i] and (j is None or i-j-1 >= A[j]['v']):
A[i], B[i] = B[i], None
hb = next((i for i in range(MERGE-1, MERGE-ZONE-1, -1) if B[i]), None)
yielder = behind(A, hb+1) if hb is not None and last == 'A' else None
nA, nB = [None]*ROAD, [None]*MERGE
for i in range(ROAD):
if not A[i]: continue
j = ahead(A, i); gap = 99 if j is None else j-i-1
if i == yielder and hb > i: gap = min(gap, hb-i-1) # let one car in
nasch(A[i], gap); k = i + A[i]['v']
if i < MERGE <= k: last = A[i]['lane']
if k >= ROAD: out += 1
else: nA[k] = A[i]
for i in range(MERGE):
if not B[i]: continue
j = ahead(B, i); gap = MERGE-1-i if j is None else j-i-1
nasch(B[i], gap); nB[i + B[i]['v']] = B[i]
A, B = nA, nB
for k in range(2): # arrivals, cars per hour
if random.random() < demand/2/3600:
wait['A' if mode == 'early' or k == 0 else 'B'] += 1
for name, lane in (('A', A), ('B', B)): # enter if the first cell is free
if wait[name] and not lane[0]:
j = ahead(lane, 0)
lane[0] = {'v': VMAX if j is None else min(VMAX, j-1), 'lane': name}
wait[name] -= 1
return out*60/minutes, wait
def capacity(v, T, length=5.0): # cars per hour per lane at speed v (m/s), time gap T (s)
return 3600*v/(v*T + length)
for mode in ('early', 'zipper'):
print(mode, run(mode, 2400))
print('2-second rule at 60 mph: %.0f cars/h, limit %.0f' % (capacity(26.8, 2.0), 3600/2.0))