feat: Unify piecewise API behind add_piecewise_formulation (sign + LP dispatch)

[FBumann]

Summary

Follows up on #602. Replaces the descriptor pattern (PiecewiseExpression, PiecewiseConstraintDescriptor, piecewise()) with a single, stateless construction layer: add_piecewise_formulation.

add_piecewise_formulation becomes the one entry point for piecewise work — equality, one-sided inequality, N-variable linking, per-entity breakpoints, unit-commitment gating, and automatic method dispatch all live behind one call. tangent_lines() survives as a low-level helper for users who want to build chord expressions manually.

The API

Each (expression, breakpoints[, sign]) tuple links a variable to its breakpoints. All tuples share interpolation weights, coupling them on the same curve piece. The optional 3rd element marks one tuple as bounded by the curve instead of pinned to it.

Equality (default)

# 2-variable: fuel = f(power)
m.add_piecewise_formulation(
    (power, [0, 30, 60, 100]),
    (fuel,  [0, 36, 84, 170]),
)

# N-variable (CHP plant — same pattern, more tuples)
m.add_piecewise_formulation(
    (power, [0, 30, 60, 100]),
    (fuel,  [0, 40, 85, 160]),
    (heat,  [0, 25, 55, 95]),
)

# Disjunctive (disconnected operating regions)
m.add_piecewise_formulation(
    (power, linopy.segments([(0, 0), (50, 80)])),
    (cost,  linopy.segments([(0, 0), (125, 200)])),
)

# Per-entity breakpoints (different curves per generator)
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")),
)

# Breakpoints from slopes
m.add_piecewise_formulation(
    (power, [0, 30, 60, 100]),
    (fuel,  linopy.breakpoints(slopes=[1.2, 1.4, 1.7], x_points=[0, 30, 60, 100], y0=0)),
)

Inequality bounds — per-tuple sign

Append "<=" or ">=" as a third tuple element to mark one expression as bounded by the curve. The other tuples stay pinned. The bounded tuple need not be first — its role is visible at the call site.

# fuel ≤ f(power).  "auto" picks the cheapest correct formulation — a pure
# LP chord formulation when the curvature matches the sign (concave + "<=",
# or convex + ">="), SOS2/incremental otherwise.
m.add_piecewise_formulation(
    (fuel,  [0, 20, 30, 35], "<="),  # bounded
    (power, [0, 10, 20, 30]),         # pinned
)

On a matching-curvature curve this dispatches to method="lp": one chord inequality per piece plus a domain bound x_min ≤ x ≤ x_max, and no auxiliary variables. Mismatched curvature (convex + "<=" / concave + ">=") describes the wrong region (not just a looser one) — method="auto" detects this and falls back to SOS2/incremental with an explanatory info log; method="lp" raises with a clear error.

Restrictions (current):

• At most one tuple may carry a non-equality sign.
• With 3+ tuples, all signs must be "==".

The N-input bounded case (e.g. z >= f(x, y)) is rejected with an error that invites a GitHub issue if anyone has a concrete use case — we don't promise a specific future API shape.

Unit-commitment gating — active

# active=1: power in [30, 100], fuel = f(power)
# active=0: power = 0, fuel = 0
m.add_piecewise_formulation(
    (power, [30, 60, 100]),
    (fuel,  [40, 90, 170]),
    active=commit,  # binary variable
)

Works with SOS2, incremental, and disjunctive. With a bounded tuple, deactivation only pushes the signed bound to 0 — set lower=0 on naturally non-negative outputs (fuel, cost, heat) and the combination forces y = 0 automatically.

Low-level: tangent_lines

# Returns a LinearExpression with one chord per piece — no variables.
# Most users should prefer add_piecewise_formulation with a bounded tuple,
# which builds on this helper and adds domain bounds + curvature checks.
t = linopy.tangent_lines(power, x_pts, y_pts)
m.add_constraints(fuel <= t)

Key changes

New public API

• add_piecewise_formulation — one entry point for all piecewise constructions. Returns a PiecewiseFormulation carrying .method (resolved formulation), .convexity (when well-defined), .variables, .constraints.
• Per-tuple sign — optional 3rd tuple element, default "==". At most one non-equality; 3+ tuples must all be equality.
• method parameter — "auto" (default), "sos2", "incremental", "lp".
• active parameter — binary gating for unit commitment.
• linopy.breakpoints() / linopy.segments() factories.
• linopy.slopes_to_points() utility.
• linopy.tangent_lines() low-level helper.

Auto-dispatch (method="auto")

1. 2-variable inequality on a matching-curvature curve → lp (chord + domain bounds, no aux vars)
2. Strictly monotonic breakpoints → incremental (deltas + binaries)
3. Otherwise → sos2 (lambdas + SOS2 constraint)
4. Disjunctive (segments) → sos2 with binary segment selection

The resolved method is logged at INFO; when LP is skipped, the log explains why (curvature, NaN layout, tuple count, or active).

Correctness properties

