Piecewise Linear Constraints#

Piecewise linear (PWL) constraints approximate nonlinear functions as connected linear pieces, allowing you to model cost curves, efficiency curves, or production functions within a linear programming framework.

Terminology used in this page:

breakpoint — an \((x, y)\) knot where the slope can change.
piece — a linear part between two adjacent breakpoints on a single connected curve. n breakpoints define n − 1 pieces.
segment — a disjoint operating region in the disjunctive formulation, built via the segments() factory. Within one segment the curve is itself piecewise-linear (made of pieces); between segments there are gaps.

Quick Start #

Equality — lock variables onto the piecewise curve:

import linopy

m = linopy.Model()
power = m.add_variables(name="power", lower=0, upper=100)
fuel = m.add_variables(name="fuel")

# fuel = f(power) on the piecewise curve defined by these breakpoints
m.add_piecewise_formulation(
    (power, [0, 30, 60, 100]),
    (fuel, [0, 36, 84, 170]),
)

Inequality — bound one expression by the curve:

# fuel <= f(power).  "auto" picks the cheapest correct formulation
# (pure LP with chord constraints when the curve's curvature matches
# the requested sign; SOS2/incremental otherwise).
m.add_piecewise_formulation(
    (fuel, [0, 20, 30, 35], "<="),  # bounded by the curve
    (power, [0, 10, 20, 30]),  # pinned to the curve
)

Each (expression, breakpoints[, sign]) tuple pairs a variable with its breakpoint values, and optionally marks it as bounded by the curve ("<=" or ">=") instead of pinned to it. All tuples share interpolation weights, so at any feasible point every variable corresponds to the same point on the piecewise curve.

API #

`add_piecewise_formulation`#

m.add_piecewise_formulation(
    (expr1, breakpoints1),  # pinned (sign defaults to "==")
    (expr2, breakpoints2, "<="),  # or with an explicit sign
    ...,
    method="auto",  # "auto", "sos2", "incremental", or "lp"
    active=None,  # binary variable to gate the constraint
    name=None,  # base name for generated variables/constraints
)

Creates auxiliary variables and constraints that enforce either a joint equality (all tuples on the curve, the default) or a one-sided bound (at most one tuple bounded by the curve, the rest pinned).

`breakpoints` and `segments`#

Two factories with distinct geometric meaning:

breakpoints() — values along a single connected curve. Linear pieces between adjacent breakpoints are interpolated continuously.
segments() — disjoint operating regions with gaps between them (e.g. forbidden zones). Builds a 2-D array consumed by the disjunctive formulation, where exactly one region is active at a time.

linopy.breakpoints([0, 50, 100])  # connected
linopy.breakpoints({"gen1": [0, 50], "gen2": [0, 80]}, dim="gen")  # per-entity
linopy.Slopes(
    [1.2, 1.4], y0=0
)  # from slopes (deferred — pairs with a sibling tuple)
linopy.segments([(0, 10), (50, 100)])  # two disjoint regions
linopy.segments({"gen1": [(0, 10)], "gen2": [(0, 80)]}, dim="gen")

Per-tuple sign — equality vs inequality #

By default each tuple’s expression is pinned to the piecewise curve. Pass a third tuple element ("<=" or ">=") to mark a single expression as bounded by the curve — it can undershoot ("<=") or overshoot (">=") the interpolated value, while every other tuple stays pinned.

# Joint equality (default): both expressions on the curve.
m.add_piecewise_formulation((y, y_pts), (x, x_pts))

# Bounded above: y <= f(x), x pinned.
m.add_piecewise_formulation((y, y_pts, "<="), (x, x_pts))

# Bounded below: y >= f(x), x pinned.
m.add_piecewise_formulation((y, y_pts, ">="), (x, x_pts))

# 3-variable equality (CHP heat/power/fuel): all three on one curve.
m.add_piecewise_formulation((power, p_pts), (fuel, f_pts), (heat, h_pts))

Restrictions (current):

At most one tuple may carry a non-equality sign — a single bounded side.
With 3 or more tuples, all signs must be "==".

Multi-bounded and N≥3-inequality use cases aren’t supported yet. If you have a concrete use case, please open an issue at PyPSA/linopy#issues so we can scope it properly.

Formulation. For methods that introduce shared interpolation weights (SOS2 and incremental — see below), only the link constraint between the weights and the bounded expression changes. Pinned tuples \(j\) keep the equality, and the bounded tuple \(b\) flips to the requested sign:

