[FeatReq] duration cost

.github/workflows/CI.yml

@@ -93,12 +93,12 @@ jobs:
         run: |
           uv run pytest
           uv run ninja coverage-xml -C build
-      - uses: codecov/codecov-action@v5
+      - uses: actions/upload-artifact@v4
         with:
-          fail_ci_if_error: true
-          token: ${{ secrets.CODECOV_TOKEN }}
-          plugins: pycoverage
-          files: build/meson-logs/coverage.xml
+          name: coverage-${{ matrix.compiler }}-${{ matrix.compiler-version }}-py${{ matrix.python-version }}
+          path: |
+            coverage.xml
+            build/meson-logs/coverage.xml
       - if: matrix.compiler == 'gcc'
         name: Install Valgrind
         run: sudo apt-get install -y valgrind

.github/workflows/CodSpeed.yml

@@ -53,9 +53,13 @@ jobs:
         run: |
           uv sync
           uv run buildtools/build_extensions.py --build_type debugoptimized --clean
-      - name: Run benchmarks
+      - name: Run benchmarks with CodSpeed
+        if: ${{ secrets.CODSPEED_TOKEN != '' }}
         # We evaluate the microbenchmarks using the release build.
         uses: CodSpeedHQ/action@v3
         with:
           token: ${{ secrets.CODSPEED_TOKEN }}
           run: uv run pytest benchmarks/ --codspeed
+      - name: Run benchmarks without CodSpeed upload
+        if: ${{ secrets.CODSPEED_TOKEN == '' }}
+        run: uv run pytest benchmarks/

notebooks/duration_constraints.ipynb

