Skip to content

Traffic

Traffic: cars on roads.

ArrivalInteraction(label, *, target, demand=1.0, width=30, height=30, phase='arrival')

Bases: LifeCycleInteraction

Arrival of new riders at random cells with random destinations, a life cycle of arrivals only.

Every step a Poisson number of riders, with mean demand, appear at random cells with random destinations.

Parameters:

Name Type Description Default
label str

The label of the arrivals.

 required
target Annotated[str, RoleSpec()]

The path of the group of riders.

 required
demand Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The mean number of new riders in a step.

 1.0
width int

The number of columns of the grid.

 30
height int

The number of rows of the grid.

 30
phase str | None

The phase the riders arrive in.

 'arrival'

Examples:

>>> from skcomplex.models.engineering.traffic import ArrivalInteraction
>>> ArrivalInteraction('arrival', target='riders').reports
('births', 'deaths', 'arrivals')
Source code in src/skcomplex/models/engineering/traffic/_arrival.py 
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    demand: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    width: int = 30,
    height: int = 30,
    phase: str | None = 'arrival',
) -> None:
    super().__init__(
        label,
        target=target,
        placement=self.place,
        values=self.destine,
        arrival=self.expect,
        phase=phase,
    )
    self.demand = demand
    self.width = width
    self.height = height

destine(agent, model)

Return a random destination on the grid for a new rider.

Parameters:

Name Type Description Default
agent RiderView | None

The parent, None for an arrival.

 required
model ModelView

The model view of the arrivals.

 required

Returns:

Name Type Description
values dict[str, object]

The coordinates of the destination by state name.
Source code in src/skcomplex/models/engineering/traffic/_arrival.py 
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
def destine(self, agent: RiderView | None, model: ModelView) -> dict[str, object]:
    """Return a random destination on the grid for a new rider.

    Args:
        agent:
            The parent, `None` for an arrival.

        model:
            The model view of the arrivals.

    Returns:
        values:
            The coordinates of the destination by state name.
    """
    to_x, to_y = model.rng.integers((self.width, self.height))
    return {'destination_x': float(to_x), 'destination_y': float(to_y)}

expect(model)

Return the mean number of new riders in a step.

Parameters:

Name Type Description Default
model ModelView

The model view of the arrivals.

 required

Returns:

Name Type Description
demand float

The demand.
Source code in src/skcomplex/models/engineering/traffic/_arrival.py 
72
73
74
75
76
77
78
79
80
81
82
83
def expect(self, model: ModelView) -> float:
    """Return the mean number of new riders in a step.

    Args:
        model:
            The model view of the arrivals.

    Returns:
        demand:
            The demand.
    """
    return model.parameter('demand')

place(agent, model)

Return a random cell of the grid for a new rider.

Parameters:

Name Type Description Default
agent RiderView | None

The parent, None for an arrival.

 required
model ModelView

The model view of the arrivals.

 required

Returns:

Name Type Description
cell tuple[int, int]

The coordinates of the cell.
Source code in src/skcomplex/models/engineering/traffic/_arrival.py 
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
def place(self, agent: RiderView | None, model: ModelView) -> tuple[int, int]:
    """Return a random cell of the grid for a new rider.

    Args:
        agent:
            The parent, `None` for an arrival.

        model:
            The model view of the arrivals.

    Returns:
        cell:
            The coordinates of the cell.
    """
    x, y = model.rng.integers((self.width, self.height))
    return int(x), int(y)

CabElement(label, visual=VisualSpec('agents', {'color_by': 'task', 'color_range': (0.0, 2.0)}), *, task=IDLE, target_x=0.0, target_y=0.0, destination_x=0.0, destination_y=0.0, rider=-1.0, seats=1.0, money=0.0, base_fare=2.0, block_fare=0.5)

Bases: Element

Cab that drives one block a step to its rider, collects the fare, and drives the rider to its destination.

The fare is a base fare plus a fare per block of the trip, paid as a transfer from rider to cab.

Parameters:

Name Type Description Default
label str

The label of the agent or group.

 required
visual VisualSpec | None

How the element is drawn, or None for default_visual_.

 VisualSpec('agents', {'color_by': 'task', 'color_range': (0.0, 2.0)})
task Annotated[float, StateSpec(domain=NN)]

0 when idle, 1 on the way to a rider, and 2 on the way to a destination.

 IDLE
target_x Annotated[float, StateSpec(domain=NN)]

The column of the cell the cab drives to.

 0.0
target_y Annotated[float, StateSpec(domain=NN)]

The row of the cell the cab drives to.

 0.0
destination_x Annotated[float, StateSpec(domain=NN)]

The column of the destination of the rider.

 0.0
destination_y Annotated[float, StateSpec(domain=NN)]

The row of the destination of the rider.

 0.0
rider Annotated[float, StateSpec()]

The id of the rider the cab is sent to.

 -1.0
seats Annotated[float, StateSpec(domain=closed_open(0, oo))]

The free seats of the cab, one when idle and zero on a trip.

 1.0
money Annotated[float, StateSpec()]

The fares the cab collected.

 0.0
base_fare Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The fare of every trip.

 2.0
block_fare Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The fare per block of a trip.

 0.5

Examples:

>>> from skcomplex.core.elements import GroupElement
>>> from skcomplex.models.engineering.traffic import CabElement
>>> GroupElement('cabs', element=CabElement, size=3).size
3
Source code in src/skcomplex/models/engineering/traffic/_cab.py 
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'task', 'color_range': (0.0, 2.0)}),
    *,
    task: Annotated[float, StateSpec(domain=NN)] = IDLE,
    target_x: Annotated[float, StateSpec(domain=NN)] = 0.0,
    target_y: Annotated[float, StateSpec(domain=NN)] = 0.0,
    destination_x: Annotated[float, StateSpec(domain=NN)] = 0.0,
    destination_y: Annotated[float, StateSpec(domain=NN)] = 0.0,
    rider: Annotated[float, StateSpec()] = -1.0,
    seats: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    money: Annotated[float, StateSpec()] = 0.0,
    base_fare: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 2.0,
    block_fare: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
) -> None:
    super().__init__(label, visual=visual)
    self.task = task
    self.target_x = target_x
    self.target_y = target_y
    self.destination_x = destination_x
    self.destination_y = destination_y
    self.rider = rider
    self.seats = seats
    self.money = money
    self.base_fare = base_fare
    self.block_fare = block_fare

drive(agent, model)

Drive one block toward the target, and pick up or drop off the rider there.

Parameters:

Name Type Description Default
agent CabView

The cab.

 required
model ModelView

The model.

 required
Source code in src/skcomplex/models/engineering/traffic/_cab.py 
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
def drive(self, agent: CabView, model: ModelView) -> None:
    """Drive one block toward the target, and pick up or drop off the rider there.

    Args:
        agent:
            The cab.

        model:
            The model.
    """
    if agent.task == IDLE:
        return
    x, y = CabElement._cell(agent.location)
    target = (int(agent.target_x), int(agent.target_y))
    if (x, y) != target:
        step = (int(np.sign(target[0] - x)), 0) if x != target[0] else (0, int(np.sign(target[1] - y)))
        agent.move_to((x + step[0], y + step[1]))
    if CabElement._cell(agent.location) != target:
        return
    if agent.task == PICKING_UP:
        self._pick_up(agent, model)
    else:
        agent.task = IDLE
        agent.seats = 1.0

CarelessBrakingInteraction(label, *, cars, probability=0.1, gap=1.0)

Bases: Interaction

Model rule of Charbonneau's traffic model: a random driver brakes now and then, and every car moves.

Every step, with probability probability, one random driver brakes by the braking of its group whatever the distance to the car ahead. Then every car moves by its speed, but no closer than gap to the car ahead.

Parameters:

Name Type Description Default
label str

The label of the interaction.

 required
cars Annotated[str, RoleSpec()]

The path of the group of cars.

 required
probability Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The probability that a random driver brakes in a step.

 0.1
gap Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The smallest distance a car keeps to the car ahead.

 1.0

Examples:

>>> from skcomplex.models.engineering.traffic import CarelessBrakingInteraction
>>> CarelessBrakingInteraction('bozo', cars='cars').probability
0.1
Source code in src/skcomplex/models/engineering/traffic/_careless_braking.py 
44
45
46
47
48
49
50
51
52
53
54
55
def __init__(
    self,
    label: str,
    *,
    cars: Annotated[str, RoleSpec()],
    probability: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    gap: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label)
    self.cars = cars
    self.probability = probability
    self.gap = gap

drive(model)

Make a random driver brake with the interaction's probability, then move every car.

Parameters:

Name Type Description Default
model ModelView

The model.

 required
