VRP Route Optimization Algorithms

Batch Recalculation on Demand-Spike Days with OR-Tools

Re-dispatch on storm and seasonal surges by re-solving a fresh demand vector while reusing the cached matrix.

This page solves one operational emergency: aggregate demand jumps mid-shift — post-holiday overflow, storm debris, a seasonal surge — and dispatch must re-plan the fleet in seconds by re-solving against the new demand, not by rebuilding the whole model from cold geometry.

When the day’s plan was computed at 5 a.m. and reality diverges by lunchtime, the naive fix is to rerun the entire pipeline: rebuild the distance matrix, re-register every dimension, solve from scratch. That is wasteful and slow — the geometry has not changed, only the load at known stops. The disciplined recalculation reuses the cached matrix, re-resolves demand through the same live-then-baseline ladder the VRP Route Optimization Algorithms pipeline already relies on, warm-starts the solver from the plan currently on the road, adds capacity only if the spike demands it, and re-solves under a hard time bound. This is the on-demand recalculation path referenced by Dynamic Threshold Tuning: the same calibrated solver, re-invoked against a moved demand vector rather than a rebuilt problem.

Environment & Data Prerequisites

The recalculation needs only the pinned solver and the standard library:

python -m pip install \
  ortools==9.11.4210 \
  pytest==8.2.0

Two inputs drive a spike-day re-solve: the baseline demand vector the morning plan was built on, and the live delta arriving from the field. Demand stays in whole kilograms as int, and the node set is held fixed so the cached matrix stays valid — a spike raises the load at existing bins rather than moving them.

Field Type Unit / range Role in recalculation
stop_id str e.g. BIN_04 cache key + audit join
lat / lon float WGS84 bounds matrix key; unchanged on a spike
baseline_kg int 0 … capacity morning demand the live plan was solved on
live_kg int or None 0 … capacity current sensor reading; None falls to baseline
spike_threshold float e.g. 1.20 aggregate ratio that triggers a re-solve

Step-by-Step Implementation

Step 1 — A matrix cache that refuses to rebuild geometry

The distance matrix is keyed by the coordinate set, so as long as the stops are unchanged a re-solve is a cache hit, not a rebuild. A build counter makes the reuse auditable — and testable.

from __future__ import annotations
from dataclasses import dataclass
import math
import hashlib
import json


@dataclass(frozen=True)
class Stop:
    stop_id: str
    lat: float
    lon: float


class DistanceMatrixCache:
    """Coordinate-keyed matrix cache. builds increments only on a genuine miss."""

    def __init__(self) -> None:
        self._cache: dict[str, list[list[int]]] = {}
        self.builds = 0

    @staticmethod
    def _key(stops: list[Stop]) -> str:
        raw = json.dumps([(s.stop_id, s.lat, s.lon) for s in stops], sort_keys=True)
        return hashlib.sha256(raw.encode()).hexdigest()

    def get(self, stops: list[Stop]) -> tuple[list[list[int]], str]:
        key = self._key(stops)
        if key not in self._cache:
            self.builds += 1                 # a cold build — the thing we avoid on a spike
            self._cache[key] = self._build(stops)
        return self._cache[key], key

    @staticmethod
    def _build(stops: list[Stop]) -> list[list[int]]:
        n = len(stops)
        d = [[0] * n for _ in range(n)]
        for i in range(n):
            for j in range(n):
                if i == j:
                    continue
                mean_lat = math.radians((stops[i].lat + stops[j].lat) / 2)
                dlat = stops[i].lat - stops[j].lat
                dlon = (stops[i].lon - stops[j].lon) * math.cos(mean_lat)
                d[i][j] = int(111_320 * math.hypot(dlat, dlon))
        return d


BASE_STOPS: list[Stop] = [
    Stop("DEPOT",  40.7128, -74.0060),
    Stop("BIN_01", 40.7180, -74.0020),
    Stop("BIN_02", 40.7205, -74.0110),
    Stop("BIN_03", 40.7090, -74.0150),
    Stop("BIN_04", 40.7060, -74.0030),
    Stop("BIN_05", 40.7220, -73.9950),
    Stop("BIN_06", 40.7040, -74.0120),
]
BASELINE_DEMAND = {"BIN_01": 120, "BIN_02": 140, "BIN_03": 130,
                   "BIN_04": 110, "BIN_05": 150, "BIN_06": 125}
