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.

[FabianHofmann]

this looks cool, well done @FBumann ! I will have a look in the afternoon!

[FabianHofmann]

first round of small requests mostly outside of the piecewise module.

[FabianHofmann]

piecewise variables and constraints cannot have the same name right?

[FBumann]

Yes, they can. As regular variables and constraints can share a name too. They do not currently though.

[FabianHofmann]

@FBumann just went through the notebook, fantastic work. mind if I push some decorative changes to the notebook?

[FBumann]

Of course not. Go ahead. You dont have to ask that going forward. If you keep a proper commit history, i can check what you did and complain if i need to :)

[FabianHofmann]

Here is the review from my agents:

Critical

• piecewise.py — _lp_eligibility (1172-1200) and _try_lp (1278-1308) duplicate every check. Two sources of truth for LP eligibility rules. Collapse into one helper returning (ok, reason) and raise from _try_lp when method == "lp".
• __init__.py:60-66 — __all__ no longer alphabetical. tangent_lines sits between breakpoints and segments; align/merge/options/read_netcdf come after slopes_to_points. The latter group was already broken on master, but the new tangent_lines insertion compounds it. Likely to trip ruff RUF022.
• No tests for EvolvingAPIWarning at all. grep -rn "EvolvingAPIWarning" test/ returns zero matches. Per-key dedup, stacklevel=3, and the "no double-warn from _tangent_lines_impl inside _add_lp" claim are all unverified. Module-global dedup set will make tests order-dependent without a reset fixture.
• NetCDF roundtrip is a single 2-tuple equality case (test_formulation_survives_netcdf, line 2048). Missing: multi-formulation models, disjunctive formulations, bounded <=/>= tuples, method="lp" (empty variable_names), and convexity != None.

High

Code quality / DRY

• piecewise.py:648-659 _validate_breakpoint_shapes repeats the BREAKPOINT_DIM error twice for bp_a/bp_b. Loop the pair.
• piecewise.py:1227-1267 _build_links has two near-identical equality-branch builders.
• piecewise.py:1038-1059 _PwlInputs two-branch construction differs in 4 fields only.
• piecewise.py:127-140 PiecewiseFormulation hand-rolls __slots__+__init__; use @dataclass(slots=True) to drop ~25 lines of boilerplate.
• piecewise.py:154-188 _user_dims and __repr__ duplicate the dim-collection / _-prefix-skip loop.

Type / API

• piecewise.py:131 method: str should be PWL_METHOD (or narrowed Literal); same for convexity.
• piecewise.py:797-805 add_piecewise_formulation(**kwargs: object) defeats static type-checking on callers and conflicts with the no-getattr-style rule. Use explicit sign: str | None = None and raise.
• io.py:1265-1269 Inconsistent fail-fast posture — convexity uses d.get("convexity") (silent None) while every other field uses d["..."] (KeyError). Pick one.

Tests

• PiecewiseFormulation.variables / .constraints / __repr__ / _user_dims essentially unexercised. Only f.method/f.convexity are read in tests.
• Auto-method success log never asserted (only the LP-skip caplog branch at line 1850 exists).
• _lp_eligibility reason strings not individually asserted. Each branch returns a distinct user-facing diagnostic; parametrize directly.
• _add_continuous resolution for method ∈ {"sos2","incremental","lp","auto"} not parametrized (only ["sos2","incremental"] at line 1732).

Medium

Code quality

• piecewise.py:62-74 Custom _emitted_evolving_warnings set is over-engineered; warnings.simplefilter("once", EvolvingAPIWarning) would suffice.
• piecewise.py:1364 _resolve_sos2_vs_incremental final raise ValueError("unknown method ...") is dead — add_piecewise_formulation validates against PWL_METHODS upstream.
• piecewise.py:1156-1169 _PwlInputs.all_bps/all_coords/all_exprs three near-identical methods; collapse.
• piecewise.py:1033 link-coord fallback to str(i) can collide with a user variable named "1". Use f"_pwl_{i}".
• model.py:506 indirect _grouped_names import → consider exposing as model._piecewise_grouped_names().
• io.py:1259-1269 move reconstruction into PiecewiseFormulation.from_dict for symmetry with to_netcdf.

Type / API

• piecewise.py:1146 link_dim = "_pwl_var" is ad-hoc; promote to PWL_LINK_DIM constant alongside BREAKPOINT_DIM/SEGMENT_DIM.
• piecewise.py:1448 piece_dim = f"{dim}_piece" recomputed in _add_continuous instead of reusing LP_PIECE_DIM (imported at line 28, used in _rename_to_pieces at 236-237). Two sources of truth.
• piecewise.py:191 _grouped_names is opaque naming → _piecewise_member_names.
• piecewise.py tangent_lines (~600-635) and add_piecewise_formulation (~900-925) docstrings are long; the EvolvingAPIWarning silencing recipe is duplicated verbatim across both sites.
• variables.py:545 bare # type: ignore should be # type: ignore[override] (regression — was [override] before).
• piecewise.py:1116-1125 _stack_along_link accepts Sequence[DataArray | xr.Dataset] but returns DataArray; either generic via TypeVar or split.