Source code in src/skcomplex/models/engineering/traffic/_careless_braking.py 
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
def drive(self, model: ModelView) -> None:
    """Make a random driver brake with the interaction's probability, then move every car.

    Args:
        model:
            The model.
    """
    cars = sorted(model.agents('cars', MotoristView), key=lambda car: car.rank)
    if not cars:
        return
    if model.rng.random() < model.parameter('probability'):
        careless = cars[int(model.rng.integers(len(cars)))]
        careless.v = max(careless.v - careless.braking, 0.0)
    gap = model.parameter('gap')
    positions = [car.x for car in cars]
    for index, car in enumerate(cars):
        moved = car.x + car.v
        car.x = min(moved, positions[index + 1] - gap) if index + 1 < len(cars) else moved

DriverElement(label, visual=VisualSpec('agents', {'color_by': 'speed'}), *, speed=0.0, rank=0.0, acceleration=1.0, min_acc=-10.0, max_acc=1.0, speed_limit=40.0, eps=0.02)

Bases: Element

Driver on a ring road who accelerates as much as allowed and stops when the speed would cause a collision.

The speed changes by the acceleration, bounded by min_acc and max_acc, is multiplied by a random factor between 1 - eps and 1 + eps, and stays between 0 and the speed limit. A speed larger than the distance to the car ahead becomes 0. The car then moves along the ring by its speed, and the ring is as long as the width of the space.

Parameters:

Name Type Description Default
label str

The label of the agent or group.

 required
visual VisualSpec | None

How the element is drawn, or None for default_visual_.

 VisualSpec('agents', {'color_by': 'speed'})
speed Annotated[float, StateSpec(domain=closed_open(0, oo))]

The speed.

 0.0
rank Annotated[float, StateSpec(domain=NN)]

The order of the car along the ring.

 0.0
acceleration Annotated[float, ParameterSpec()]

The acceleration every driver chooses.

 1.0
min_acc Annotated[float, ParameterSpec()]

The largest braking, as a negative acceleration.

 -10.0
max_acc Annotated[float, ParameterSpec()]

The largest acceleration.

 1.0
speed_limit Annotated[float, ParameterSpec(domain=open(0, oo))]

The largest speed.

 40.0
eps Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The relative amplitude of the random error of speed.

 0.02

Examples:

>>> from skcomplex.core.elements import GroupElement
>>> from skcomplex.models.engineering.traffic import DriverElement
>>> GroupElement('drivers', element=DriverElement, size=2, location=[(0.0, 0.5), (10.0, 0.5)]).size
2
Source code in src/skcomplex/models/engineering/traffic/_driver.py 
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'speed'}),
    *,
    speed: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    rank: Annotated[float, StateSpec(domain=NN)] = 0.0,
    acceleration: Annotated[float, ParameterSpec()] = 1.0,
    min_acc: Annotated[float, ParameterSpec()] = -10.0,
    max_acc: Annotated[float, ParameterSpec()] = 1.0,
    speed_limit: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 40.0,
    eps: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.02,
) -> None:
    super().__init__(label, visual=visual)
    self.speed = speed
    self.rank = rank
    self.acceleration = acceleration
    self.min_acc = min_acc
    self.max_acc = max_acc
    self.speed_limit = speed_limit
    self.eps = eps

drive(agent, model)

Choose the speed from the distance to the car ahead and move along the ring.

Parameters:

Name Type Description Default
agent DriverView

The driver.

 required
model ModelView

The model.

 required
Source code in src/skcomplex/models/engineering/traffic/_driver.py 
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
def drive(self, agent: DriverView, model: ModelView) -> None:
    """Choose the speed from the distance to the car ahead and move along the ring.

    Args:
        agent:
            The driver.

        model:
            The model.
    """
    cars = model.agents(agent.group, DriverView)
    if model.dimensions is None:
        error_msg = 'Drivers move along a ring road, a continuous space with an extent.'
        raise SimulationError(error_msg)
    length = model.dimensions[0]
    ahead = DriverElement._car_ahead(agent, model, len(cars))
    x, y = DriverElement._position(agent)
    distance = (DriverElement._position(ahead)[0] - x) % length if ahead is not None else length
    acceleration = min(max(agent.acceleration, agent.min_acc), agent.max_acc)
    speed = (agent.speed + acceleration) * model.rng.uniform(1 - agent.eps, 1 + agent.eps)
    speed = min(max(speed, 0.0), agent.speed_limit)
    agent.speed = 0.0 if speed > distance else speed
    agent.move_to(((x + agent.speed) % length, y))

GasKineticTrafficElement(label, visual=None, *, rho=20.0, Q=1800.0, V0=110.0, T=1.8 / 3600, tau=35 / 3600, rho_max=160.0, gamma=1.2, A0=0.008, dA=0.02, rho_c=0.27 * 160, d_rho=0.05 * 160)

Bases: Element

Cell of a road with the traffic density rho and the flow Q of the gas-kinetic-based traffic model.

The density and the flow are conserved by conservation interactions, and the flow relaxes to a dynamic equilibrium speed that brakes on the density and the speed read ahead at the interaction point. Lengths are in kilometres and times in hours, so T and tau are fractions of an hour.

Parameters:

Name Type Description Default
label str

The label of the cell or the road.

 required
visual VisualSpec | None

How the cell is drawn, or None for default_visual_.

 None
rho Annotated[float, StateSpec(domain=open(0, oo))]

The initial density in vehicles per kilometre.

 20.0
Q Annotated[float, StateSpec()]

The initial flow in vehicles per hour.

 1800.0
V0 Annotated[float, ParameterSpec(domain=open(0, oo))]

The desired speed.

 110.0
T Annotated[float, ParameterSpec(domain=open(0, oo))]

The safe time headway.

 1.8 / 3600
tau Annotated[float, ParameterSpec(domain=open(0, oo))]

The relaxation time of the speed.

 35 / 3600
rho_max Annotated[float, ParameterSpec(domain=open(0, oo))]

The maximum density.

 160.0
gamma Annotated[float, ParameterSpec(domain=open(0, oo))]

The anticipation factor, the interaction distance in units of the safe distance 1 / rho_max + T V.

 1.2
A0 Annotated[float, ParameterSpec(domain=open(0, oo))]

The variance prefactor of free traffic.

 0.008
dA Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The increase of the variance prefactor from free to congested traffic, halved.

 0.02
rho_c Annotated[float, ParameterSpec(domain=open(0, oo))]

The density of the transition of the variance prefactor.

 0.27 * 160
d_rho Annotated[float, ParameterSpec(domain=open(0, oo))]

The width of the transition of the variance prefactor.

 0.05 * 160

Examples:

>>> from skcomplex.models.engineering.traffic import GasKineticTrafficElement
>>> GasKineticTrafficElement('cell').gamma
1.2
Source code in src/skcomplex/models/engineering/traffic/_gas_kinetic_traffic.py 
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    rho: Annotated[float, StateSpec(domain=RealSet.open(0, oo))] = 20.0,
    Q: Annotated[float, StateSpec()] = 1800.0,
    V0: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 110.0,
    T: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.8 / 3600,
    tau: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 35 / 3600,
    rho_max: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 160.0,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.2,
    A0: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.008,
    dA: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.02,
    rho_c: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.27 * 160,
    d_rho: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.05 * 160,
) -> None:
    super().__init__(label, visual=visual)
    self.rho = rho
    self.Q = Q
    self.V0 = V0
    self.T = T
    self.tau = tau
    self.rho_max = rho_max
    self.gamma = gamma
    self.A0 = A0
    self.dA = dA
    self.rho_c = rho_c
    self.d_rho = d_rho

LangevinDrivingInteraction(label, *, target, gamma=3.0, u0=1.0, rho_max=1.0, D1=0.05, road_length=100.0)

Bases: 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.

Parameters:

Name Type Description Default
label str

The label of the interaction.

 required
target Annotated[str, RoleSpec()]

The path of the group of cars, in the order they drive.

 required
gamma Annotated[object, ParameterSpec(domain=open(0, oo))]

The relaxation rate.

 3.0
u0 Annotated[object, ParameterSpec(domain=open(0, oo))]

The desired speed on an empty road.

 1.0
rho_max Annotated[object, ParameterSpec(domain=open(0, oo))]

The jam density.

 1.0
D1 Annotated[object, ParameterSpec(domain=closed_open(0, oo))]

The strength of the noise.

 0.05
road_length Annotated[object, ParameterSpec(domain=open(0, oo))]

The length of the ring road.

 100.0

Examples:

>>> from skcomplex.models.engineering.traffic import LangevinDrivingInteraction
>>> LangevinDrivingInteraction('drive', target='cars', gamma=2.0).state
'v'
Source code in src/skcomplex/models/engineering/traffic/_langevin_driving.py 
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
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