CAPACITY_KG = 600

Step 2 — Re-resolve demand through the live-then-baseline ladder

A spike does not mean every sensor reported. Demand is resolved per stop: use the live reading when present and non-negative, otherwise fall back to the morning baseline — never silently substitute zero, which would skip an overflowing bin. This mirrors the fallback ladder formalized in Fallback Routing Logic.

def resolve_demand(stops: list[Stop],
                   baseline: dict[str, int],
                   live: dict[str, int | None]) -> list[int]:
    """Index-aligned demand vector: live reading, else baseline. Depot is 0."""
    vec = []
    for s in stops:
        if s.stop_id == "DEPOT":
            vec.append(0)
            continue
        reading = live.get(s.stop_id)
        if reading is not None and reading >= 0:
            vec.append(int(reading))
        else:
            vec.append(int(baseline[s.stop_id]))
    return vec

Step 3 — Detect the spike trigger

The trigger is an aggregate ratio: total re-resolved demand over the baseline the live plan was solved on. Below the threshold, the morning plan stands and no re-solve fires; at or above it, recalculation is warranted.

def detect_spike(baseline_vec: list[int], live_vec: list[int],
                 threshold: float = 1.20) -> dict:
    agg_base = sum(baseline_vec)
    agg_live = sum(live_vec)
    ratio = agg_live / agg_base if agg_base else float("inf")
    return {"triggered": ratio >= threshold,
            "ratio": round(ratio, 3),
            "aggregate_baseline_kg": agg_base,
            "aggregate_live_kg": agg_live}

The number of trucks the spike needs is the ceiling of total demand over per-vehicle capacity — this is what turns a load surge into an added vehicle rather than an overloaded run:

Vspike=max⁡ ⁣(Vbase, ⌈∑idiQ⌉)V_{\text{spike}} = \max\!\left(V_{\text{base}},\ \left\lceil \frac{\sum_i d_i}{Q} \right\rceil\right)

where did_i is the re-resolved demand at stop ii and QQ is the per-vehicle capacity in kilograms.

Step 4 — Reuse the matrix and warm-start the re-solve

The re-solve pulls the matrix from cache (no rebuild), sizes the fleet from the formula above, and warm-starts from the plan currently on the road. ReadAssignmentFromRoutes seeds the search with the live routes; added trucks get empty routes. A hard time_limit bounds the recalculation so a spike-day re-solve cannot itself become an outage.

from ortools.constraint_solver import pywrapcp, routing_enums_pb2
import math as _math


def _params(time_limit_s: int) -> pywrapcp.DefaultRoutingSearchParameters:
    p = pywrapcp.DefaultRoutingSearchParameters()
    p.first_solution_strategy = routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC
    p.local_search_metaheuristic = routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH
    p.time_limit.FromSeconds(time_limit_s)
    return p


def _model(stops, dist, demand_vec, num_vehicles, capacity):
    manager = pywrapcp.RoutingIndexManager(len(stops), num_vehicles, 0)
    routing = pywrapcp.RoutingModel(manager)
    tc = routing.RegisterTransitCallback(
        lambda a, b: dist[manager.IndexToNode(a)][manager.IndexToNode(b)])
    routing.SetArcCostEvaluatorOfAllVehicles(tc)
    dcb = routing.RegisterUnaryTransitCallback(
        lambda a: demand_vec[manager.IndexToNode(a)])
    routing.AddDimensionWithVehicleCapacity(
        dcb, 0, [capacity] * num_vehicles, True, "Payload")
    return routing, manager


def _routes(routing, manager, sol, stops):
    out = []
    for v in range(manager.GetNumberOfVehicles()):
        idx, r = routing.Start(v), []
        while not routing.IsEnd(idx):
            n = manager.IndexToNode(idx)
            if n != 0:
                r.append(stops[n].stop_id)
            idx = sol.Value(routing.NextVar(idx))
        out.append(r)
    return out