@@ -81,7 +81,11 @@
    "outputs": [],
    "source": [
     "m = pyvrp.Model()\n",
-    "m.add_vehicle_type(2, unit_distance_cost=0, unit_duration_cost=1)\n",
+    "m.add_vehicle_type(\n",
+    "    2,\n",
+    "    unit_distance_cost=0,\n",
+    "    duration_cost=pyvrp.PiecewiseLinearFunction([], [(0, 1)]),\n",
+    ")\n",
     "\n",
     "m.add_depot(\n",
     "    location=m.add_location(x=COORDS[0][0], y=COORDS[0][1]),\n",
@@ -214,7 +218,11 @@
    "outputs": [],
    "source": [
     "m = pyvrp.Model()\n",
-    "m.add_vehicle_type(2, unit_distance_cost=0, unit_duration_cost=1)\n",
+    "m.add_vehicle_type(\n",
+    "    2,\n",
+    "    unit_distance_cost=0,\n",
+    "    duration_cost=pyvrp.PiecewiseLinearFunction([], [(0, 1)]),\n",
+    ")\n",
     "\n",
     "m.add_depot(\n",
     "    location=m.add_location(x=COORDS[0][0], y=COORDS[0][1]),\n",
@@ -269,7 +277,10 @@
     "PyVRP also supports this.\n",
     "\n",
     "Let's consider a small example where each vehicle *must* start at $t=0$, and must return to the depot by $t=30$.\n",
-    "Additionally, a regular shift takes at most 20 time units, but a maximum overtime of 10 units is allowed, at additional cost."
+    "Additionally, a regular shift takes at most 20 time units, but a maximum overtime of 10 units is allowed, at additional cost.\n",
+    "\n",
+    "The duration cost is expressed as a piecewise linear function of route duration.\n",
+    "Here, the cost increases linearly at rate 1 for durations up to the nominal shift of 20, and then at the higher rate of 6 (= 1 + 5 overtime surcharge) for any overtime beyond that."
    ]
   },
   {
@@ -287,8 +298,8 @@
     "    shift_duration=20,  # nominal shift duration\n",
     "    max_overtime=10,  # maximum allowed overtime\n",
     "    unit_distance_cost=0,\n",
-    "    unit_duration_cost=1,\n",
-    "    unit_overtime_cost=5,  # additional unit overtime cost\n",
+    "    # slope=1 for d <= 20 (regular shift), slope=6 (=1+5 surcharge) for d > 20\n",
+    "    duration_cost=pyvrp.PiecewiseLinearFunction([20], [(0, 1), (-100, 6)]),\n",
     ")\n",
     "\n",
     "m.add_depot(\n",

pyvrp/Model.py

@@ -10,6 +10,7 @@
     ClientGroup,
     Depot,
     Location,
+    PiecewiseLinearFunction,
     ProblemData,
     Solution,
     VehicleType,
@@ -382,14 +383,15 @@ def add_vehicle_type(
         shift_duration: int = np.iinfo(np.int64).max,
         max_distance: int = np.iinfo(np.int64).max,
         unit_distance_cost: int = 1,
-        unit_duration_cost: int = 0,
         profile: Profile | None = None,
         start_late: int | None = None,
         initial_load: int | list[int] = [],
         reload_depots: list[Depot] = [],
         max_reloads: int = np.iinfo(np.uint64).max,
         max_overtime: int = 0,
-        unit_overtime_cost: int = 0,
+        duration_cost: PiecewiseLinearFunction = PiecewiseLinearFunction(
+            [(np.int64(0), np.int64(0)), (np.int64(1), np.int64(0))]
+        ),
         *,
         name: str = "",
     ) -> VehicleType:
@@ -453,14 +455,13 @@ def add_vehicle_type(
             shift_duration=shift_duration,
             max_distance=max_distance,
             unit_distance_cost=unit_distance_cost,
-            unit_duration_cost=unit_duration_cost,
             profile=profile_idx,
             start_late=start_late,
             initial_load=init_load,
             reload_depots=reloads,
             max_reloads=max_reloads,
             max_overtime=max_overtime,
-            unit_overtime_cost=unit_overtime_cost,
+            duration_cost=duration_cost,
             name=name,
         )
 

pyvrp/_pyvrp.pyi

@@ -76,6 +76,7 @@ class PiecewiseLinearFunction:
     @property
     def segments(self) -> list[tuple[np.int64, np.int64]]: ...
     def is_monotonically_increasing(self) -> bool: ...
+    def is_non_negative(self, lb: np.int64) -> bool: ...
     def __eq__(self, other: object) -> bool: ...
     def __getstate__(self) -> tuple: ...
     def __setstate__(self, state: tuple, /) -> None: ...
@@ -170,14 +171,13 @@ class VehicleType:
     max_distance: int
     fixed_cost: int
     unit_distance_cost: int
-    unit_duration_cost: int
     profile: int
     start_late: int
     initial_load: list[int]
     reload_depots: list[int]
     max_reloads: int
     max_overtime: int
-    unit_overtime_cost: int
+    duration_cost: PiecewiseLinearFunction
     max_duration: int
     name: str
     def __init__(
@@ -192,14 +192,13 @@ class VehicleType:
         shift_duration: int = ...,
         max_distance: int = ...,
         unit_distance_cost: int = 1,
-        unit_duration_cost: int = 0,
         profile: int = 0,
         start_late: int | None = None,
         initial_load: list[int] = [],
         reload_depots: list[int] = [],
         max_reloads: int = ...,
         max_overtime: int = 0,
-        unit_overtime_cost: int = 0,
+        duration_cost: PiecewiseLinearFunction = ...,
         *,
         name: str = "",
     ) -> None: ...
@@ -217,14 +216,13 @@ class VehicleType:
         shift_duration: int | None = None,
         max_distance: int | None = None,
         unit_distance_cost: int | None = None,
-        unit_duration_cost: int | None = None,
         profile: int | None = None,
         start_late: int | None = None,
         initial_load: list[int] | None = None,
         reload_depots: list[int] | None = None,
         max_reloads: int | None = None,
         max_overtime: int | None = None,
-        unit_overtime_cost: int | None = None,
+        duration_cost: PiecewiseLinearFunction | None = None,
         *,
         name: str | None = None,
     ) -> VehicleType: ...

pyvrp/cpp/PiecewiseLinearFunction.h