MotoristElement(label, visual=VisualSpec('agents', {'position': 'x', 'color_by': 'v'}), *, x=0.0, v=0.0, rank=0.0, safe_distance=5.0, braking=3.0, acceleration=1.0, speed_limit=10.0)

Bases: Element

Car of Charbonneau's traffic model, which brakes when it is closer than a safe distance to the car ahead.

Parameters:

Name Type Description Default
label str

The label of the agent or group.

 required
visual VisualSpec | None

How the element is drawn, or None for default_visual_.

 VisualSpec('agents', {'position': 'x', 'color_by': 'v'})
x Annotated[float, StateSpec(role=position)]

The position along the road.

 0.0
v Annotated[float, StateSpec(domain=closed_open(0, oo), role=velocity)]

The speed.

 0.0
rank Annotated[float, StateSpec(domain=NN)]

The order of the car along the road, 0 for the last car.

 0.0
safe_distance Annotated[float, ParameterSpec(domain=open(0, oo))]

The distance to the car ahead below which the driver brakes.

 5.0
braking Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The speed lost when braking.

 3.0
acceleration Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The speed gained when accelerating.

 1.0
speed_limit Annotated[float, ParameterSpec(domain=open(0, oo))]

The largest speed.

 10.0

Examples:

>>> from skcomplex.core.elements import GroupElement
>>> from skcomplex.models.engineering.traffic import MotoristElement
>>> GroupElement('cars', element=MotoristElement, size=3, values={'x': [0.0, 10.0, 20.0]}).size
3
Source code in src/skcomplex/models/engineering/traffic/_motorist.py 
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'position': 'x', 'color_by': 'v'}),
    *,
    x: Annotated[float, StateSpec(role='position')] = 0.0,
    v: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo), role='velocity')] = 0.0,
    rank: Annotated[float, StateSpec(domain=NN)] = 0.0,
    safe_distance: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 5.0,
    braking: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 3.0,
    acceleration: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    speed_limit: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 10.0,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.v = v
    self.rank = rank
    self.safe_distance = safe_distance
    self.braking = braking
    self.acceleration = acceleration
    self.speed_limit = speed_limit

adjust(agent, model)

Brake when the car ahead is closer than the safe distance and accelerate otherwise.

The lead car, which has no car ahead, keeps its distance to the car behind it instead.

Parameters:

Name Type Description Default
agent MotoristView

The car.

 required
model ModelView

The model.

 required
Source code in src/skcomplex/models/engineering/traffic/_motorist.py 
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
def adjust(self, agent: MotoristView, model: ModelView) -> None:
    """Brake when the car ahead is closer than the safe distance and accelerate otherwise.

    The lead car, which has no car ahead, keeps its distance to the car behind it instead.

    Args:
        agent:
            The car.

        model:
            The model.
    """
    cars = model.agents(agent.group, MotoristView)
    ahead = next((car for car in cars if int(car.rank) == int(agent.rank) + 1), None)
    if ahead is None:
        behind = next((car for car in cars if int(car.rank) == int(agent.rank) - 1), None)
        distance = agent.x - behind.x if behind is not None else np.inf
    else:
        distance = ahead.x - agent.x
    change = -agent.braking if distance < agent.safe_distance else agent.acceleration
    agent.v = float(np.clip(agent.v + change, 0.0, agent.speed_limit))

RiderElement(label, visual=VisualSpec('agents', {'color_by': 'waited'}), *, destination_x=0.0, destination_y=0.0, waited=0.0, seats=1.0, money=100.0)

Bases: Element

Rider who waits at a cell for a cab to a destination, and pays the fare when the cab picks it up.

Parameters:

Name Type Description Default
label str

The label of the agent or group.

 required
visual VisualSpec | None

How the element is drawn, or None for default_visual_.

 VisualSpec('agents', {'color_by': 'waited'})
destination_x Annotated[float, StateSpec(domain=NN)]

The column of the destination.

 0.0
destination_y Annotated[float, StateSpec(domain=NN)]

The row of the destination.

 0.0
waited Annotated[float, StateSpec(domain=closed_open(0, oo))]

The number of steps the rider has waited.

 0.0
seats Annotated[float, StateSpec(domain=closed_open(0, oo))]

The seats the rider still requests, one until a cab is on its way and zero after.

 1.0
money Annotated[float, StateSpec()]

The money of the rider.

 100.0

Examples:

>>> from skcomplex.core.elements import GroupElement
>>> from skcomplex.models.engineering.traffic import RiderElement
>>> GroupElement('riders', element=RiderElement, size=2, values={'money': 50.0}).size
2
Source code in src/skcomplex/models/engineering/traffic/_rider.py 
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = VisualSpec('agents', {'color_by': 'waited'}),
    *,
    destination_x: Annotated[float, StateSpec(domain=NN)] = 0.0,
    destination_y: Annotated[float, StateSpec(domain=NN)] = 0.0,
    waited: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
    seats: Annotated[float, StateSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    money: Annotated[float, StateSpec()] = 100.0,
) -> None:
    super().__init__(label, visual=visual)
    self.destination_x = destination_x
    self.destination_y = destination_y
    self.waited = waited
    self.seats = seats
    self.money = money

RoadElement(label, visual=None, *, velocity=-1, slowdown=0.3)

Bases: Element

Cell of the Nagel-Schreckenberg road: the velocity of the car on it, or -1 when it is empty.

The local rules of the road make cars accelerate by one up to the largest velocity, brake to the gap ahead, slow down by one with probability slowdown, and move forward by their velocity, every step.

Parameters:

Name Type Description Default
label str

The label of the cell or the lattice group.

 required
visual VisualSpec | None

How the element is drawn, or None for default_visual_.

 None
velocity Annotated[float, StateSpec(domain=ZZ)]

The initial velocity, or -1 for an empty cell.

 -1
slowdown Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The probability of random braking.

 0.3

Examples:

>>> from skcomplex.models.engineering.traffic import RoadElement
>>> RoadElement('cell').slowdown
0.3
Source code in src/skcomplex/models/engineering/traffic/_traffic.py 
65
66
67
68
69
70
71
72
73
74
75
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    velocity: Annotated[float, StateSpec(domain=ZZ)] = -1,
    slowdown: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.3,
) -> None:
    super().__init__(label, visual=visual)
    self.velocity = velocity
    self.slowdown = slowdown

TrafficMomentsElement(label, visual=None, *, rho=0.5, w=0.25, gamma=1.0, u0=1.0, rho_max=1.0, D1=0.2)

Bases: Element

Density rho and flow w of Langevin traffic on a road, the first two moments of its speeds.

Source code in src/skcomplex/models/engineering/traffic/_traffic_moments.py 
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    rho: Annotated[object, StateSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
    w: Annotated[object, StateSpec()] = 0.25,
    gamma: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 1.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.2,
) -> None:
    super().__init__(label, visual=visual)
    self.rho = rho
    self.w = w
    self.gamma = gamma
    self.u0 = u0
    self.rho_max = rho_max
    self.D1 = D1

charbonneau_traffic(n_cars=100, *, probability=0.1, safe_distance=5.0, speed_limit=10.0, random_state=0)

Return Charbonneau's traffic model, where random braking of single drivers sets off traffic jams.

Cars on a one-way road start at rest with random spacings between 3 and 17. Every step, each driver brakes by 3 when the car ahead is closer than the safe distance and accelerates by 1 otherwise, within the speed limit. A random driver brakes with a probability every step, and then every car moves by its speed.

Parameters:

Name Type Description Default
n_cars int

The number of cars.

 100
probability float

The probability that a random driver brakes in a step.

 0.1
safe_distance float

The distance to the car ahead below which drivers brake.

 5.0
speed_limit float

The largest speed.

 10.0
random_state int | Generator | None

The seed or generator of the initial spacings.

 0

Returns:

Name Type Description
system System

The model, with the group cars of positions x, speeds v, and ranks rank.

Examples:

