Langevin car-following toward a density-dependent speed
Cars on a ring road whose speeds relax with noise towards a desired speed set by the local density.
- Name:
langevin_traffic - Field:
engineering/traffic - Formalism: stochastic differential equations
- References: Bertin 2021
Factory
@register_model(
field='engineering/traffic',
tags=('traffic', 'sde', 'langevin', 'car following'),
references=('Bertin 2021',),
defaults={'until': 200.0, 'dt': 1.0, 'solver': EulerMaruyamaSolver(dt=0.01)},
)
def langevin_traffic(
n_cars: int = 50,
length: float = 100.0,
*,
gamma: float = 3.0,
u0: float = 1.0,
rho_max: float = 1.0,
D1: float = 0.05,
) -> System:
"""Return cars on a ring road whose speeds relax with noise towards a desired speed set by the local density.
Car `i` follows `dv_i = -gamma (v_i - u(rho_i)) dt + sqrt(2 D(rho_i)) dW` with the desired speed
`u(rho) = u0 (1 - rho / rho_max)` and the noise intensity `D(rho) = D1 rho`, where the local density `rho_i` is the
inverse of the headway to the car ahead. Cars start evenly spaced at the speed `u(rho)` of the mean density
`rho = n_cars / length`. The mean speed stays at `u(rho)`, so the mean flow follows the fundamental diagram
`rho u(rho)`, and the speeds fluctuate with the variance `D1 rho / gamma`. The uniform flow is linearly stable while
`u0 rho^2 / rho_max < gamma / 2`, as for the optimal velocity model. Strong noise lets cars pass each other, and the
headway of a car that has passed the car ahead wraps around the ring.
Args:
n_cars:
The number of cars.
length:
The length of the ring road.
gamma:
The relaxation rate of the speed.
u0:
The desired speed on an empty road.
rho_max:
The jam density, at which the desired speed vanishes.
D1:
The noise intensity per unit density.
Returns:
system:
The model, with the group `cars` of positions `cars.x` and speeds `cars.v`, car `i + 1` ahead of car `i`.
Examples:
>>> from skcomplex.models.engineering.traffic import langevin_traffic
>>> float(langevin_traffic(n_cars=10, length=20.0).elements[0].values['v'][0])
0.5
"""
density = n_cars / length
positions = np.arange(n_cars) * length / n_cars
cars = GroupElement(
'cars',
element=CarElement,
size=n_cars,
values={'x': positions, 'v': np.full(n_cars, u0 * (1 - density / rho_max))},
)
driving = LangevinDrivingInteraction(
'drive',
target='cars',
gamma=gamma,
u0=u0,
rho_max=rho_max,
D1=D1,
road_length=length,
)
return System(elements=[cars], interactions=[driving])
Components
class CarElement(Element):
"""Car on a road whose position `x` advances at its speed `v`."""
def __init__(
self,
label: str,
visual: VisualSpec | None = None,
*,
x: Annotated[object, StateSpec(role='position')] = 0.0,
v: Annotated[object, StateSpec(role='velocity')] = 0.0,
) -> None:
super().__init__(label, visual=visual)
self.x = x
self.v = v
def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
"""Return the advance of the position at the speed."""
return [RateLaw(s.x, s.v)]
class LangevinDrivingInteraction(RelaxationInteraction):
"""Relaxation of the speed of every car of a ring road towards the speed of its local density, with noise.
The speed follows `dv = -gamma (v - u(rho)) dt + sqrt(2 D1 rho) dW`, with the desired speed
`u(rho) = u0 (1 - rho / rho_max)` of the density `rho`, the inverse of the headway to the next car, an
Ornstein-Uhlenbeck relaxation over the time `1 / gamma` with a noise that grows with the density.
Args:
label:
The label of the interaction.
target:
The path of the group of cars, in the order they drive.
gamma:
The relaxation rate.
u0:
The desired speed on an empty road.
rho_max:
The jam density.
D1:
The strength of the noise.
road_length:
The length of the ring road.
Examples:
>>> from skcomplex.models.engineering.traffic import LangevinDrivingInteraction
>>> LangevinDrivingInteraction('drive', target='cars', gamma=2.0).state
'v'
"""
def __init__(
self,
label: str,
*,
target: Annotated[str, RoleSpec()],
gamma: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 3.0,
u0: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
rho_max: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
D1: Annotated[object, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.05,
road_length: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 100.0,
) -> None:
super().__init__(label, target=target, state='v')
self.gamma = gamma
self.u0 = u0
self.rho_max = rho_max
self.D1 = D1
self.road_length = road_length
def _get_goal(self, namespace: SimpleNamespace) -> object:
"""Return the desired speed of the local density of every car."""
return LangevinDrivingInteraction._target_speed(
LangevinDrivingInteraction._local_density(namespace.target.x, namespace.road_length),
namespace,
)
def _get_time(self, namespace: SimpleNamespace) -> object:
"""Return the relaxation time `1 / gamma`."""
return 1 / namespace.gamma
def _get_noise(self, namespace: SimpleNamespace) -> object:
"""Return the amplitude `sqrt(2 D1 rho)` of the noise of the local density."""
return sqrt(
2 * namespace.D1 * LangevinDrivingInteraction._local_density(namespace.target.x, namespace.road_length),
)
@staticmethod
def _target_speed(density: Expression, s: SimpleNamespace) -> Expression:
"""Return the desired speed `u(rho) = u0 (1 - rho / rho_max)` at a density."""
return s.u0 * (1 - density / s.rho_max)
@staticmethod
def _local_density(x: Expression, road_length: Expression) -> Expression:
"""Return the density seen by every car of a ring road, the inverse of its headway."""
return 1 / OptimalVelocityInteraction._headway(x, road_length)
Equations
cars.x[i]: diff(x_cars(t), t) == v_cars(t)
cars.v[i]: diff(v_cars(t), t) == -(u0_drive*(1/(rho_max_drive*mod(x_cars__shift0 - x_cars(t), road_length_drive)) - 1) + v_cars(t))*gamma_drive