Tests

• add_piecewise_formulation unexpected-kwargs TypeError path untested.
• _breakpoints_from_slopes 3D slopes path, dict y0 missing key, per-entity length mismatch — direct tests missing.
• _validate_breakpoint_shapes first-arg branch (mirror of the bp_b test at 1360) untested.
• Public breakpoints / segments factories: empty list, single point, all-NaN, ±inf, and duplicate x (zero dx → div-by-zero in tangent_lines / _detect_convexity) all uncovered.
• LP _add_lp piece_count == 1 shortcut (piece_mask all-True) not tested.
• _broadcast_points "Cannot broadcast: missing dim" branch has no negative test.
• SOS2/incremental masking branch of bounded path with NaN not directly verified.

Low

• piecewise.py:229-231 _strip_nan returns list[float]; use arr[~np.isnan(arr)].tolist() (faster, type-consistent).
• piecewise.py:1187, 1287 two assert ... is not None # narrowed by is_equality check comments — drop the trailing rationale comments.
• piecewise.py:550-551 _tangent_lines_impl calls _coerce_breaks before the type check; reverse for fail-fast.
• piecewise.py:216-217 _repr_summary uses len(pwl.variables) / len(pwl.constraints) (constructs views) — len(pwl.variable_names) is O(1).
• piecewise.py:1239 _PwlLinks.eq_bp = stacked_bp is redundant in the equality case.
• constants.py:58,60 PWL_METHODS / PWL_CONVEXITIES could be frozenset for immutability.
• constants.py:84-107 EvolvingAPIWarning ~24-line docstring is verbose; two silencing examples could be one.
• model.py:808 add_piecewise_formulation = add_piecewise_formulation mypy-invisible binding survives refactor (pre-existing).
• test_piecewise_feasibility.py:274-356 TestHandComputedAnchors six near-identical methods — collapse via parametrization.
• test_piecewise_constraints.py:1823 manual for sign in signs loop instead of pytest.parametrize.
• test_piecewise_constraints.py:1933-2034 from linopy.piecewise import _detect_convexity repeated in every method — hoist to module scope.
• test_formulation_survives_netcdf:2048 excludes _model from the slot comparison without affirming the back-reference is rebound — that back-reference is what makes .variables / .constraints work after read.
• TestTangentLines:370-421 covers var/linexpr/dataarray/constraint paths but no BreaksLike size-mismatch or NaN tests.

Positive observations

• *pairs variadic API at piecewise.py:797 is a clean upgrade over the old descriptor >=-overload trickery.
• Per-tuple (expr, breakpoints, "<=") sense convention reads naturally and is well-checked.
• EvolvingAPIWarning design (subclass FutureWarning at constants.py:84, per-key dedup at piecewise.py:66-74 with stacklevel=3, message-prefix targeting) is well thought out.
• PiecewiseFormulation return value with .variables / .constraints views is a real improvement over returning a bare Constraint.
• API asymmetry to consider: add_piecewise_formulation is a Model method (model.py:808) but tangent_lines is module-level (piecewise.py:573). Either expose model.tangent_lines(...) for discoverability or document the split.

---

Recommended action plan

• Collapse _lp_eligibility / _try_lp duplication.
• Restore alphabetical __all__ sort (consider also fixing the pre-existing ordering).
• Add EvolvingAPIWarning tests with a fixture clearing _emitted_evolving_warnings.
• Expand netCDF roundtrip tests (multi-formulation, bounded, LP, convexity).
• Tighten types: method: PWL_METHOD; replace **kwargs: object with explicit sign: param.
• Make io.py:1269 convexity loading consistent with other field reads (drop .get).
• Cover PiecewiseFormulation.variables / .constraints / __repr__, parametrize _lp_eligibility reasons, assert the auto-resolution log.
• DRY refactors in _validate_breakpoint_shapes (648-659), _build_links (1227-1267), _PwlInputs (1038-1059).

Defer to follow-up cleanup PR (Medium + Low): dataclass conversion, constants extraction (PWL_LINK_DIM, reuse LP_PIECE_DIM), docstring trimming, test parametrization cleanups, _strip_nan speedup, comment removal.

I would say nothing of this is really critical, but I think the recommended points should be easy to address. perhaps @FBumann you want to go through this and see what of it you find helpful. I would say nothing of this is a requirement.

[FBumann]

@FabianHofmann Seems reasonable.

A few comments on things i disagree on:

1. kwargs: Good catch, but those should just go entirely. They are mere an artifact from iterations during dev...
2. _build_links and _PwlInputs: I dont see how we could compress this meaningfully.

[FBumann]

Great quality improvements @FabianHofmann 🎉