>>> from skcomplex.models.engineering.traffic import charbonneau_traffic
>>> charbonneau_traffic(5).elements[0].check().prepare().size_
5
Source code in src/skcomplex/models/engineering/traffic/_drivers.py 
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'jams', 'agents', 'avalanches', 'stochastic'),
    references=('Charbonneau 2017',),
    defaults={'until': 100},
)
def charbonneau_traffic(
    n_cars: int = 100,
    *,
    probability: float = 0.1,
    safe_distance: float = 5.0,
    speed_limit: float = 10.0,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return Charbonneau's traffic model, where random braking of single drivers sets off traffic jams.

    Cars on a one-way road start at rest with random spacings between 3 and 17. Every step, each driver brakes by 3
    when the car ahead is closer than the safe distance and accelerates by 1 otherwise, within the speed limit. A random
    driver brakes with a probability every step, and then every car moves by its speed.

    Args:
        n_cars:
            The number of cars.

        probability:
            The probability that a random driver brakes in a step.

        safe_distance:
            The distance to the car ahead below which drivers brake.

        speed_limit:
            The largest speed.

        random_state:
            The seed or generator of the initial spacings.

    Returns:
        system:
            The model, with the group `cars` of positions `x`, speeds `v`, and ranks `rank`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import charbonneau_traffic
        >>> charbonneau_traffic(5).elements[0].check().prepare().size_
        5
    """
    rng = np.random.default_rng(random_state)
    positions = np.cumsum(rng.uniform(3.0, 17.0, n_cars))
    cars = GroupElement(
        'cars',
        element=MotoristElement,
        size=n_cars,
        values={
            'x': positions,
            'v': np.zeros(n_cars),
            'rank': np.arange(n_cars, dtype=float),
            'safe_distance': safe_distance,
            'speed_limit': speed_limit,
        },
    )
    return System(
        elements=[cars],
        interactions=[CarelessBrakingInteraction('braking', cars='cars', probability=probability)],
    )

delayed_car_following(n_cars=8, *, spacing=10.0, u=1.0, amplitude=0.1, omega=0.5, reaction=1.0, T=1.0, power=0.0)

Return a chain of cars that react to the speed difference to the car ahead after a reaction time.

Car j + 1 follows car j as x''_{j+1}(t) = λ Δv(t - T) / Δx(t - T)^power, where Δv and Δx are the speed and position differences to the car ahead. The lead car drives at u (1 - amplitude sin(ω t)). Oscillations grow along the chain when ω < 2 λ sin(ω T) for the linear model, power = 0.

Parameters:

Name Type Description Default
n_cars int

The number of cars, the lead car included.

 8
spacing float

The initial distance between neighbouring cars.

 10.0
u float

The mean speed of the lead car.

 1.0
amplitude float

The relative amplitude of the speed oscillation of the lead car.

 0.1
omega float

The angular frequency of the speed oscillation of the lead car.

 0.5
reaction float

The reaction constant λ.

 1.0
T float

The reaction time.

 1.0
power float

The power of the distance in the reaction, 0 for the linear model.

 0.0

Returns:

Name Type Description
system System

The model, with the lead car car_0 and the following cars car_<j>, each with position x and speed v.

Examples:

>>> from skcomplex.models.engineering.traffic import delayed_car_following
>>> delayed_car_following(n_cars=3).elements[2].x
-20.0
Source code in src/skcomplex/models/engineering/traffic/_car_following.py 
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'delay', 'car following', 'instability'),
    references=('Gros 2015', 'Chandler, Herman, and Montroll 1958'),
    defaults={'until': 60.0, 'dt': 0.5},
)
def delayed_car_following(
    n_cars: int = 8,
    *,
    spacing: float = 10.0,
    u: float = 1.0,
    amplitude: float = 0.1,
    omega: float = 0.5,
    reaction: float = 1.0,
    T: float = 1.0,
    power: float = 0.0,
) -> System:
    """Return a chain of cars that react to the speed difference to the car ahead after a reaction time.

    Car `j + 1` follows car `j` as `x''_{j+1}(t) = λ Δv(t - T) / Δx(t - T)^power`, where `Δv` and `Δx` are the speed
    and position differences to the car ahead. The lead car drives at `u (1 - amplitude sin(ω t))`. Oscillations grow
    along the chain when `ω < 2 λ sin(ω T)` for the linear model, `power = 0`.

    Args:
        n_cars:
            The number of cars, the lead car included.

        spacing:
            The initial distance between neighbouring cars.

        u:
            The mean speed of the lead car.

        amplitude:
            The relative amplitude of the speed oscillation of the lead car.

        omega:
            The angular frequency of the speed oscillation of the lead car.

        reaction:
            The reaction constant `λ`.

        T:
            The reaction time.

        power:
            The power of the distance in the reaction, 0 for the linear model.

    Returns:
        system:
            The model, with the lead car `car_0` and the following cars `car_<j>`, each with position `x` and speed `v`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import delayed_car_following
        >>> delayed_car_following(n_cars=3).elements[2].x
        -20.0
    """
    lead = LeadCarElement('car_0', x=0.0, v=u, u=u, amplitude=amplitude, omega=omega)
    followers = [CarElement(f'car_{index}', x=-index * spacing, v=u) for index in range(1, n_cars)]
    interactions = [
        DelayedFollowingInteraction(
            f'follow_{index}',
            target=f'car_{index}',
            leader=f'car_{index - 1}',
            reaction=reaction,
            T=T,
            power=power,
        )
        for index in range(1, n_cars)
    ]
    return System(elements=[lead, *followers], interactions=interactions)

gkt_traffic(n=200, *, length=10.0, density=35.0, amplitude=2.0, V0=110.0, T=1.8 / 3600, tau=35 / 3600, rho_max=160.0, gamma=1.2)

Return the gas-kinetic-based traffic model, where drivers brake on the traffic at an interaction point ahead.

The density and the flow Q = rho V on a ring road follow rho_t + Q_x = 0 and Q_t + (Q**2 / rho + rho A(rho) V**2)_x = rho (Ve* - V) / tau, where the dynamic equilibrium speed

Ve* = V0 (1 - (theta + theta_a) / (2 A(rho_max)) (rho_a T / (1 - rho_a / rho_max))**2 B(dV))

reads the density rho_a and the speed V_a at the interaction point x + gamma (1 / rho_max + T V), with the variances theta = A(rho) V**2, the variance prefactor A(rho) = A0 + dA (tanh((rho - rho_c) / d_rho) + 1), the dimensionless speed difference dV = (V - V_a) / sqrt(theta + theta_a), and the Boltzmann factor B(d) = 2 (d phi(d) + (1 + d**2) Phi(d)) of the standard normal density phi and distribution Phi. Lengths are in kilometres and times in hours, and the defaults are the calibrated values of Treiber, Hennecke, and Helbing. The flow starts in equilibrium with their dipole perturbation of the density, of height amplitude. Homogeneous traffic is unstable between about 24 and 51 vehicles per kilometre, where the perturbation grows into stop-and-go waves whose downstream fronts move against the traffic.

Parameters:

Name Type Description Default
n int

The number of cells.

 200
length float

The length of the ring road in kilometres.

 10.0
density float

The mean density in vehicles per kilometre.

 35.0
amplitude float

The height of the dipole perturbation in vehicles per kilometre.

 2.0
V0 float

The desired speed in kilometres per hour.

 110.0
T float

The safe time headway in hours.

 1.8 / 3600
tau float

The relaxation time of the speed in hours.

 35 / 3600
rho_max float

The maximum density in vehicles per kilometre.

 160.0
gamma float

The anticipation factor.

 1.2

Returns:

Name Type Description
system System

The model, with the density road.rho and the flow road.Q.

Examples:

>>> from skcomplex.models.engineering.traffic import gkt_traffic
>>> gkt_traffic(10).elements[0].values['gamma']
1.2
Source code in src/skcomplex/models/engineering/traffic/_continuum.py 
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'field', 'conservation law', 'finite volume', 'non-local', 'stop-and-go', 'pde'),
    references=('Treiber, Hennecke, and Helbing 1999', 'arXiv:cond-mat/0012229'),
    defaults={'until': 1.0, 'dt': 1 / 60, 'solver': SSPRKSolver(cfl=0.4)},
)
def gkt_traffic(
    n: int = 200,
    *,
    length: float = 10.0,
    density: float = 35.0,
    amplitude: float = 2.0,
    V0: float = 110.0,
    T: float = 1.8 / 3600,
    tau: float = 35 / 3600,
    rho_max: float = 160.0,
    gamma: float = 1.2,
) -> System:
    """Return the gas-kinetic-based traffic model, where drivers brake on the traffic at an interaction point ahead.

    The density and the flow `Q = rho V` on a ring road follow `rho_t + Q_x = 0` and
    `Q_t + (Q**2 / rho + rho A(rho) V**2)_x = rho (Ve* - V) / tau`, where the dynamic equilibrium speed

        Ve* = V0 (1 - (theta + theta_a) / (2 A(rho_max)) (rho_a T / (1 - rho_a / rho_max))**2 B(dV))

    reads the density `rho_a` and the speed `V_a` at the interaction point `x + gamma (1 / rho_max + T V)`, with the
    variances `theta = A(rho) V**2`, the variance prefactor `A(rho) = A0 + dA (tanh((rho - rho_c) / d_rho) + 1)`, the
    dimensionless speed difference `dV = (V - V_a) / sqrt(theta + theta_a)`, and the Boltzmann factor
    `B(d) = 2 (d phi(d) + (1 + d**2) Phi(d))` of the standard normal density `phi` and distribution `Phi`. Lengths are
    in kilometres and times in hours, and the defaults are the calibrated values of Treiber, Hennecke, and Helbing. The
    flow starts in equilibrium with their dipole perturbation of the density, of height `amplitude`. Homogeneous
    traffic is unstable between about 24 and 51 vehicles per kilometre, where the perturbation grows into stop-and-go
    waves whose downstream fronts move against the traffic.

    Args:
        n:
            The number of cells.

        length:
            The length of the ring road in kilometres.

        density:
            The mean density in vehicles per kilometre.

        amplitude:
            The height of the dipole perturbation in vehicles per kilometre.

        V0:
            The desired speed in kilometres per hour.

        T:
            The safe time headway in hours.

        tau:
            The relaxation time of the speed in hours.

        rho_max:
            The maximum density in vehicles per kilometre.

        gamma:
            The anticipation factor.

    Returns:
        system:
            The model, with the density `road.rho` and the flow `road.Q`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import gkt_traffic
        >>> gkt_traffic(10).elements[0].values['gamma']
        1.2
    """
    centers = (np.arange(n) + 0.5) * length / n
    start, width, wake = 0.3 * length, 0.20125, 0.805
    rho = density + amplitude * (
        np.cosh((centers - start) / width) ** -2 - width / wake * np.cosh((centers - start - 1.00625) / wake) ** -2
    )
    values = {'V0': V0, 'T': T, 'tau': tau, 'rho_max': rho_max, 'gamma': gamma}
    cell = GasKineticTrafficElement('cell', V0=V0, T=T, tau=tau, rho_max=rho_max, gamma=gamma)
    road = FieldElement(
        'road',
        element=GasKineticTrafficElement,
        space=LatticeSpace(n, spacing=length / n),
        values={'rho': rho, 'Q': rho * _gkt_equilibrium_speed(rho, cell), **values},
        discretization='finite-volume',
    )
    vehicles = ConservationInteraction('density', target='road', state='rho', flux=lambda s: s.target.Q)
    momentum = ConservationInteraction(
        'flow',
        target='road',
        state='Q',
        flux=_get_gkt_momentum_flux,
        source=_get_gkt_relaxation,
    )
    return System(elements=[road], interactions=[vehicles, momentum])

intelligent_driver_model(n_cars=20, length=500.0, *, v0=30.0, T=1.5, a=1.0, b=1.5, s0=2.0, delta=4.0, car_length=5.0, perturbation=1.0)

Return the intelligent driver model of cars on a ring road.

Car i accelerates as v' = a (1 - (v / v0)^δ - (s* / s)^2) with the gap s to the car ahead and the desired gap s* = s0 + v T + v Δv / (2 sqrt(a b)), where Δv is the speed difference to the car ahead. Cars start at rest. The cars are one group, and the car ahead of car i is member i + 1 along the ring.

Parameters:

Name Type Description Default
n_cars int

The number of cars.

 20
length float

The length of the ring road.

 500.0
v0 float

The desired speed.

 30.0
T float

The safe time headway.

 1.5
a float

The largest acceleration.

 1.0
b float

The comfortable deceleration.

 1.5
s0 float

The jam distance.

 2.0
delta float

The acceleration exponent.

 4.0
car_length float

The length of a car, subtracted from the headway to give the gap.

 5.0
perturbation float

The forward shift of the first car from uniform spacing.

 1.0

Returns:

Name Type Description
system 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 intelligent_driver_model
>>> intelligent_driver_model(n_cars=4).interactions[0].road_length
500.0
Source code in src/skcomplex/models/engineering/traffic/_car_following.py 
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'jams', 'ode', 'car following'),
    references=('Treiber, Hennecke, and Helbing 2000', 'Chowdhury, Santen, and Schadschneider 2000'),
    defaults={'until': 100.0, 'dt': 1.0},
)
def intelligent_driver_model(
    n_cars: int = 20,
    length: float = 500.0,
    *,
    v0: float = 30.0,
    T: float = 1.5,
    a: float = 1.0,
    b: float = 1.5,
    s0: float = 2.0,
    delta: float = 4.0,
    car_length: float = 5.0,
    perturbation: float = 1.0,
) -> System:
    """Return the intelligent driver model of cars on a ring road.

    Car `i` accelerates as `v' = a (1 - (v / v0)^δ - (s* / s)^2)` with the gap `s` to the car ahead and the desired
    gap `s* = s0 + v T + v Δv / (2 sqrt(a b))`, where `Δv` is the speed difference to the car ahead. Cars start at rest.
    The cars are one group, and the car ahead of car `i` is member `i + 1` along the ring.

    Args:
        n_cars:
            The number of cars.

        length:
            The length of the ring road.

        v0:
            The desired speed.

        T:
            The safe time headway.

        a:
            The largest acceleration.

        b:
            The comfortable deceleration.

        s0:
            The jam distance.

        delta:
            The acceleration exponent.

        car_length:
            The length of a car, subtracted from the headway to give the gap.

        perturbation:
            The forward shift of the first car from uniform spacing.

    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 intelligent_driver_model
        >>> intelligent_driver_model(n_cars=4).interactions[0].road_length
        500.0
    """
    cars = GroupElement(
        'cars',
        element=CarElement,
        size=n_cars,
        values={'x': _ring_positions(n_cars, length, perturbation), 'v': 0.0},
    )
    following = IntelligentDriverInteraction(
        'follow',
        target='cars',
        v0=v0,
        T=T,
        a=a,
        b=b,
        s0=s0,
        delta=delta,
        length=car_length,
        road_length=length,
    )
    return System(elements=[cars], interactions=[following])

kerner_konhauser(n=40, *, length=10.0, density=0.25, bump=0.1, V0=1.0, rho_max=1.0, tau=1.0, c0=0.45, mu=0.2)

Return the Kerner-Konhäuser model, where viscous traffic flow forms dense moving clusters.

The density and the flow q = rho v on a ring road follow rho_t + q_x = 0 and q_t + (q^2 / rho + c0^2 rho)_x = (rho V(rho) - q) / τ + μ v_xx, with the equilibrium speed V(rho) = V0 (1 / (1 + exp((rho / rho_max - 0.25) / 0.06)) - 3.72e-6). The flow starts in equilibrium with a bump of density, which grows into a cluster at intermediate densities.

Parameters:

Name Type Description Default
n int

The number of cells.

 40
length float

The length of the ring road.

 10.0
density float

The mean density.

 0.25
bump float

The relative height of the density bump.

 0.1
V0 float

The scale of the equilibrium speed.

 1.0
rho_max float

The jam density.

 1.0
tau float

The relaxation time of the speed.

 1.0
c0 float

The speed of sound of the traffic pressure c0^2 rho.

 0.45
mu float

The viscosity.

 0.2

Returns:

Name Type Description
system System

The model, with the density road.rho and the flow road.q.

Examples:

>>> from skcomplex.models.engineering.traffic import kerner_konhauser
>>> kerner_konhauser(10).elements[0].values['mu']
0.2
Source code in src/skcomplex/models/engineering/traffic/_continuum.py 
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'field', 'conservation law', 'finite volume', 'jams', 'pde'),
    references=('Kerner and Konhäuser 1993', 'Chowdhury, Santen, and Schadschneider 2000'),
    defaults={'until': 100.0, 'dt': 0.2, 'solver': SSPRKSolver(cfl=0.15)},
)
def kerner_konhauser(
    n: int = 40,
    *,
    length: float = 10.0,
    density: float = 0.25,
    bump: float = 0.1,
    V0: float = 1.0,
    rho_max: float = 1.0,
    tau: float = 1.0,
    c0: float = 0.45,
    mu: float = 0.2,
) -> System:
    """Return the Kerner-Konhäuser model, where viscous traffic flow forms dense moving clusters.

    The density and the flow `q = rho v` on a ring road follow `rho_t + q_x = 0` and
    `q_t + (q^2 / rho + c0^2 rho)_x = (rho V(rho) - q) / τ + μ v_xx`, with the equilibrium speed
    `V(rho) = V0 (1 / (1 + exp((rho / rho_max - 0.25) / 0.06)) - 3.72e-6)`. The flow starts in equilibrium with a
    bump of density, which grows into a cluster at intermediate densities.

    Args:
        n:
            The number of cells.

        length:
            The length of the ring road.

        density:
            The mean density.

        bump:
            The relative height of the density bump.

        V0:
            The scale of the equilibrium speed.

        rho_max:
            The jam density.

        tau:
            The relaxation time of the speed.

        c0:
            The speed of sound of the traffic pressure `c0^2 rho`.

        mu:
            The viscosity.

    Returns:
        system:
            The model, with the density `road.rho` and the flow `road.q`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import kerner_konhauser
        >>> kerner_konhauser(10).elements[0].values['mu']
        0.2
    """
    rho = _density_bump(n, density, bump)
    values = {
        'rho': rho,
        'q': rho * _kerner_konhauser_speed(rho, V0, rho_max),
        'V0': V0,
        'rho_max': rho_max,
        'tau': tau,
        'c0': c0,
        'mu': mu,
    }
    space = LatticeSpace(n, spacing=length / n)
    road = FieldElement(
        'road',
        element=KernerKonhauserElement,
        space=space,
        values=values,
        discretization='finite-volume',
    )
    density = ConservationInteraction('density', target='road', state='rho', flux=lambda s: s.target.q)
    flow = ConservationInteraction(
        'flow',
        target='road',
        state='q',
        flux=lambda s: s.target.q**2 / s.target.rho + s.target.c0**2 * s.target.rho,
    )
    relaxation = RelaxationInteraction(
        'relaxation',
        target='road',
        state='q',
        goal=lambda s: s.target.rho * _get_kerner_konhauser_speed(s.target),
        time=lambda s: s.target.tau,
    )
    viscosity = RateInteraction(
        'viscosity',
        target='road',
        state='q',
        rate=lambda s: s.target.mu * laplacian(s.target.q / s.target.rho),
    )
    return System(elements=[road], interactions=[density, flow, relaxation, viscosity])

langevin_traffic(n_cars=50, length=100.0, *, gamma=3.0, u0=1.0, rho_max=1.0, D1=0.05)

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.

Parameters:

Name Type Description Default
n_cars int

The number of cars.

 50
length float

The length of the ring road.

 100.0
gamma float

The relaxation rate of the speed.

 3.0
u0 float

The desired speed on an empty road.

 1.0
rho_max float

The jam density, at which the desired speed vanishes.

 1.0
D1 float

The noise intensity per unit density.

 0.05

Returns:

Name Type Description
system 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
Source code in src/skcomplex/models/engineering/traffic/_langevin.py 
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
@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])

lighthill_whitham_richards(n=200, *, upstream=0.3, downstream=0.8, vmax=1.0, rho_max=1.0)

Return the Lighthill-Whitham-Richards model, the conservation of cars with the Greenshields speed law.

The density follows rho_t + (rho V(rho))_x = 0 with V(rho) = vmax (1 - rho / rho_max) on a road of length 1 with open edges. Free traffic at density upstream runs into a queue at density downstream in the second half of the road, and the shock between them moves at vmax (1 - (upstream + downstream) / rho_max).

Parameters:

Name Type Description Default
n int

The number of cells.

 200
upstream float

The density in the first half of the road.

 0.3
downstream float

The density in the second half of the road.

 0.8
vmax float

The speed on an empty road.

 1.0
rho_max float

The jam density.

 1.0

Returns:

Name Type Description
system System

The model, with the density road.rho.

Examples:

>>> from skcomplex.models.engineering.traffic import lighthill_whitham_richards
>>> float(lighthill_whitham_richards(10).elements[0].values['rho'][-1])
0.8
Source code in src/skcomplex/models/engineering/traffic/_continuum.py 
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'field', 'conservation law', 'finite volume', 'shocks', 'pde'),
    references=('Lighthill and Whitham 1955', 'Richards 1956', 'Chowdhury, Santen, and Schadschneider 2000'),
    defaults={'until': 1.0, 'dt': 0.05, 'solver': SSPRKSolver()},
)
def lighthill_whitham_richards(
    n: int = 200,
    *,
    upstream: float = 0.3,
    downstream: float = 0.8,
    vmax: float = 1.0,
    rho_max: float = 1.0,
) -> System:
    """Return the Lighthill-Whitham-Richards model, the conservation of cars with the Greenshields speed law.

    The density follows `rho_t + (rho V(rho))_x = 0` with `V(rho) = vmax (1 - rho / rho_max)` on a road of length 1
    with open edges. Free traffic at density `upstream` runs into a queue at density `downstream` in the second half of
    the road, and the shock between them moves at `vmax (1 - (upstream + downstream) / rho_max)`.

    Args:
        n:
            The number of cells.

        upstream:
            The density in the first half of the road.

        downstream:
            The density in the second half of the road.

        vmax:
            The speed on an empty road.

        rho_max:
            The jam density.

    Returns:
        system:
            The model, with the density `road.rho`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import lighthill_whitham_richards
        >>> float(lighthill_whitham_richards(10).elements[0].values['rho'][-1])
        0.8
    """
    density = np.where(np.arange(n) < n // 2, upstream, downstream).astype(float)
    values = {'rho': density, 'vmax': vmax, 'rho_max': rho_max}
    space = LatticeSpace(n, spacing=1 / n)
    road = FieldElement(
        'road',
        element=LighthillWhithamRichardsElement,
        space=space,
        boundary={'rho': 'zero-flux'},
        values=values,
        discretization='finite-volume',
    )
    flux = _get_lighthill_whitham_richards_flux
    density = ConservationInteraction('density', target='road', state='rho', flux=flux, riemann='godunov')
    return System(elements=[road], interactions=[density])

nagel_schreckenberg(length=1000, *, density=0.2, slowdown=0.3, vmax=5, random_state=0)

Return the Nagel-Schreckenberg traffic model on a ring road, where jams emerge above a critical density.

Parameters:

Name Type Description Default
length int

The number of cells of the road.

 1000
density float

The fraction of cells with a car.

 0.2
slowdown float

The probability of random braking.

 0.3
vmax int

The largest velocity.

 5
random_state int | Generator | None

The seed or generator of the initial positions of the cars.

 0

Returns:

Name Type Description
system System

The model, with the velocities road.velocity, -1 on empty cells.

Examples:

>>> from skcomplex.models.engineering.traffic import nagel_schreckenberg
>>> int((nagel_schreckenberg(100, density=0.1).elements[0].values['velocity'] >= 0).sum())
10
Source code in src/skcomplex/models/engineering/traffic/_traffic.py 
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
@register_model(
    field='engineering/traffic',
    tags=('lattice', 'traffic', 'jams', 'stochastic'),
    references=('Nagel and Schreckenberg 1992', 'Charbonneau 2017'),
    defaults={'until': 200, 'dt': 1},
)
def nagel_schreckenberg(
    length: int = 1000,
    *,
    density: float = 0.2,
    slowdown: float = 0.3,
    vmax: int = 5,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return the Nagel-Schreckenberg traffic model on a ring road, where jams emerge above a critical density.

    Args:
        length:
            The number of cells of the road.

        density:
            The fraction of cells with a car.

        slowdown:
            The probability of random braking.

        vmax:
            The largest velocity.

        random_state:
            The seed or generator of the initial positions of the cars.

    Returns:
        system:
            The model, with the velocities `road.velocity`, -1 on empty cells.

    Examples:
        >>> from skcomplex.models.engineering.traffic import nagel_schreckenberg
        >>> int((nagel_schreckenberg(100, density=0.1).elements[0].values['velocity'] >= 0).sum())
        10
    """
    rng = np.random.default_rng(random_state)
    velocity = np.full(length, -1.0)
    velocity[rng.choice(length, size=round(density * length), replace=False)] = 0.0
    road = LatticeElement(
        'road',
        element=RoadElement,
        space=LatticeSpace(length),
        values={'velocity': velocity, 'slowdown': slowdown},
    )
    ahead = LatticeNeighborhood(tuple((distance,) for distance in range(vmax + 1)))
    behind = LatticeNeighborhood(tuple((-distance,) for distance in range(vmax + 1)))
    return System(
        elements=[road],
        interactions=[
            LocalRuleInteraction(
                'acceleration',
                target='road',
                state='velocity',
                rule=_accelerate,
                neighborhood=ahead,
                stochastic=True,
            ),
            LocalRuleInteraction('movement', target='road', state='velocity', rule=_move, neighborhood=behind),
        ],
    )

optimal_velocity(n_cars=20, length=40.0, *, a=1.0, V0=1.0, hc=2.0, perturbation=0.1)

Return the optimal velocity model of cars on a ring road, where uniform flow breaks into stop-and-go waves.

Car i relaxes its speed towards the optimal velocity of its headway h, v_i' = a (V(h) - v_i) with V(h) = V0 (tanh(h - hc) + tanh(hc)). Uniform flow with headway b = length / n_cars is unstable when V'(b) > a / 2, and then a jam forms that moves against the direction of travel. The cars are one group, and the headway of car i is the distance along the ring to member i + 1.

Parameters:

Name Type Description Default
n_cars int

The number of cars.

 20
length float

The length of the ring road.

 40.0
a float

The sensitivity of drivers.

 1.0
V0 float

The scale of the optimal velocity, which tends to V0 (1 + tanh(hc)) on an empty road.

 1.0
hc float

The headway at which the optimal velocity changes fastest.

 2.0
perturbation float

The forward shift of the first car from uniform spacing.

 0.1

Returns:

Name Type Description
system 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 optimal_velocity
>>> system = optimal_velocity(n_cars=5)
>>> system.elements[0].size, len(system.interactions)
(5, 1)
Source code in src/skcomplex/models/engineering/traffic/_car_following.py 
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'jams', 'ode', 'car following'),
    references=('Bando et al. 1995', 'Chowdhury, Santen, and Schadschneider 2000'),
    defaults={'until': 200.0, 'dt': 1.0},
)
def optimal_velocity(
    n_cars: int = 20,
    length: float = 40.0,
    *,
    a: float = 1.0,
    V0: float = 1.0,
    hc: float = 2.0,
    perturbation: float = 0.1,
) -> System:
    """Return the optimal velocity model of cars on a ring road, where uniform flow breaks into stop-and-go waves.

    Car `i` relaxes its speed towards the optimal velocity of its headway `h`, `v_i' = a (V(h) - v_i)` with
    `V(h) = V0 (tanh(h - hc) + tanh(hc))`. Uniform flow with headway `b = length / n_cars` is unstable when
    `V'(b) > a / 2`, and then a jam forms that moves against the direction of travel. The cars are one group, and the
    headway of car `i` is the distance along the ring to member `i + 1`.

    Args:
        n_cars:
            The number of cars.

        length:
            The length of the ring road.

        a:
            The sensitivity of drivers.

        V0:
            The scale of the optimal velocity, which tends to `V0 (1 + tanh(hc))` on an empty road.

        hc:
            The headway at which the optimal velocity changes fastest.

        perturbation:
            The forward shift of the first car from uniform spacing.

    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 optimal_velocity
        >>> system = optimal_velocity(n_cars=5)
        >>> system.elements[0].size, len(system.interactions)
        (5, 1)
    """
    headway = length / n_cars
    speed = V0 * (np.tanh(headway - hc) + np.tanh(hc))
    cars = GroupElement(
        'cars',
        element=CarElement,
        size=n_cars,
        values={'x': _ring_positions(n_cars, length, perturbation), 'v': float(speed)},
    )
    following = OptimalVelocityInteraction('follow', target='cars', a=a, V0=V0, hc=hc, road_length=length)
    return System(elements=[cars], interactions=[following])