@@ -51,7 +51,7 @@ namespace pyvrp
 template <typename Dom, typename Co> class PiecewiseLinearFunction
 {
 public:
-    using Segment = std::pair<Dom, Dom>;
+    using Segment = std::pair<Co, Co>;
     using Point = std::pair<Dom, Co>;
 
 private:
@@ -86,7 +86,13 @@ template <typename Dom, typename Co> class PiecewiseLinearFunction
      */
     [[nodiscard]] bool isMonotonicallyIncreasing() const;
 
+    /**
+     * Returns whether this function is non-negative for all :math:`x \ge lb`.
+     */
+    [[nodiscard]] bool isNonNegative(Dom lb) const;
+
     bool operator==(PiecewiseLinearFunction const &other) const = default;
+    auto operator<=>(PiecewiseLinearFunction const &other) const = default;
 };
 
 template <typename Dom, typename Co>
@@ -112,8 +118,8 @@ PiecewiseLinearFunction<Dom, Co>::PiecewiseLinearFunction(
         if (curr.first == next.first)  // nothing to do; we have a jump at this
             continue;                  // point, but not a new segment.
 
-        auto const dy = next.second - curr.second;
-        auto const dx = next.first - curr.first;
+        auto const dy = static_cast<int64_t>(next.second - curr.second);
+        auto const dx = static_cast<int64_t>(next.first - curr.first);
 
         if (dx < 0)
             throw std::invalid_argument("Points must be non-decreasing in x.");
@@ -124,8 +130,9 @@ PiecewiseLinearFunction<Dom, Co>::PiecewiseLinearFunction(
         if (idx != 0)  // breakpoints separate segments
             breakpoints_.push_back(curr.first);
 
-        auto const slope = dy / dx;
-        auto const intercept = curr.second - slope * curr.first;
+        auto const slope = static_cast<Co>(dy / dx);
+        auto const intercept
+            = curr.second - slope * static_cast<Co>(curr.first);
         segments_.emplace_back(intercept, slope);
     }
 
@@ -158,7 +165,7 @@ Co PiecewiseLinearFunction<Dom, Co>::operator()(Dom x) const
         std::upper_bound(breakpoints_.begin(), breakpoints_.end(), x));
 
     auto const [intercept, slope] = segments_[idx];
-    return static_cast<Co>(intercept + slope * x);
+    return intercept + slope * static_cast<Co>(x);
 }
 
 template <typename Dom, typename Co>
@@ -187,15 +194,54 @@ bool PiecewiseLinearFunction<Dom, Co>::isMonotonicallyIncreasing() const
         auto const [prevIntercept, prevSlope] = segments_[idx];
         auto const [nextIntercept, nextSlope] = segments_[idx + 1];
 
-        auto const left = prevIntercept + prevSlope * breakpoint;
-        auto const right = nextIntercept + nextSlope * breakpoint;
+        auto const left
+            = prevIntercept + prevSlope * static_cast<Co>(breakpoint);
+        auto const right
+            = nextIntercept + nextSlope * static_cast<Co>(breakpoint);
 
         if (right < left)
             return false;
     }
 
     return true;
 }
+
+template <typename Dom, typename Co>
+bool PiecewiseLinearFunction<Dom, Co>::isNonNegative(Dom lb) const
+{
+    // The last segment extends to +inf. A negative slope would drive the
+    // function to -inf, violating non-negativity.
+    if (segments_.back().second < 0)
+        return false;
+
+    // Check the value at the lower bound of the domain.
+    if ((*this)(lb) < 0)
+        return false;
+
+    // Check all breakpoints strictly greater than lb. Since each segment is
+    // linear, its minimum over any interval is at one of the endpoints; those
+    // endpoints are lb (already checked) and the breakpoints below.
+    for (size_t idx = 0; idx != breakpoints_.size(); ++idx)
+    {
+        auto const breakpoint = breakpoints_[idx];
+
+        if (breakpoint <= lb)
+            continue;
+
+        auto const [prevIntercept, prevSlope] = segments_[idx];
+        auto const [nextIntercept, nextSlope] = segments_[idx + 1];
+
+        auto const left
+            = prevIntercept + prevSlope * static_cast<Co>(breakpoint);
+        auto const right
+            = nextIntercept + nextSlope * static_cast<Co>(breakpoint);
+
+        if (left < 0 || right < 0)
+            return false;
+    }
+
+    return true;
+}
 }  // namespace pyvrp
 
 #endif  // PYVRP_PIECEWISELINEARFUNCTION_H