\[ \begin{align}\begin{aligned}&e_j = \sum_{i=0}^{n} \lambda_i \, B_{j,i} \quad \text{(pinned, } j \ne b \text{)}\\&e_b \ \text{sign}\ \sum_{i=0}^{n} \lambda_i \, B_{b,i} \quad \text{(bounded)}\end{aligned}\end{align} \]

Internally this shows up as a stacked *_link equality covering the pinned tuples plus a separate signed *_output_link for the bounded tuple. The method="lp" path encodes the same one-sided semantics without weights — see the LP section below.

Geometry. For 2 variables with sign="<=" on a concave curve \(f\), the feasible region is the hypograph of \(f\) on its domain:

\[\{ (x, y) \ :\ x_0 \le x \le x_n,\ y \le f(x) \}.\]

For convex \(f\) with sign=">=" it is the epigraph. Mismatched sign + curvature (convex + "<=", or concave + ">=") describes a non-convex region — method="auto" falls back to SOS2/incremental and method="lp" raises.

Choice of bounded tuple. The bounded tuple should correspond to a quantity with a mechanism for below-curve operation — typically a controllable dissipation path: heat rejection via cooling tower (also called thermal curtailment), electrical curtailment, or emissions after post-treatment. Marking a consumption-side variable such as fuel intake as bounded yields a valid but loose formulation: the characteristic curve fixes fuel draw at a given load, so "<=" on fuel admits operating points the plant cannot physically realise. An objective that rewards lower fuel may then find a non-physical optimum — safe only when no such objective pressure exists.

When is a one-sided bound wanted?

For continuous curves, the main reason to reach for "<=" / ">=" is to unlock the LP chord formulation — no SOS2, no binaries, just pure LP. On a convex/concave curve with a matching sign, the chord inequalities are as tight as SOS2, so you get the same optimum with a cheaper model. Inequality formulations also tighten the LP relaxation of SOS2/incremental, which can reduce branch-and-bound work even when LP itself is not applicable.

For disjunctive curves (segments(...)), the per-tuple sign is a first-class tool in its own right: disconnected operating regions with a bounded output, always exact regardless of segment curvature (see the disjunctive section below).

If the curvature doesn’t match the sign (convex + "<=", or concave + ">="), LP is not applicable — method="auto" falls back to SOS2/incremental with the signed link, which gives a valid but much more expensive model. In that case prefer "==" unless you genuinely need the one-sided semantics. See the Piecewise inequalities — per-tuple sign notebook for a full walkthrough.

Warning

With a bounded tuple and active=0, the output is only forced to 0 on the signed side — the complementary bound still comes from the output variable’s own lower/upper bound. In the common case of non-negative outputs (fuel, cost, heat), set lower=0 on that variable: combined with the y ≤ 0 constraint from deactivation, this forces y = 0 automatically. See the docstring for the full recipe.

Breakpoint Construction #

From lists #

The simplest form — pass Python lists directly in the tuple:

m.add_piecewise_formulation(
    (power, [0, 30, 60, 100]),
    (fuel, [0, 36, 84, 170]),
)

With the `breakpoints()` factory #

Equivalent, but explicit about the DataArray construction:

m.add_piecewise_formulation(
    (power, linopy.breakpoints([0, 30, 60, 100])),
    (fuel, linopy.breakpoints([0, 36, 84, 170])),
)

From slopes #

When you know marginal costs (slopes) rather than absolute values, wrap them in linopy.Slopes. The x grid is borrowed from the sibling tuple — no need to repeat it:

m.add_piecewise_formulation(
    (power, [0, 50, 100, 150]),
    (cost, linopy.Slopes([1.1, 1.5, 1.9], y0=0)),
)
# cost breakpoints: [0, 55, 130, 225]

For standalone resolution outside of add_piecewise_formulation, call linopy.Slopes.to_breakpoints() with an explicit x grid:

bp = linopy.Slopes([1.1, 1.5, 1.9], y0=0).to_breakpoints([0, 50, 100, 150])

Per-entity breakpoints #

Different generators can have different curves. Pass a dict to breakpoints() with entity names as keys:

m.add_piecewise_formulation(
    (
        power,
        linopy.breakpoints(
            {"gas": [0, 30, 60, 100], "coal": [0, 50, 100, 150]}, dim="gen"
        ),
    ),
    (
        fuel,
        linopy.breakpoints(
            {"gas": [0, 40, 90, 180], "coal": [0, 55, 130, 225]}, dim="gen"
        ),
    ),
)

Ragged lengths are NaN-padded automatically. Breakpoints are auto-broadcast over remaining dimensions (e.g. time).

Disjunctive segments #

For disconnected operating regions (e.g. forbidden zones), use segments():

m.add_piecewise_formulation(
    (power, linopy.segments([(0, 0), (50, 80)])),
    (cost, linopy.segments([(0, 0), (125, 200)])),
)