payne(n=200, *, density=0.3, bump=0.2, vmax=1.0, rho_max=1.0, tau=0.5)

Return the Payne model, where the speed of traffic relaxes to the Greenshields speed and anticipates density.

The density and the flow q = rho v on a ring road of length 1 follow rho_t + q_x = 0 and q_t + (q^2 / rho + P(rho))_x = (rho V(rho) - q) / τ, with V(rho) = vmax (1 - rho / rho_max) and the anticipation pressure P(rho) = vmax rho / (2 τ rho_max). This is the speed equation v_t + v v_x = (V - v) / τ + V'(rho) rho_x / (2 τ rho) in conservation form. The flow starts in equilibrium with a bump of density.

Parameters:

Name Type Description Default
n int

The number of cells.

 200
density float

The mean density.

 0.3
bump float

The relative height of the density bump.

 0.2
vmax float

The speed on an empty road.

 1.0
rho_max float

The jam density.

 1.0
tau float

The relaxation time of the speed.

 0.5

Returns:

Name Type Description
system System

The model, with the density road.rho and the flow road.q.

Examples:

>>> from skcomplex.models.engineering.traffic import payne
>>> payne(10).elements[0].values['tau']
0.5
Source code in src/skcomplex/models/engineering/traffic/_continuum.py 
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'field', 'conservation law', 'finite volume', 'pde'),
    references=('Payne 1971', 'Chowdhury, Santen, and Schadschneider 2000'),
    defaults={'until': 1.0, 'dt': 0.05, 'solver': SSPRKSolver()},
)
def payne(
    n: int = 200,
    *,
    density: float = 0.3,
    bump: float = 0.2,
    vmax: float = 1.0,
    rho_max: float = 1.0,
    tau: float = 0.5,
) -> System:
    """Return the Payne model, where the speed of traffic relaxes to the Greenshields speed and anticipates density.

    The density and the flow `q = rho v` on a ring road of length 1 follow `rho_t + q_x = 0` and
    `q_t + (q^2 / rho + P(rho))_x = (rho V(rho) - q) / τ`, with `V(rho) = vmax (1 - rho / rho_max)` and the
    anticipation pressure `P(rho) = vmax rho / (2 τ rho_max)`. This is the speed equation
    `v_t + v v_x = (V - v) / τ + V'(rho) rho_x / (2 τ rho)` in conservation form. The flow starts in equilibrium
    with a bump of density.

    Args:
        n:
            The number of cells.

        density:
            The mean density.

        bump:
            The relative height of the density bump.

        vmax:
            The speed on an empty road.

        rho_max:
            The jam density.

        tau:
            The relaxation time of the speed.

    Returns:
        system:
            The model, with the density `road.rho` and the flow `road.q`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import payne
        >>> payne(10).elements[0].values['tau']
        0.5
    """
    rho = _density_bump(n, density, bump)
    values = {'rho': rho, 'q': rho * vmax * (1 - rho / rho_max), 'vmax': vmax, 'rho_max': rho_max, 'tau': tau}
    road = FieldElement(
        'road',
        element=PayneElement,
        space=LatticeSpace(n, spacing=1 / n),
        values=values,
        discretization='finite-volume',
    )
    density = ConservationInteraction('density', target='road', state='rho', flux=lambda s: s.target.q)
    flow = ConservationInteraction(
        'flow',
        target='road',
        state='q',
        flux=lambda s: s.target.q**2 / s.target.rho
        + s.target.vmax * s.target.rho / (2 * s.target.tau * s.target.rho_max),
    )
    relaxation = RelaxationInteraction(
        'relaxation',
        target='road',
        state='q',
        goal=lambda s: s.target.rho * s.target.vmax * (1 - s.target.rho / s.target.rho_max),
        time=lambda s: s.target.tau,
    )
    return System(elements=[road], interactions=[density, flow, relaxation])

