Piecewise Linear Constraints

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.

Throughout this page: a breakpoint is a knot where the slope can change; a piece is the linear part between adjacent breakpoints; a segment is a disjoint operating region in the disjunctive formulation. Per-tuple breakpoint arrays are paired by index — the i-th entries across all tuples together define one knot.

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) on the same heat-rate curve as above.
m.add_piecewise_formulation(
    (fuel, [0, 36, 84, 170], ">="),
    (power, [0, 30, 60, 100]),
)

Over-fuelling is physically admissible but wasteful, so minimising fuel pulls the operating point onto the curve. method="auto" picks the cheapest correct formulation: pure LP (chord constraints) here, since convex + ">=" is LP-applicable; SOS2/incremental otherwise.

Each (expression, breakpoints[, sign]) tuple pairs a variable with its breakpoint values. The optional sign (default "==") is "<=" or ">=" to mark that expression as bounded by the curve. With every sign "==", all tuples land on the same point of the piecewise curve — see Per-tuple sign below for the geometry of the inequality cases.

API#

add_piecewise_formulation#

m.add_piecewise_formulation(
    (expr1, breakpoints1),  # sign defaults to "==" (equality role)
    (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
)

Adds constraints — and, depending on the resolved method, auxiliary variables — for either an all-equality joint (every tuple at the same point on the curve, the default) or a one-sided bound on a single tuple. The pure-LP path adds chord and domain constraints only; SOS2, incremental, and disjunctive also add interpolation weights and/or binaries (see Formulation Methods below).

Breakpoint Construction#

Each tuple’s breakpoints come from breakpoints() (a single connected curve) or segments() (disjoint operating bands). Slopes can stand in for breakpoints() when per-piece slopes are the natural input — it resolves to a breakpoints array.

breakpoints() — a connected curve#

Values along a single connected piecewise curve — the default case for efficiency curves, heat rates, and cost curves.

The simplest form passes a Python list directly in the tuple:

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

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])),
)

Per-entity curves. 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).

Specifying by slopes. linopy.Slopes resolves to a breakpoint array from per-piece slopes plus an initial y0, instead of from absolute y-values — useful when slopes are the natural input (e.g. marginal costs). The x grid is borrowed from the sibling tuple, so the y breakpoints don’t have to be computed by hand:

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

For standalone resolution outside 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])

segments() — disjoint operating bands#

For equipment with disconnected operating bands. Each segment is one band’s (range, curve); a binary picks exactly one per operating point, with continuous interpolation within the chosen band.

# Stepped pump with two speed bands.
m.add_piecewise_formulation(
    (flow, linopy.segments([(5, 25), (40, 100)])),
    (power, linopy.segments([(1, 7), (15, 50)])),
)

Bounded tuples ("<=" / ">=") are supported on disjunctive curves too.

For a single on/off gate on one continuous curve, prefer active=... (see Active parameter (unit commitment)) — using a degenerate (0, 0) segment to encode “off” mixes the disjunctive concept with on/off logic.

N-variable linking#

Independent of the building block used, any number of variables can be linked 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. Sign restrictions apply (see Per-tuple sign below) — for N ≥ 3 tuples all signs must be "==".

Per-tuple sign — equality vs inequality#

Roles and restrictions#

Each tuple’s optional third element is a sign:

  • "==" (default) — equality role: the tuple enters as an equality.

  • "<=" / ">=" — bounded: the expression undershoots / overshoots the curve.

Note

Current restrictions.

  • 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.

Geometry#

What the formulation actually constrains depends on the tuple count and signs:

  • All-equality (default). Shared interpolation weights put the joint \((e_1, \ldots, e_N)\) exactly on the curve.

  • One bounded + one equality (2 tuples). The joint \((x, y)\) lies in the hypograph ("<=" on a concave \(f\)) or epigraph (">=" on a convex \(f\)):

    \[\{ (x, y) \ :\ x_{\min} \le x \le x_{\max},\ y \le f(x) \} \qquad \text{(hypograph)}\]

    The equality axis is just confined to its breakpoint domain \([x_{\min}, x_{\max}]\) — a single coordinate can’t locate a curve point. Same projection in LP (enforced directly) and SOS2/incremental (enforced via the weight link).

  • Mismatched sign + curvature (convex + "<=", or concave + ">=") describes a non-convex region — method="auto" falls back to SOS2/incremental and method="lp" raises.

# All-equality: joint (x, y) on the curve.
m.add_piecewise_formulation((y, y_pts), (x, x_pts))

# Bounded: joint (x, y) in the hypograph — y ≤ f(x), x ∈ [x_min, x_max].
m.add_piecewise_formulation((y, y_pts, "<="), (x, x_pts))

# Bounded: joint (x, y) in the epigraph — y ≥ f(x), x ∈ [x_min, x_max].
m.add_piecewise_formulation((y, y_pts, ">="), (x, x_pts))

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

Choice of bounded tuple and sign#

Pick the sign matching the physically admissible direction for that expression:

  • "<=" for a quantity with a controllable dissipation path — heat rejection via cooling tower (thermal curtailment), electrical curtailment, emissions after post-treatment — so undershooting the curve is realisable.

  • ">=" for an input whose over-supply is admissible but wasteful — fuel, raw materials — so overshooting the curve is realisable (objective pressure then pulls the operating point onto the curve).

The wrong direction ("<=" on fuel, ">=" on a non-curtailable output) yields a valid but loose formulation that admits operating points the plant cannot physically realise; an objective rewarding the wrong direction 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 Creating Piecewise Inequality Bounds notebook for a full walkthrough.

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) → SOS2 applied per segment with binary segment selection (the disjunctive formulation in the table below).

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

Continuous + binary + SOS2

active= supported

No

Yes

Yes

Yes

Solver requirement

Any LP solver

MIP-capable

SOS2-capable (or MIP via Big-M reformulation)

SOS2 + MIP (or MIP via Big-M reformulation)

Note

Disjunctive formulations report method="sos2" (the underlying per-segment encoding uses SOS2), but the table treats them as a separate column because the per-segment binaries change the auxiliary-variable structure and solver requirements.

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_{\min} \le x \le x_{\max}\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\). The domain bound uses \(x_{\min}\) and \(x_{\max}\) rather than the first/last breakpoint so that descending x grids work too — strictly-monotonic breakpoints are accepted in either order. 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 Creating Piecewise Inequality Bounds). 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.

Chord expressions as a building block#

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

chord = linopy.tangent_lines(x, x_pts, y_pts)
m.add_constraints(y <= chord + slack)

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 equality-signed tuples keep the equality above.

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" (gating needs a binary). Use method="auto", or Chord expressions as a building block for manual gating.

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.

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.

Tutorials#

See Also#