• _detect_convexity is direction-invariant — ascending and descending x give the same label for the same graph.
• NaN-padded per-entity breakpoints are masked everywhere, including the LP chord constraints (padded pieces don't create spurious y ≤ 0 constraints).
• method="lp" validates curvature + sign + tuple count + active absence upfront; mismatch raises with a clear pointer to method="auto".
• PiecewiseFormulation.method and .convexity persist across netCDF round-trip.

Architecture (internal)

• _PwlInputs dataclass categorizes the user's tuples into pinned_* and bounded_* slots — explicit role assignment, no positional convention. Carries through the whole dispatch chain.
• _PwlLinks is the stacked link representation consumed by SOS2/incremental/disjunctive builders; _build_links() produces it from _PwlInputs.
• _try_lp() / _resolve_sos2_vs_incremental() flatten the dispatch — _add_continuous is ~10 lines of real logic.
• User tuple order is preserved end-to-end; the bounded one is found via _PwlInputs, not by assuming position 0.

Removed / deprecated

• PiecewiseExpression, PiecewiseConstraintDescriptor, linopy.piecewise() — descriptor pattern gone.
• Old add_piecewise_constraints shape that dispatched on Variable.__le__/__ge__/__eq__ return types.
• Variable operator overloads (__le__, __ge__, __eq__) simplified — no more piecewise dispatch.

Design principles

• API surface stays minimal: Model, Variables, Constraints, Expression, plus the single PiecewiseFormulation return type. No user-facing intermediate state.
• Piecewise is a construction layer that produces regular linopy objects.
• One entry point for equality and inequality — per-tuple sign reads as the relation it encodes ((y, y_pts, "<=") ⇒ y ≤ f(...)).
• Variable names as link coords: the internal _pwl_var dimension shows power, fuel, heat instead of 0, 1, 2.

Alternative designs considered

A) Descriptor pattern (PR #602, previous master)

pw = linopy.piecewise(x, x_pts, y_pts)
m.add_piecewise_constraints(pw == y)

Rejected: Introduces PiecewiseExpression and PiecewiseConstraintDescriptor as user-facing state. Breaks the "only Variables, Constraints, Expressions" principle. Structurally limited to 2-variable x → y.

B) Dict of expressions + shared breakpoints DataArray

m.add_piecewise_formulation(exprs={"power": power, "fuel": fuel}, breakpoints=bp)

Considered: Supports N-variable. But requires string keys to match breakpoint coordinates — an indirection layer.

C) Global sign= keyword (initial implementation, refactored in #664)

m.add_piecewise_formulation((fuel, y_pts), (power, x_pts), sign="<=")

The original implementation used a single formulation-level sign= kwarg with a positional convention ("the first tuple is the bounded one when sign != '=='"). #664 refactored this to per-tuple sign because the role belongs to a specific tuple — encoding it globally + positionally was a wart that required a paragraph of docs to teach. The internal dispatch was rewritten from positional access to explicit categorization (_PwlInputs with bounded_* and pinned_* slots).

D) Chosen: per-tuple sign on the tuple it applies to

m.add_piecewise_formulation(
    (fuel,  y_pts, "<="),  # bounded (role visible at the call site)
    (power, x_pts),         # pinned
)

Future work this enables

• 2-D triangulated piecewise (z ≥ f(x, y)). The reserved N≥3-inequality slot reads naturally as the epigraph of a bivariate function: (z, z_pts, ">="), (x, x_pts), (y, y_pts), triangulation=tris. Vertex coordinates stay as plain 1-D DataArrays; the simplex adjacency rides along as one shared sibling kwarg. _PwlInputs would grow a single optional field; the formulation builders for the triangulated path are new but consume the same input carrier.
• Multi-bounded relations (if a use case appears) — the "at most one bounded tuple" check is a single line; the dispatch chain already handles per-tuple roles.

Defered Features

• More polished slopes api, potentially with a new class

Out of scope / follow-ups

Directions the design enables or considered but deliberately not shipping in this PR:

• 2-D triangulated piecewise (z ≥ f(x, y)). The reserved N≥3-inequality slot reads naturally as the epigraph of a bivariate function: (z, z_pts, ">="), (x, x_pts), (y, y_pts), triangulation=tris. Vertex coordinates stay as plain 1-D DataArrays; the simplex adjacency rides along as one shared sibling kwarg. _PwlInputs would grow a single optional field; the formulation builders for the triangulated path are new but consume the same input carrier.
• Multi-bounded relations (if a use case appears) — the "at most one bounded tuple" check is a single line; the dispatch chain already handles per-tuple roles.
• More polished slopes API, potentially with a new class. breakpoints(slopes=..., x_points=..., y0=...) is shipped today. What's deferred is letting slopes inherit the x grid from a sibling tuple (e.g. (fuel, Slopes([1.2, 1.4, 1.7], y0=0)) paired with a reference tuple that carries the x grid). Two viable shapes were discussed — extending breakpoints(slopes=...) to allow omitted x_points, or introducing a small frozen-dataclass Slopes value type — each with tradeoffs around type-system visibility and dispatch clarity. EvolvingAPIWarning covers either direction as a non-breaking follow-up.