ridesharing(n_cabs=30, *, demand=1.0, width=30, height=30, base_fare=2.0, block_fare=0.5, random_state=0)

Return riders who request trips on a city grid and cabs that a stable matching sends to them.

Every step new riders appear at random cells with random destinations, and the dispatch matches the waiting riders with the idle cabs by a stable matching in which riders prefer nearer cabs and cabs prefer riders who waited longer. Every cab drives one block a step to its rider, collects the fare, and drives the rider to its destination. While the demand stays below the trips the fleet can serve, riders wait about as long as a cab takes to reach them. As the demand approaches that capacity, every cab is busy, the queue of waiting riders grows, and the waiting time rises steeply, as in a queue near saturation.

Parameters:

Name Type Description Default
n_cabs int

The number of cabs.

 30
demand float

The mean number of new riders in a step.

 1.0
width int

The number of columns of the grid.

 30
height int

The number of rows of the grid.

 30
base_fare float

The fare of every trip.

 2.0
block_fare float

The fare per block of a trip.

 0.5
random_state int | Generator | None

The seed or generator of the cells of the cabs.

 0

Returns:

Name Type Description
system System

The model, with the groups riders and cabs, the arrivals arrival, and the matching dispatch.

Examples:

>>> from skcomplex.models.engineering.traffic import ridesharing
>>> [group.size for group in ridesharing(5).elements]
[1, 5]
Source code in src/skcomplex/models/engineering/traffic/_ridesharing.py 
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
@register_model(
    field='engineering/traffic',
    tags=('ridesharing', 'agents', 'matching', 'queues', 'grid'),
    references=('Carmichael and Hadzikadic 2019 ch. 7',),
    defaults={'until': 200},
    observables=_ridesharing_observables,
)
def ridesharing(
    n_cabs: int = 30,
    *,
    demand: float = 1.0,
    width: int = 30,
    height: int = 30,
    base_fare: float = 2.0,
    block_fare: float = 0.5,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return riders who request trips on a city grid and cabs that a stable matching sends to them.

    Every step new riders appear at random cells with random destinations, and the dispatch matches the waiting
    riders with the idle cabs by a stable matching in which riders prefer nearer cabs and cabs prefer riders who waited
    longer. Every cab drives one block a step to its rider, collects the fare, and drives the rider to its destination.
    While the demand stays below the trips the fleet can serve, riders wait about as long as a cab takes to reach them.
    As the demand approaches that capacity, every cab is busy, the queue of waiting riders grows, and the waiting time
    rises steeply, as in a queue near saturation.

    Args:
        n_cabs:
            The number of cabs.

        demand:
            The mean number of new riders in a step.

        width:
            The number of columns of the grid.

        height:
            The number of rows of the grid.

        base_fare:
            The fare of every trip.

        block_fare:
            The fare per block of a trip.

        random_state:
            The seed or generator of the cells of the cabs.

    Returns:
        system:
            The model, with the groups `riders` and `cabs`, the arrivals `arrival`, and the matching `dispatch`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import ridesharing
        >>> [group.size for group in ridesharing(5).elements]
        [1, 5]
    """
    rng = np.random.default_rng(random_state)
    cells = rng.integers((width, height), size=(n_cabs + 1, 2))
    riders = GroupElement(
        'riders',
        element=RiderElement,
        size=1,
        values={'destination_x': float(rng.integers(width)), 'destination_y': float(rng.integers(height))},
        variable=True,
        location=[(int(cells[0, 0]), int(cells[0, 1]))],
    )
    cabs = GroupElement(
        'cabs',
        element=CabElement,
        size=n_cabs,
        values={'base_fare': base_fare, 'block_fare': block_fare},
        location=[(int(x), int(y)) for x, y in cells[1:]],
    )
    return System(
        elements=[riders, cabs],
        interactions=[
            ArrivalInteraction('arrival', target='riders', demand=demand, width=width, height=height),
            MapInteraction(
                'waiting',
                target='riders',
                state='waited',
                function=lambda agent, model: 1.0,
                mode='increment',
                form='agent',
                phase='drive',
            ),
            MatchingInteraction(
                'dispatch',
                requesters='riders',
                providers='cabs',
                rule=StableMatchingSpec(requester_preference=_nearness, provider_preference=_patience),
                quantity='seats',
                phase='dispatch',
                on_match=_dispatch,
            ),
        ],
        space=LatticeSpace((width, height), boundary='fixed', neighborhood='von_neumann'),
        time=DiscreteTime(
            schedule=PhasesScheduleSpec([PhaseSpec('arrival'), PhaseSpec('dispatch'), PhaseSpec('drive')]),
        ),
    )

ring_road_drivers(n_cars=30, length=1000.0, *, eps=0.02, speed_limit=40.0)

Return drivers on a circular highway, where small random errors of speed set off traffic jams.

The cars start at rest, equally spaced. Every driver accelerates by 1 per step up to the speed limit, the speed takes a random relative error of up to eps, and a speed that would reach the car ahead becomes 0.

Parameters:

Name Type Description Default
n_cars int

The number of cars.

 30
length float

The length of the ring road.

 1000.0
eps float

The relative amplitude of the random error of speed.

 0.02
speed_limit float

The largest speed.

 40.0

Returns:

Name Type Description
system System

The model, with the group drivers of speeds speed on a ring road of length length.

Examples:

>>> from skcomplex.models.engineering.traffic import ring_road_drivers
>>> ring_road_drivers(10).space.check().prepare().extent_[0]
1000.0
Source code in src/skcomplex/models/engineering/traffic/_drivers.py 
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'jams', 'agents', 'stochastic'),
    references=('Downey 2018', 'Resnick 1994'),
    defaults={'until': 100},
)
def ring_road_drivers(
    n_cars: int = 30,
    length: float = 1000.0,
    *,
    eps: float = 0.02,
    speed_limit: float = 40.0,
) -> System:
    """Return drivers on a circular highway, where small random errors of speed set off traffic jams.

    The cars start at rest, equally spaced. Every driver accelerates by 1 per step up to the speed limit, the speed
    takes a random relative error of up to `eps`, and a speed that would reach the car ahead becomes 0.

    Args:
        n_cars:
            The number of cars.

        length:
            The length of the ring road.

        eps:
            The relative amplitude of the random error of speed.

        speed_limit:
            The largest speed.

    Returns:
        system:
            The model, with the group `drivers` of speeds `speed` on a ring road of length `length`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import ring_road_drivers
        >>> ring_road_drivers(10).space.check().prepare().extent_[0]
        1000.0
    """
    positions = np.linspace(0.0, length, n_cars, endpoint=False)
    drivers = GroupElement(
        'drivers',
        element=DriverElement,
        size=n_cars,
        values={'rank': np.arange(n_cars, dtype=float), 'eps': eps, 'speed_limit': speed_limit},
        location=[(float(position), ROAD_WIDTH / 2) for position in positions],
    )
    return System(elements=[drivers], space=ContinuousSpace(extent=(length, ROAD_WIDTH), boundary='periodic'))

traffic_density_equation(n=100, *, length=50.0, density=0.5, gamma=1.0, u0=1.0, rho_max=1.0, D1=0.2, perturbation=0.01, random_state=0)

Return the density and flux equations of Langevin traffic on a ring road, whose uniform flow breaks into jams.

The density rho and the flux w of cars that relax to the desired speed u(rho) = u0 (1 - rho / rho_max) at the rate gamma, with the diffusion coefficient D(rho) = D1 rho, follow the moment equations rho_t + w_x = (D rho_x)_x and w_t + (2 u w - rho u^2)_x = -gamma (w - rho u) + (D w_x)_x, closed by expanding the velocity distribution around u(rho). A long density wave of the uniform flow grows when the effective diffusion coefficient D(rho) + ((rho u^2)' - (rho u)'^2) / gamma = D(rho) - rho^2 u'(rho)^2 / gamma is negative, so the flow is unstable above the density D1 gamma rho_max^2 / u0^2, and the perturbation grows into dense jams. The uniform density starts with independent Gaussian perturbations of relative size perturbation in every cell, at the flux rho u(rho).

Parameters:

Name Type Description Default
n int

The number of cells.

 100
length float

The length of the ring road.

 50.0
density float

The mean density.

 0.5
gamma float

The relaxation rate of the speed.

 1.0
u0 float

The desired speed on an empty road.

 1.0
rho_max float

The jam density, at which the desired speed vanishes.

 1.0
D1 float

The diffusion coefficient per unit density.

 0.2
perturbation float

The relative standard deviation of the initial density perturbations.

 0.01
random_state int | Generator | None

The seed or generator of the initial perturbations.

 0

Returns:

Name Type Description
system System

The model, with the density road.rho and the flux road.w.

Examples:

>>> from skcomplex.models.engineering.traffic import traffic_density_equation
>>> traffic_density_equation(10).elements[0].values['D1']
0.2
Source code in src/skcomplex/models/engineering/traffic/_langevin.py 
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
@register_model(
    field='engineering/traffic',
    tags=('traffic', 'field', 'jams', 'instability', 'pde'),
    references=('Bertin 2021',),
    defaults={'until': 400.0, 'dt': 4.0, 'solver': RK4Solver(dt=0.1)},
)
def traffic_density_equation(
    n: int = 100,
    *,
    length: float = 50.0,
    density: float = 0.5,
    gamma: float = 1.0,
    u0: float = 1.0,
    rho_max: float = 1.0,
    D1: float = 0.2,
    perturbation: float = 0.01,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return the density and flux equations of Langevin traffic on a ring road, whose uniform flow breaks into jams.

    The density `rho` and the flux `w` of cars that relax to the desired speed `u(rho) = u0 (1 - rho / rho_max)` at
    the rate `gamma`, with the diffusion coefficient `D(rho) = D1 rho`, follow the moment equations
    `rho_t + w_x = (D rho_x)_x` and `w_t + (2 u w - rho u^2)_x = -gamma (w - rho u) + (D w_x)_x`, closed by expanding
    the velocity distribution around `u(rho)`. A long density wave of the uniform flow grows when the effective
    diffusion coefficient `D(rho) + ((rho u^2)' - (rho u)'^2) / gamma = D(rho) - rho^2 u'(rho)^2 / gamma` is negative,
    so the flow is unstable above the density `D1 gamma rho_max^2 / u0^2`, and the perturbation grows into dense jams.
    The uniform density starts with independent Gaussian perturbations of relative size
    `perturbation` in every cell, at the flux `rho u(rho)`.

    Args:
        n:
            The number of cells.

        length:
            The length of the ring road.

        density:
            The mean density.

        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 diffusion coefficient per unit density.

        perturbation:
            The relative standard deviation of the initial density perturbations.

        random_state:
            The seed or generator of the initial perturbations.

    Returns:
        system:
            The model, with the density `road.rho` and the flux `road.w`.

    Examples:
        >>> from skcomplex.models.engineering.traffic import traffic_density_equation
        >>> traffic_density_equation(10).elements[0].values['D1']
        0.2
    """
    rng = np.random.default_rng(random_state)
    rho = density * (1 + perturbation * rng.standard_normal(n))
    values = {
        'rho': rho,
        'w': rho * u0 * (1 - rho / rho_max),
        'gamma': gamma,
        'u0': u0,
        'rho_max': rho_max,
        'D1': D1,
    }
    space = LatticeSpace(n, spacing=length / n)
    road = FieldElement('road', element=TrafficMomentsElement, space=space, values=values)
    return System(elements=[road])