def warm_resolve(stops, dist, demand_vec, num_vehicles, capacity,
                 prior_routes: list[list[str]], time_limit_s: int = 3):
    routing, manager = _model(stops, dist, demand_vec, num_vehicles, capacity)
    id_to_node = {s.stop_id: i for i, s in enumerate(stops)}
    seed = [[id_to_node[sid] for sid in r] for r in prior_routes]
    while len(seed) < num_vehicles:          # added trucks start empty
        seed.append([])
    routing.CloseModelWithParameters(_params(time_limit_s))
    initial = routing.ReadAssignmentFromRoutes(seed, True)
    sol = routing.SolveFromAssignmentWithParameters(initial, _params(time_limit_s))
    if sol is None:
        return None
    return {"objective": sol.ObjectiveValue(),
            "routes": _routes(routing, manager, sol, stops),
            "vehicles": num_vehicles}

Step 5 — Emit the delta-dispatch audit record

Every recalculation writes one structured record naming the trigger, the capacity change, and the objective delta, so an auditor can reconstruct why the road plan changed at 12:14. The JSONFormatter is the site-wide structured-logging pattern, with fields distinct to a re-dispatch event.

import logging


class JSONFormatter(logging.Formatter):
    def format(self, record: logging.LogRecord) -> str:
        payload = {
            "event": record.getMessage(),
            "trigger_reason": getattr(record, "trigger_reason", None),
            "spike_ratio": getattr(record, "spike_ratio", None),
            "vehicles_before": getattr(record, "vehicles_before", None),
            "vehicles_after": getattr(record, "vehicles_after", None),
            "delta_objective": getattr(record, "delta_objective", None),
            "matrix_cache_reused": getattr(record, "matrix_cache_reused", None),
        }
        return json.dumps(payload, separators=(",", ":"))


logger = logging.getLogger("waste_ops.recalc")
if not logger.handlers:
    _h = logging.StreamHandler()
    _h.setFormatter(JSONFormatter())
    logger.addHandler(_h)
logger.setLevel(logging.INFO)


def recalculate(cache: DistanceMatrixCache, stops, baseline, live,
                base_plan: dict, capacity: int, threshold: float = 1.20) -> dict:
    dist, key = cache.get(stops)              # reuse — no rebuild on a spike
    baseline_vec = resolve_demand(stops, baseline, {})   # all-baseline reference
    live_vec = resolve_demand(stops, baseline, live)
    spike = detect_spike(baseline_vec, live_vec, threshold)
    if not spike["triggered"]:
        return {"recalculated": False, **spike}

    needed = _math.ceil(sum(live_vec) / capacity)
    vehicles_after = max(base_plan["vehicles"], needed)
    plan = warm_resolve(stops, dist, live_vec, vehicles_after, capacity,
                        base_plan["routes"])
    record = {
        "recalculated": True,
        "trigger_reason": f"aggregate demand ratio {spike['ratio']} >= {threshold}",
        "spike_ratio": spike["ratio"],
        "vehicles_before": base_plan["vehicles"],
        "vehicles_after": vehicles_after,
        "delta_objective": plan["objective"] - base_plan["objective"],
        "matrix_cache_reused": True,
        "routes": plan["routes"],
    }
    logger.info("demand_spike_recalculated", extra=record)
    return record

Compliance Output

A triggered recalculation prints one delta-dispatch record. On a 1.6× storm-debris spike the six-bin baseline (two trucks) needs a third, and the matrix is served from cache:

{
  "recalculated": true,
  "trigger_reason": "aggregate demand ratio 1.6 >= 1.2",
  "spike_ratio": 1.6,
  "vehicles_before": 2,
  "vehicles_after": 3,
  "delta_objective": 977,
  "matrix_cache_reused": true,
  "routes": [["BIN_02"], ["BIN_04", "BIN_06", "BIN_03"], ["BIN_01", "BIN_05"]]
}

Each field carries an operational or audit purpose:

  • trigger_reason — the human-readable justification for a mid-shift re-dispatch, so the change is defensible rather than mysterious; it pairs with the threshold logic in Dynamic Threshold Tuning.
  • spike_ratio — the aggregate demand multiple that crossed the threshold; the primary signal distinguishing a genuine surge from sensor noise.
  • vehicles_before / vehicles_after — the capacity change, proving the spike was absorbed by adding a truck rather than overloading the existing fleet past its registered GVWR.
  • delta_objective — the added routing cost the surge imposed, reconciled later against fuel and hours-of-service accounting.
  • matrix_cache_reused — the record that geometry was not rebuilt, which is both the performance guarantee and the evidence that the re-solve used the same audited distances as the morning plan.