The disjunctive formulation is selected automatically when breakpoints have a segment dimension. A bounded tuple ("<=" / ">=") also works here.

N-variable linking #

Link any number of variables through shared breakpoints (joint equality):

m.add_piecewise_formulation(
    (power, [0, 30, 60, 100]),
    (fuel, [0, 40, 85, 160]),
    (heat, [0, 25, 55, 95]),
)

All variables are symmetric here; every feasible point is the same λ-weighted combination of breakpoints across all three. With 3 or more tuples, only "==" signs are accepted — bounding one expression by a multi-input curve isn’t supported yet; see the per-tuple sign section above for the issue link.

Formulation Methods #

Pass method="auto" (the default) and linopy picks the cheapest correct formulation based on sign, curvature and breakpoint layout:

2-variable inequality on a convex/concave curve → lp (chord lines, no auxiliary variables)
All breakpoints monotonic → incremental
Otherwise → sos2
Disjunctive (segments) → always sos2 with binary segment selection

The resolved choice is exposed on the returned PiecewiseFormulation via .method (and .convexity when well-defined). An INFO-level log line explains the resolution whenever method="auto" is in play.

At-a-glance comparison:

Property	`lp`	`incremental`	`sos2`	Disjunctive
Curve layout	Connected	Connected	Connected	Disconnected
Supported per-tuple sign	one `<=` or `>=` (required)	all `==` or one `<=`/`>=`	all `==` or one `<=`/`>=`	all `==` or one `<=`/`>=`
Number of tuples	Exactly 2	≥ 2 (3+ requires all `==`)	≥ 2 (3+ requires all `==`)	≥ 2 (3+ requires all `==`)
Breakpoint order	Strictly monotonic	Strictly monotonic	Any	Any (per segment)
Curvature requirement	Concave (`<=`) or convex (`>=`)	None	None	None
Auxiliary variables	None	Continuous + binary	Continuous + SOS2	Binary + SOS2
`active=` supported	No	Yes	Yes	Yes
Solver requirement	Any LP solver	MIP-capable	SOS2-capable	SOS2 + MIP

LP (chord-line) Formulation #

For 2-variable inequality on a convex or concave curve. Adds one chord inequality per piece plus a domain bound — no auxiliary variables and no MIP relaxation:

\[ \begin{align}\begin{aligned}&y \ \text{sign}\ m_k \cdot x + c_k \quad \forall\ \text{pieces } k\\&x_0 \le x \le x_n\end{aligned}\end{align} \]

where \(m_k = (y_{k+1} - y_k)/(x_{k+1} - x_k)\) and \(c_k = y_k - m_k\, x_k\). For concave \(f\) with sign="<=", the intersection of all chord inequalities equals the hypograph of \(f\) on its domain.

The LP dispatch requires curvature and sign to match: sign="<=" needs concave (or linear); sign=">=" needs convex (or linear). A mismatch is not just a loose bound — it describes the wrong region (see the Piecewise inequalities — per-tuple sign). method="auto" detects this and falls back; method="lp" raises.

# y <= f(x) on a concave f — auto picks LP
m.add_piecewise_formulation((y, yp, "<="), (x, xp))

# Or explicitly:
m.add_piecewise_formulation((y, yp, "<="), (x, xp), method="lp")

Not supported with method="lp": all-equality, more than 2 tuples, and active. method="auto" falls back to SOS2/incremental in all three cases.

The underlying chord expressions are also exposed as a standalone helper, linopy.tangent_lines(x, x_pts, y_pts), which returns the per-piece chord as a LinearExpression with no variables created. Use it directly if you want to compose the chord bound with other constraints by hand, without the domain bound that method="lp" adds automatically.

Incremental (Delta) Formulation #

The default MIP encoding when method="auto" is in play and breakpoints are strictly monotonic — produces a tight MIP directly, without going through an SOS2 layer. Uses fill-fraction variables \(\delta_i\) with binary indicators \(z_i\):

\[ \begin{align}\begin{aligned}&\delta_i \in [0, 1], \quad z_i \in \{0, 1\}\\&\delta_{i+1} \le \delta_i, \quad z_{i+1} \le \delta_i, \quad \delta_i \le z_i\\&e_j = B_{j,0} + \sum_{i=1}^{n} \delta_i \, (B_{j,i} - B_{j,i-1})\end{aligned}\end{align} \]

With a bounded tuple, the link to that tuple’s expression flips to the requested sign while the pinned tuples keep the equality above (see the Per-tuple sign section’s Formulation block).

m.add_piecewise_formulation((power, xp), (fuel, yp), method="incremental")