Notebook examples

examples/piecewise-linear-constraints.ipynb (main tutorial):

#	Section	Feature
1	Getting started	Basic 2-variable equality + method inspection via pwf
2	Picking a method	auto vs sos2 vs incremental give the same optimum
3	Disjunctive segments	segments() for gaps in operation
4	Inequality bounds	Per-tuple "<=" (auto-dispatches LP)
5	Unit commitment	active parameter
6	N-variable linking	CHP plant (3 variables, joint equality)
7	Per-entity breakpoints	Fleet with different curves

examples/piecewise-inequality-bounds.ipynb (deep-dive on per-tuple sign / LP):

• Setup with a concave reference curve
• Three methods (LP, SOS2, incremental) verified to produce the same feasible region
• 2-D hypograph feasibility plot
• Curvature × sign rule table — the "wrong region" cases and why
• Auto-dispatch fallback demonstration (non-convex curve, mismatched sign, explicit method="lp" raise)

Formulation math lives in the reference page, not in the notebook.

Test plan

• pytest test/test_piecewise_constraints.py test/test_piecewise_feasibility.py — 518 passed
• Full suite (excl. solver/remote) — 1578 passed, 7 skipped
• Both notebooks execute end-to-end against Gurobi
• Ruff lint + format clean, mypy clean on the full repo
• Round-trip tested: PiecewiseFormulation.method and .convexity persist across to_netcdf / read_netcdf
• Review

🤖 Generated with Claude Code

[FBumann]

@FabianHofmann This is pretty urgent to me, as the current piecewise API is already on master and should NOT be included in the next release in my opinion.
But i understand that the CSR PR is more important atm.

The current code can be reorganized a bit, but id like to hear your thoughts about the API first.

[FBumann]

@coroa I thought a lot about the inequality in non lp piecewise relationships.
I added a notebook and plot about it in #662 (5992788). Essentially, the sign needs to be "translated" to behave like tangent lines. And we need to make a decision about the convention, because givin all options to the user makes a pretty bad UX whie "allowing" cases which probably noone wants.
Here is a visualization for feasible regions for a 2-variable case with SOS2/lambda/delta formulations based on which sign is used for which variable:

Compared to tangent line one:

We can see that we need to fix one of the variables to the slope! TO be precise, all but one variable needs to be fixed!

[FBumann]

@coroa What do you think about the new api?
I think its quite good and clear now. The testing is also already quite good imo.

One thing we should discuss is the notebooks and docs. If you have the time to read it, feedback and potentially changes would be great!
The one who implemented it usually has a different angle on the docs than the reviewer/user.
https://linopy--638.org.readthedocs.build/en/638/piecewise-linear-constraints.html

[FBumann]

@coroa One API detail im not sure about yet, is where to but the sign.
It could be per tuple as a third optional parameter:

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

with defaulting to '==' where not specified.

I decided against this api for simplicity, as using more than i single equality produces mathematically weird bounds as far as i understand it (only relevant for the 3 var case). ALso, for the 2 variable LP case, it doestn make any sense to have more than a single sign. So there is a lot of API surface that would lead to errors we need to catch.

We could however also just disallow the 3 variable case with inequalities for now, as i don`t understand the use case that well. This would allow the arguably more obvious notation of having the sing right next to the expression...

With the EvolvingAPIWarning however, im fine with keeping it and adjusting it if needed.

Looking forward, this API seems better suited for triangular bounding, where we could actually allow multiple inequalities. But the implementation of this is not part of this PR and quite complex i think, so we can just change it as soon as this potential feature lands

[FBumann]

One thing I didn't cover yet, is a dedicated slopes object. Should we add it? Or can we do it in a follow up?

[FBumann]

@coroa One API detail im not sure about yet, is where to but the sign. It could be per tuple as a third optional parameter:

I thought a bit more about it, and its much clearer having the sign next to the variable/factor. THis also extends naturally to triangulation etc, and feels less "magic". I refactored it in #664. This also made the docs cleaner.

I also decided to disallow inequalities for 3+ variables, as im not sure about the mathematical formulation and would rather ship the PR as is without further discussion on the not yet needed 3+ != case.

[FBumann]

Refactor: Disambiguate usage of "segment" by renaming to "piece" in most contexts

Both terms are used in the PWL literature — Wikipedia and the Northwestern optimization wiki use them interchangeably; JuMP's PiecewiseLinearOpt.jl prefers "pieces", Pyomo's API leans "segments". So the rename isn't about correctness - its about readability.
I argue strongly in favor of Pieces, while keeping Segments as a collection of curves, where each curve/segment consist of several pieces (disjunctive). This is clear and mathematically sound.

I added a clear definition into the docs

[FBumann]

@coroa @FabianHofmann I think everything is sorted out now and the api is both intuitive, well documented, tested and future proof.

Please take a look