The re-resolved live readings themselves arrive pre-validated from the Telematics & Sensor Data Ingestion pipeline, so a spike is never triggered by a corrupt fix.

Verification

The tests prove the two properties that define a correct spike-day recalculation: a genuine surge adds capacity, and the distance matrix is reused rather than rebuilt.

def _base_plan(cache):
    dist, _ = cache.get(BASE_STOPS)
    vec = resolve_demand(BASE_STOPS, BASELINE_DEMAND, {})
    # cold morning solve at two vehicles
    routing, manager = _model(BASE_STOPS, dist, vec, 2, CAPACITY_KG)
    sol = routing.SolveWithParameters(_params(5))
    return {"objective": sol.ObjectiveValue(),
            "routes": _routes(routing, manager, sol, BASE_STOPS),
            "vehicles": 2}


def test_spike_adds_capacity_and_reuses_matrix():
    cache = DistanceMatrixCache()
    base = _base_plan(cache)
    assert cache.builds == 1                 # exactly one cold build so far

    live = {sid: int(kg * 1.6) for sid, kg in BASELINE_DEMAND.items()}
    result = recalculate(cache, BASE_STOPS, BASELINE_DEMAND, live, base, CAPACITY_KG)

    assert result["recalculated"] is True
    assert result["vehicles_after"] > result["vehicles_before"]   # capacity added
    assert result["matrix_cache_reused"] is True
    assert cache.builds == 1                 # NOT rebuilt during recalculation


def test_small_fluctuation_does_not_trigger():
    cache = DistanceMatrixCache()
    base = _base_plan(cache)
    live = {sid: int(kg * 1.05) for sid, kg in BASELINE_DEMAND.items()}
    result = recalculate(cache, BASE_STOPS, BASELINE_DEMAND, live, base, CAPACITY_KG)
    assert result["recalculated"] is False   # 1.05 < 1.20 threshold
    assert cache.builds == 1

Common Errors

Every recalculation rebuilds the distance matrix and the re-solve takes as long as a cold start. Root cause: constructing the matrix inside the re-solve instead of reading it from a keyed cache, so the O(n2)O(n^2) geometry build runs again even though no stop moved. Fix: route every solve through DistanceMatrixCache.get, key it on the coordinate set, and assert cache.builds stays flat across a spike, exactly as the verification does.

A spike-day re-solve runs for minutes and blows the dispatch deadline. Root cause: the warm re-solve was called without a time_limit, so guided local search keeps refining indefinitely — a recalculation meant to take three seconds becomes its own incident. Fix: always bound the re-solve with p.time_limit.FromSeconds(...) as _params does, and treat the warm-started first solution as a valid fallback if refinement has not converged in time.

ReadAssignmentFromRoutes returns None or the solve raises on a mismatched vector. Root cause: the demand vector length no longer matches the node count — usually because debris stops were appended to the stop list but the demand vector (or the warm-start routes) was not extended in lockstep. Fix: build the demand vector from the same stops list every time via resolve_demand, and if you genuinely add stops, extend the matrix with only the new rows and columns and rebuild the warm-start routes rather than reusing the stale plan.

FAQ

Why reuse the matrix instead of just rebuilding — it is only a few milliseconds here?

At six stops the rebuild is trivial, but a metropolitan network is thousands of stops and an n2n^2 matrix is tens of millions of cells; rebuilding it mid-shift under a spike is exactly when you have the least time to spare. Reuse also preserves auditability: the re-dispatched plan is provably scored against the same distances as the morning plan, so the delta_objective is a true comparison rather than an artifact of a regenerated matrix.

Does warm-starting bias the solver toward the old plan?

It seeds the search, it does not pin it. ReadAssignmentFromRoutes gives guided local search a feasible starting point built from the routes currently on the road, which is far better than a cold first solution when only demand has shifted. The solver is free to move stops between trucks — including onto the newly added vehicle — and the bounded re-solve converges faster because it starts near a good plan.

Up: Dynamic Threshold Tuning