Limitation: breakpoint sequences must be strictly monotonic.

SOS2 (Convex Combination)#

Fallback when breakpoints aren’t strictly monotonic (the only case method="auto" does not pick incremental for a connected curve). Introduces interpolation weights \(\lambda_i\) with an SOS2 adjacency constraint:

\[ \begin{align}\begin{aligned}&\sum_{i=0}^{n} \lambda_i = 1, \qquad \text{SOS2}(\lambda_0, \ldots, \lambda_n)\\&e_j = \sum_{i=0}^{n} \lambda_i \, B_{j,i} \quad \text{for each expression } j\end{aligned}\end{align} \]

The SOS2 constraint ensures at most two adjacent \(\lambda_i\) are non-zero, so every expression is interpolated within the same piece. With a bounded tuple, the bounded link flips to the requested sign as above.

Note

Solvers with native SOS2 support handle the adjacency constraint via branch-and-bound. Solvers without it see the Big-M reformulation linopy applies (controlled by reformulate_sos= on solve) — see SOS Reformulation for the reformulated MIP form, which is the model those solvers actually receive. When breakpoints are monotonic, prefer method="incremental" (or just "auto"): it builds a similar MIP encoding directly and does not depend on solver SOS2 support or the reformulation step.

m.add_piecewise_formulation((power, xp), (fuel, yp), method="sos2")

Disjunctive (Disaggregated Convex Combination)#

For disconnected segments (gaps between operating regions). Binary indicators \(z_k\) select exactly one segment; SOS2 applies within it:

\[ \begin{align}\begin{aligned}&z_k \in \{0, 1\}, \quad \sum_{k} z_k = 1\\&\sum_{i} \lambda_{k,i} = z_k, \qquad e_j = \sum_{k} \sum_{i} \lambda_{k,i} \, B_{j,k,i}\end{aligned}\end{align} \]

No big-M constants are needed, giving a tight LP relaxation.

Disjunctive + bounded tuple. A per-tuple "<=" / ">=" works here too, applied to the bounded tuple exactly as for the continuous methods. Because the disjunctive machinery already carries a per-segment binary, there is no curvature requirement on the segments — inequality is always exact on the hypograph (or epigraph) of the active segment, whatever its slope pattern. This makes disjunctive plus a bounded tuple a first-class tool for “bounded output on disconnected operating regions” that method="lp" cannot handle.

Advanced Features #

Active parameter (unit commitment)#

The active parameter gates the piecewise function with a binary variable. When active=0, all auxiliary variables (and thus the linked expressions) are forced to zero:

commit = m.add_variables(name="commit", binary=True, coords=[time])
m.add_piecewise_formulation(
    (power, [30, 60, 100]),
    (fuel, [40, 90, 170]),
    active=commit,
)

commit=1: power operates in [30, 100], fuel = f(power)
commit=0: power = 0, fuel = 0

Not supported with method="lp".

Note

With a bounded tuple, deactivation only pushes the signed bound to 0 — the complementary side comes from the output variable’s own lower/upper bound. Set lower=0 on naturally non-negative outputs (fuel, cost, heat) to pin the output to zero on deactivation. See the per-tuple sign section above for details.

Auto-broadcasting #

Breakpoints are automatically broadcast to match expression dimensions — you don’t need expand_dims:

time = pd.Index([1, 2, 3], name="time")
x = m.add_variables(name="x", lower=0, upper=100, coords=[time])
y = m.add_variables(name="y", coords=[time])

# 1D breakpoints auto-expand to match x's time dimension
m.add_piecewise_formulation((x, [0, 50, 100]), (y, [0, 70, 150]))

NaN masking #

Trailing NaN values in breakpoints mask the corresponding lambda / delta variables (and, for LP, the corresponding chord constraints). This is useful for per-entity breakpoints with ragged lengths:

# gen1 has 3 breakpoints, gen2 has 2 (NaN-padded)
bp = linopy.breakpoints({"gen1": [0, 50, 100], "gen2": [0, 80]}, dim="gen")

Interior NaN values (gaps in the middle) are not supported and raise an error.

Inspecting generated objects #

The returned PiecewiseFormulation exposes .variables and .constraints as live views into the model — use them to introspect exactly what was generated, rather than relying on documented name conventions:

f = m.add_piecewise_formulation((y, y_pts, "<="), (x, x_pts))
print(f)  # method, convexity, vars/cons summary

The comparison table above describes the kind of auxiliary objects each method creates (continuous + SOS2, binary + SOS2, none, …); exact name suffixes are an implementation detail and may evolve.

Piecewise Linear Constraints

Contents