diff --git a/.github/workflows/benchmark-smoke.yml b/.github/workflows/benchmark-smoke.yml index 96e75ae1..1765fe58 100644 --- a/.github/workflows/benchmark-smoke.yml +++ b/.github/workflows/benchmark-smoke.yml @@ -1,7 +1,7 @@ name: Benchmark smoke -# Builds every spec and fires every phase once (--benchmark-disable): -# a "did a refactor break a spec?" check, not timing. +# Builds every spec and fires every phase once (--benchmark-disable) under both +# arithmetic conventions: a "did a refactor break a spec?" check, not timing. on: push: @@ -15,8 +15,14 @@ concurrency: jobs: smoke: - name: Benchmark smoke (quick) + name: Benchmark smoke (${{ matrix.semantics }}) runs-on: ubuntu-latest + strategy: + fail-fast: false + matrix: + semantics: [legacy, v1] + env: + LINOPY_BENCH_SEMANTICS: ${{ matrix.semantics }} steps: - uses: actions/checkout@v7 diff --git a/.github/workflows/semantics-report.yml b/.github/workflows/semantics-report.yml new file mode 100644 index 00000000..a71232de --- /dev/null +++ b/.github/workflows/semantics-report.yml @@ -0,0 +1,121 @@ +name: Semantics report (v1 vs legacy) + +on: + push: + branches: [ master ] + pull_request: + branches: [ '*' ] + +concurrency: + group: ${{ github.workflow }}-${{ github.ref }} + cancel-in-progress: true + +jobs: + report: + name: Build cost β€” v1 vs legacy + runs-on: ubuntu-latest + permissions: + contents: write + pull-requests: write + steps: + - uses: actions/checkout@v7 + with: + fetch-depth: 0 + - name: Set up Python 3.13 + uses: actions/setup-python@v6 + with: + python-version: "3.13" + - name: Install package and benchmark dependencies + run: | + python -m pip install uv + uv pip install --system -e ".[dev,benchmarks]" + - name: Build under each semantics (memray peak + time) + run: | + pytest benchmarks/drivers/test_build.py \ + --benchmark-only --benchmark-memory --benchmark-json=legacy.json + LINOPY_BENCH_SEMANTICS=v1 pytest benchmarks/drivers/test_build.py \ + --benchmark-only --benchmark-memory --benchmark-json=v1.json + - name: Render plots + table + run: | + mkdir -p semantics-report + benchmem plot legacy.json v1.json --columns peak --view scatter -o semantics-report/peak-relative.html + benchmem plot legacy.json v1.json --columns peak --view compare -o semantics-report/peak-absolute.html + benchmem plot legacy.json v1.json --columns time --view scatter -o semantics-report/time-relative.html + benchmem plot legacy.json v1.json --columns peak --view scatter -o semantics-report/peak-relative.svg + benchmem plot legacy.json v1.json --columns time --view scatter -o semantics-report/time-relative.svg + benchmem compare legacy.json v1.json --columns time,peak --stat mean --sort change \ + --csv semantics-report/comparison.csv \ + | sed 's/\x1b\[[0-9;]*m//g' > semantics-report/comparison.txt + run_url="$GITHUB_SERVER_URL/$GITHUB_REPOSITORY/actions/runs/$GITHUB_RUN_ID" + { + echo '## Build cost β€” v1 vs legacy (peak + time, mean)' + echo '' + echo '_v1 build cost relative to legacy on this commit β€” not a comparison against master (that is CodSpeed)._' + echo '' + echo '```' + cat semantics-report/comparison.txt + echo '```' + echo '' + echo "πŸ“Š Plots + CSV: the **semantics-report-v1-vs-legacy** artifact on [this run]($run_url)." + } >> "$GITHUB_STEP_SUMMARY" + - name: Upload report + uses: actions/upload-artifact@v4 + with: + name: semantics-report-v1-vs-legacy + path: semantics-report/ + - name: Publish plot image + compose comment + if: >- + github.event_name == 'pull_request' && + github.event.pull_request.head.repo.full_name == github.repository + env: + GH_TOKEN: ${{ github.token }} + PR: ${{ github.event.pull_request.number }} + HEAD_SHA: ${{ github.event.pull_request.head.sha }} + ASSETS_BRANCH: ci-assets-semantics-report-v1-vs-legacy-plot-images + run: | + remote="https://x-access-token:${GH_TOKEN}@github.com/${GITHUB_REPOSITORY}.git" + tmp="$(mktemp -d)" + git clone --depth=1 --branch "$ASSETS_BRANCH" "$remote" "$tmp" 2>/dev/null || { + git clone --depth=1 "$remote" "$tmp" + git -C "$tmp" checkout --orphan "$ASSETS_BRANCH" + git -C "$tmp" rm -rf . >/dev/null 2>&1 || true + } + mkdir -p "$tmp/pr-$PR" + cp semantics-report/peak-relative.svg semantics-report/time-relative.svg "$tmp/pr-$PR/" + git -C "$tmp" config user.name "github-actions[bot]" + git -C "$tmp" config user.email "41898282+github-actions[bot]@users.noreply.github.com" + git -C "$tmp" add "pr-$PR" + git -C "$tmp" commit -m "semantics-report plot for pr-$PR" || echo "no image change" + git -C "$tmp" push origin "$ASSETS_BRANCH" + base="https://raw.githubusercontent.com/${GITHUB_REPOSITORY}/${ASSETS_BRANCH}/pr-$PR" + peak_img="$base/peak-relative.svg?$HEAD_SHA" + time_img="$base/time-relative.svg?$HEAD_SHA" + run_url="$GITHUB_SERVER_URL/$GITHUB_REPOSITORY/actions/runs/$GITHUB_RUN_ID" + { + echo '## Build cost β€” v1 vs legacy' + echo '' + echo '_v1 build peak & time relative to legacy, on this commit β€” not a comparison against master (that is CodSpeed)._' + echo '' + echo '| peak β€” v1 / legacy | time β€” v1 / legacy |' + echo '| :---: | :---: |' + echo "| \"peak | \"time |" + echo '' + echo '
Full table (time + peak, mean)' + echo '' + echo '```' + cat semantics-report/comparison.txt + echo '```' + echo '
' + echo '' + echo "πŸ“Š Interactive plots + CSV: download the **semantics-report-v1-vs-legacy** artifact from [this run]($run_url)." + echo '' + echo 'Report-only Β· not a gate Β· refreshed on every push Β· obsolete once legacy is dropped.' + } > report.md + - name: Sticky PR comment + if: >- + github.event_name == 'pull_request' && + github.event.pull_request.head.repo.full_name == github.repository + uses: marocchino/sticky-pull-request-comment@0ea0beb66eb9baf113663a64ec522f60e49231c0 + with: + header: semantics-report + path: report.md diff --git a/arithmetics-design/convention.md b/arithmetics-design/convention.md new file mode 100644 index 00000000..2c685ae2 --- /dev/null +++ b/arithmetics-design/convention.md @@ -0,0 +1,264 @@ +# The v1 convention + +The strict ("v1") convention for linopy. Goals and rollout plan: +[`goals.md`](goals.md). The bugs it fixes are catalogued in [#714]. + +An object-scope statement, then thirteen sections in three groups: absence +(Β§1–§7), coordinate alignment (Β§8–§11), then constraints and reductions +(Β§12–§13). + +## Object scope + +The convention governs every operation a linopy object takes part in, +whatever the other operand is β€” a DataArray, a pandas Series or DataFrame, a +numpy array, a list, or a scalar. A non-linopy operand is converted to a +labelled array (`as_dataarray` in `linopy.common`) and from there behaves +exactly like the constant-only expression holding the same values and +coordinates: `x + arr` builds what `x + arr_expr` builds, for every operator +and in either operand position. No rule below depends on what type an +operand arrived as. How an operand *gets* its labels β€” its own coordinates, +or pairing by size when it has none β€” is the alignment group's first rule. + +(The lone exception is type-decided: an expression is never a valid divisor, +so `x / arr` works where `x / arr_expr` raises `TypeError`.) + +## Absence + +Absence β€” a labelled slot the model does not cover β€” is the richer half of the +convention. The sections below say what it is (Β§1–§3), how it arises (Β§4–§5), +and how it flows through arithmetic and is resolved (Β§6–§7). + +### Β§1. Absence is a first-class state + +A *slot* β€” one labelled position β€” is either present or *absent*. An absent +slot is one the model does not cover. Absence is a state in its own right, +never a stand-in for a number: an absent variable is not a variable fixed to +zero ([#712]). + +### Β§2. Encoding absence + +The *marker* is how an absent slot is stored: `NaN` in floating-point fields +(`coeffs`, `const`, numeric constants), and `-1` in integer label fields (a +variable's `labels`, an expression's `vars`, which cannot hold a NaN). The two +encodings are one concept β€” an absent slot, whatever the dtype. + +Within a single slot, the markers move together: `const.isnull()` at a slot +implies *every* term at that slot has `coeffs = NaN` and `vars = -1`. Operators +that introduce absence at a slot also absorb any live terms there, so the +storage never carries a half-absent row. A term at a *present* slot may still +carry `vars = -1` after `fillna(value)` revives the slot β€” that's a *dead +term*, inert at the solver layer, and only meaningful as storage book-keeping. + +### Β§3. Testing absence + +`isnull()` is the one predicate for absence. It reads the marker β€” `NaN` or +`-1`, whichever the field uses β€” and reports absence slot by slot. Every rule +that speaks of an "absent slot" means exactly what `isnull()` reports; the +caller never inspects the raw marker. + +### Β§4. Creating absence + +Absence enters a model only through named operations: `mask=` at construction +marks slots absent up front; `.where(cond)` masks slots in place, keeping +shape; `.reindex()`, `.reindex_like()`, `.shift()`, and `.unstack()` +restructure a coordinate and leave the new positions absent. Operations that +merely move or select existing data β€” `.roll()`, `.sel()`, `.isel()` β€” never +introduce it. + +### Β§5. User-supplied NaN raises + +A NaN in a user-supplied constant raises `ValueError`. linopy trusts NaN only +from its own structural operations (Β§4), which genuinely mark absence. A NaN in +user data is ambiguous β€” a deliberate "absent", or a data error β€” so linopy +refuses to guess and asks the caller to resolve it with `fillna()`. This +replaces today's silent per-operator fills, which guessed a different value for +every operator ([#713]). To mark slots absent, use the mechanisms of Β§4 β€” a +bare NaN in a constant is not one of them. + +The alternative β€” reading user NaN as "absent" instead of raising β€” was +discussed in [#627] and closed: ambiguous overload of a numeric value +defeats goal #1, since a data-error NaN is silently re-labelled as +intentional absence. + +### Β§6. Absence propagates through every operator + +Every operator carries absence through unchanged: a slot absent in any operand +is absent in the result. `shifted * 3` is absent; `shifted + 5` is absent; +`x + shifted` is absent wherever `shifted` is β€” even though `x` itself is fine +there. + +linopy never fills an absent slot on the user's behalf, because the right fill +depends on intent it cannot see: 0 for a sum, 1 for a product, or "leave this +out" entirely. Because every operator propagates the same way, the algebraic +laws of Β§10 carry over to absent slots untouched β€” absence absorbs, so every +grouping of an expression agrees. And `shifted * 3` staying absent, rather than +collapsing to `0`, is what preserves the absent-vs-zero distinction of Β§1. + +### Β§7. Resolving absence + +Because Β§6 never fills, turning an absent slot into a value is the caller's +explicit act, never linopy's. `fillna(value)` fills an expression's absent +slots; `.fillna(...)` fills a constant before it enters the arithmetic; +`fill_value=` on a named method fills as part of the call. Filling at the call +site documents the intent: `x + y.shift(time=1).fillna(0)` says "treat the +missing earlier step as zero" exactly where it matters. + +## Coordinate alignment + +linopy's operands are xarray objects, so the convention starts from xarray's +alignment model (goal 4): coordinates align by *label*, never by position; +non-shared dimensions broadcast; a mismatch on a shared dimension is resolved +by an explicit *join*. + +Operands that carry coordinates β€” a DataArray, a pandas Series or DataFrame β€” +align by them, under the rules below. *Unlabeled* operands β€” numpy arrays, +lists, polars Series β€” carry no labels to align by, so they pair with the +linopy operand's dimensions by size: each axis adopts the dimension (and the +coordinates) whose length matches, and the rules below apply from there. The +pairing must be determined by the sizes alone. A length-4 array meeting a +variable with dims `(a: 4, time: 5)` pairs with `a`; meeting a variable with +dims `(a: 4, b: 4)` it could pair with either, so the operation raises β€” as +it does when no dimension matches. The same goes for a 4Γ—4 array against +`(a: 4, b: 4)`: sizes cannot tell `(a, b)` from `(b, a)`. To name the +dimensions, wrap the array in a DataArray. + +A scalar broadcasts over every dimension and so needs no pairing. A 0-d +array is treated as a scalar; a Python `list` is read as a numpy array +(it carries values, not labels). Implemented in `linopy.alignment` +([#736]). + +### Β§8. Shared dimensions must carry the same labels + +If two operands share a dimension, their coordinate labels must be the same +*set*, or the operator raises `ValueError`. Order is immaterial: the same +labels in a different order are the same coordinate and align by label (a +reindex), following "by label, never by position" above β€” only a difference in +the label set raises. + +This is close to xarray's `arithmetic_join="exact"` β€” deliberately stricter +than xarray's own default (`inner`) β€” but order-independent, where xarray's +`exact` would reject a pure reorder. An inner join silently drops the +non-overlapping labels, and in an optimization model a dropped coordinate is a +dropped term or constraint: a silent wrong answer. Matching on the label set +surfaces a real mismatch where it happens. (The [pyoframe] library uses the +same model.) + +Because the rule is identical for every operator, the operator-alignment split +([#708]) β€” `*` aligning by label while `+`, `-`, `/` go by position β€” +disappears. + +### Β§9. Non-shared dimensions broadcast freely + +A dimension present in only one operand broadcasts over the other, with no +restriction β€” for both expressions and constants. Only *shared* dimensions are +subject to Β§8. + +### Β§10. Mismatches resolve via an explicit join + +When coordinates genuinely differ, Β§8 raises β€” and the caller says how to +resolve it. Several primitives bring operands into agreement: + +- `.sel()` / `.isel()` cut operands down to a shared subset β€” often the + clearest fix. +- The named methods β€” `.add` `.sub` `.mul` `.div` `.le` `.ge` `.eq` β€” take a + `join=` argument: `exact`, `inner`, `outer`, `left`, `right`, or `override`. + `override` is the old positional behavior β€” still available, but now opt-in + and named rather than triggered by a size coincidence. It still requires the + shared dimensions to match in size: a genuine size mismatch raises rather + than relabelling mismatched data, so reach for a label join (`inner` / + `outer` / `left` / `right`) when the sizes really differ. +- `.reindex()` / `.reindex_like()` conform an operand to a target index + (extending past the original creates absent positions β€” Β§4). +- `.assign_coords()` relabels an operand outright. Unlike `join="override"`, + which checks the shared dims match in size, this is an *unguarded* relabel: + it renames positions whether or not they correspond, so the "made explicit" + here is the caller's responsibility, not a safety check. +- `linopy.align()` pre-aligns several operands at once. + +Because no operator silently drops coordinates, the associativity break +([#711]) cannot occur: the operation that used to drop coordinates now raises. +Every standard algebraic law β€” commutativity, associativity, distributivity, +the identities β€” holds for same-coordinate operands. + +### Β§11. Auxiliary-coordinate conflicts raise + +Auxiliary (non-dimension) coordinates are user-attached metadata: a coord +defined on some dimension but not itself a dimension, like a `B(A)` group +label on dimension `A`. linopy *validates* them (the conflict-raise rule +below) and *propagates* them through arithmetic unchanged, but never +*computes* with them β€” they describe the data, they don't enter the math. + +When two operands carry an aux coord with the same name and values agree, +the coord propagates to the result. When only one operand carries the +coord, it propagates from that operand unchanged β€” asymmetric presence is +not a conflict. When the values *do* disagree (same name on both sides, +different values), the operator raises β€” `xarray` silently drops the +conflict, which is the [#295] bug. The caller resolves it explicitly with +`.drop_vars(name)` (remove the coord) or `.assign_coords(name=...)` +(relabel one side). + +**Stacked MultiIndex dimensions.** A stacked MultiIndex dim (e.g. PyPSA's +`(period, timestep)` `snapshot`) stores its *levels* as auxiliary +coordinates β€” `period` and `timestep` are non-dimension coords on +`snapshot` β€” and its elements are *level combinations* (one tuple per +position). An operand whose *dimension* names one of those levels β€” a +per-`period` weighting meeting a `snapshot`-indexed expression β€” is a +same-name conflict between a dimension and an auxiliary coordinate, and it +raises like any other conflict of this section. There is no implicit +projection; write it explicitly by selecting with the dimension's level +values: + + weights.sel(period=expr.indexes["snapshot"].get_level_values("period")) + +An input that reconstructs the *entire* MultiIndex (all levels, every +combination) is not a conflict β€” it is the same coordinate spelled +differently, and aligns by tuple under Β§8, in any order. Order is +immaterial here exactly as for a plain dimension: the same tuples in a +different order are reordered to match, not rejected. + +(Legacy projects implicitly and warns β€” scenario B of the [#732]/[#737] +discussion; the implicit projection is removed at 1.0.) + +## Constraints and reductions + +Two kinds of operation build on the rules above without being binary operators: +the comparisons that form constraints, and the reductions that collapse a +dimension. + +### Β§12. Constraints follow the same rules + +A constraint is built by comparing two sides with `<=`, `>=`, or `==` β€” and a +comparison is an operator like any other. It aligns its sides by Β§8 and carries +absence by Β§6, exactly as `+`, `-`, `*`, and `/` do. So algebraically equal +forms build the same constraint: `x - a <= 0` and `x <= a` agree, where today +they do not ([#707]). + +Each slot becomes one constraint row. An absent slot yields no row β€” absence +propagated into a comparison drops the constraint there, the same outcome as +masking it. + +### Β§13. Reductions skip absent slots + +Reductions collapse a dimension rather than combining two operands, so the +NaN propagation of Β§6 does not apply: they *skip* absent slots instead. `sum` +(including `groupby.sum`) adds only the present terms, and the sum of none is +the zero expression. The objective totals its terms the way `sum` does. + +Further reductions (`mean`, `resample`, `coarsen`) are not in linopy yet; they +are added under v1 only ([#703]) β€” as new operations with no legacy behaviour β€” +and follow this same skip-absent rule. + + +[pyoframe]: https://github.com/Bravos-Power/pyoframe +[#732]: https://github.com/PyPSA/linopy/pull/732 +[#737]: https://github.com/PyPSA/linopy/pull/737 +[#736]: https://github.com/PyPSA/linopy/issues/736 +[#714]: https://github.com/PyPSA/linopy/issues/714 +[#703]: https://github.com/PyPSA/linopy/issues/703 +[#713]: https://github.com/PyPSA/linopy/issues/713 +[#712]: https://github.com/PyPSA/linopy/issues/712 +[#711]: https://github.com/PyPSA/linopy/issues/711 +[#708]: https://github.com/PyPSA/linopy/issues/708 +[#707]: https://github.com/PyPSA/linopy/issues/707 +[#627]: https://github.com/PyPSA/linopy/issues/627 +[#295]: https://github.com/PyPSA/linopy/issues/295 diff --git a/arithmetics-design/docs-plan.md b/arithmetics-design/docs-plan.md new file mode 100644 index 00000000..9c55750e --- /dev/null +++ b/arithmetics-design/docs-plan.md @@ -0,0 +1,32 @@ +# Docs plan β€” user-facing migration guide + +Early-stage outline for the v1 migration docs. Not the guide itself β€” +the three pieces it needs to cover, written when someone picks it up. + +## Three audiences, one migration + +- **Downstream library maintainers** (PyPSA, pypsa-eur, calliope, …) β€” + carry the bulk of the migration work: opt their codebases into v1, + fix the raises, ship a release that no longer warns under legacy. +- **Direct users of linopy** β€” write linopy code themselves and need + to know what changes for their own call sites. +- **End users of downstream libraries** β€” never touch linopy directly, + but may see a `LinopySemanticsWarning` in CI logs and need a + pointer to "this is upstream; your maintainer will handle it". + +## Three things to cover + +1. **Why v1 exists.** One paragraph: legacy silently mishandled NaN, + coord mismatches, and absent variables. The bug catalogue in #714 + has the case-by-case detail. + +2. **What's changing and when.** The rollout timeline: + - v1 ships opt-in via `linopy.options['semantics'] = 'v1'`. + - v1 becomes the default in a later minor release (date TBD). + - Legacy removed at 1.0. + +3. **How to migrate.** What downstream maintainers do to flip their + codebase: opt in on a branch, run tests, fix the raises. The + legacy warning text already names the rule and the fix per site, + so the guide is mostly the high-level recipe plus a pointer to + the spec (`arithmetics-design/convention.md`) for the rule list. diff --git a/arithmetics-design/goals.md b/arithmetics-design/goals.md new file mode 100644 index 00000000..ba434836 --- /dev/null +++ b/arithmetics-design/goals.md @@ -0,0 +1,47 @@ +# The v1 convention β€” design & transitioning goals + +Goals for linopy's strict ("v1") convention. The bugs that motivate +it are catalogued in [#714]; the convention itself is in +[`convention.md`](convention.md). + +## Design goals + +The convention serves four goals, in priority order: + +1. **No silent wrong answers.** Every bug in the catalogue ([#714]) returns a + plausible result with no error. The overriding goal: a mismatch linopy + cannot resolve unambiguously must raise, not get guessed. Where the library + cannot decide, the caller does β€” with an explicit join, `.sel()`, or + `fill_value=`. +2. **Preserve the algebraic laws.** Commutativity, associativity, + distributivity, the identities. Optimization code builds expressions by + rearranging terms, and the convention must keep that safe. +3. **Absence is first-class.** A variable can be genuinely absent at a slot β€” + masked out, or shifted past the edge. The data model needs an explicit + marker for that absence, kept distinct from a zero term, so absent-vs-zero + is never a silent guess. +4. **Least surprise.** linopy is built on xarray and its users know xarray. The + convention should behave the way xarray already taught them β€” align by + label, broadcast non-shared dimensions, resolve mismatches with a named + join β€” not invent linopy-specific rules. Auxiliary coordinates the user + attached are the user's; linopy validates and carries them through, + never silently dropped or rewritten. + +## Transitioning goals + +1. **Non-breaking.** Existing code keeps working β€” legacy stays available and + unchanged until it is removed at linopy 1.0. +2. **Actionable warnings.** Warn every legacy user about behaviour changes β€” + what changes under v1, and how to fix it β€” aiming for 100% coverage. +3. **No silent change.** Opting into v1 never silently changes a model β€” every + difference is either raised, or was warned about in legacy mode. + +**Schedule:** + +1. Introduce v1 as opt-in β€” warn about behaviour changes on legacy, raise if + opted into v1. +2. Make v1 the default, allow opt-out. +3. linopy 1.0 β€” drop the legacy convention entirely. + + +[#714]: https://github.com/PyPSA/linopy/issues/714 diff --git a/arithmetics-design/legacy-removal.md b/arithmetics-design/legacy-removal.md new file mode 100644 index 00000000..5d563044 --- /dev/null +++ b/arithmetics-design/legacy-removal.md @@ -0,0 +1,165 @@ +# Legacy removal checklist (for linopy 1.0) + +The v1 convention ships alongside legacy from 0.x onward and replaces it +entirely at 1.0. This file enumerates everything to delete when that +release happens, in dependency order. Most edits are mechanical; +`grep "LEGACY: remove at 1.0"` finds every inline marker comment in the +source tree. + +## Implementation + +### `linopy/config.py` + +- Drop `LEGACY_SEMANTICS`, `V1_SEMANTICS`, `VALID_SEMANTICS` constants. +- Drop `LEGACY_SEMANTICS_MESSAGE`. +- Drop the `LinopySemanticsWarning` class. +- Remove the `semantics` key from `options`. The option no longer exists; + callers don't need to opt in. +- Remove the v1/legacy validation branch in `set_value`. + +### `linopy/semantics.py` + +- Delete `is_v1()` (always true). Inline `True` at the four import sites + in `expressions.py`/`variables.py` or, better, delete the import and + the now-dead `else` branches alongside it. +- Drop the legacy-warn branch from `check_user_nan_scalar` / + `check_user_nan_array` β€” both become a single `raise ValueError(...)`. +- Delete `dim_coords_differ`: only the legacy `_align_constant` default + path uses it. +- `merge_shared_user_coords_differ`, `conflicting_aux_coord`, + `absorb_absence` stay β€” they're v1 enforcement helpers. + +### `linopy/expressions.py` + +- `_add_constant`: delete `_add_constant_legacy`; inline + `_add_constant_v1` into the now-trivial dispatcher (or rename it to + `_add_constant`). +- `_apply_constant_op`: same treatment with `_apply_constant_op_legacy`. +- `_align_constant`: delete the `else` branch under `if join is None:` + that handles the legacy size-aware default (`other.sizes == self.const.sizes` + positional + `reindex_like` left-join paths). The explicit-join code + below stays. +- `to_constraint`: drop the `if is_v1(): ... return ...` wrapper and + keep its body; delete the legacy auto-mask fallthrough that follows + (the `rhs_nan_mask` plumbing plus the `rhs.reindex_like(..., + fill_value=np.nan)` pad). +- `LinearExpression.isnull`: drop the legacy `(self.vars == -1).all(...) + & self.const.isnull()` branch β€” `self.const.isnull()` is the v1 answer. +- `merge`: + - Drop the `if differ: warn(...)` line and the `if aux_conflict: + warn(...)` line β€” these are the Β§8 / Β§11 legacy warns. The raises + above stay. + - The `skipna = not is_v1()` simplifies to `skipna = False` (v1's + propagation rule). + - The trailing `if is_v1(): ds = absorb_absence(ds)` becomes + unconditional. +- Drop the `LinopySemanticsWarning` / `LEGACY_SEMANTICS_MESSAGE` imports + from `expressions.py`. + +### `linopy/variables.py` + +- `Variable.to_linexpr`: drop the `else` branch (legacy + `reindex_like(fill_value=0).fillna(0)`); make the v1 `reindex_like( + fill_value=NaN) β†’ .where(~absent)` path unconditional. The + `const = NaN`/`0` assign also becomes unconditional. +- Drop the `from linopy.semantics import is_v1` import. + +### `linopy/alignment.py` + +- `_enforce_implicit_projections`: drop the legacy `warn_legacy(...)` + branch β€” the v1 raise for partial-level / coverage-gap projections + becomes unconditional. The projection machinery itself + (`_project_onto_multiindex_levels`, `_LevelProjection`) stays: + full-coverage full-level projections remain legal under v1 (they are + the same coordinate spelled differently, Β§8). +- `_dims_for_unlabeled_operand`: drop the legacy positional-pairing + fallback (the `warn_legacy(...)` branches plus the `return + list(candidates)`); the v1 size-pairing β€” the `is_v1()` block that + raises on ambiguity / no-match β€” becomes the whole function. The + `as_constant` / `_pair_axes_by_size` helpers stay (v1-clean). + +### `linopy/piecewise.py` / `linopy/sos_reformulation.py` + +Nothing to remove; these are v1-clean (the `drop=True` / `assign_coords` +fixes from Slice P are correct for both semantics). + +## Tests + +### `test/conftest.py` + +- Drop the `LEGACY_SEMANTICS`/`V1_SEMANTICS`/`VALID_SEMANTICS` imports + and the `LinopySemanticsWarning` import. +- Drop the `legacy` / `v1` marker registration in `pytest_configure`. +- Delete the autouse `semantics` fixture entirely (no more parameterization, + no more warning suppression). + +### `test/test_legacy_violations.py` + +- Delete the file. Everything in it either documents legacy behaviour + (gone) or tests v1 raises (covered by the per-module test files we + add tests to alongside the implementation). + + Before deleting, move any v1 tests that don't have a per-Β§ home into + the appropriate module: + - `TestExactAlignmentConstant`, `TestExactAlignmentMerge`, + `TestBroadcastNonSharedDim` β†’ `test_linear_expression.py`. + - `TestConstraintRHS` β†’ `test_constraints.py`. + - The rest are small enough to fold in alongside related tests. + +### `test/test_convention.py` + +- Delete the file (it tests the `options["semantics"]` framework, which + is gone). + +### Marker stripping + +`grep -rn "@pytest.mark.legacy\|@pytest.mark.v1\|pytestmark = pytest.mark.legacy" +test/` finds every marker: + +- `@pytest.mark.legacy` decorators β€” delete the decorator (the test + body is documenting old behaviour; deleting the whole test is + usually right). Spot-check before each delete; a few "legacy" + marks turned out to gate on legacy auto-mask semantics and the + test itself stays valid under v1. +- `@pytest.mark.v1` decorators β€” strip the decorator (the test stays). +- `pytestmark = pytest.mark.legacy` at module level β€” was only used + while `piecewise.py` was non-v1-aware; removed in Slice P. Verify + none remain. + +## Documentation + +### `arithmetics-design/goals.md` + +- Drop the entire "Transitioning goals" section (the three transitioning + goals + the schedule are about the legacy bridge, which is gone). +- Goal #4's mention of `LinopySemanticsWarning` (if any) goes. + +### `arithmetics-design/convention.md` + +- Drop the `## Legacy` framing if it exists. Each Β§ already describes + v1 directly; any "where today this does X" asides referencing legacy + behaviour can go. +- Update the intro: "The strict ('v1') convention" β†’ "The arithmetic + convention" (drops the v1 framing now that there's only one). + +### This file (`legacy-removal.md`) + +- Delete after the 1.0 release ships. + +## Order of operations + +A safe sequence (each step compiles and tests pass): + +1. Delete legacy test infrastructure (`test/conftest.py` fixture, + `test/test_convention.py`, `test/test_legacy_violations.py` after + moving v1 tests out). +2. Strip `@pytest.mark.legacy` decorators (the tests fail under v1 + anyway once the legacy paths are gone β€” delete or update each). +3. Delete legacy implementation branches in `expressions.py` / + `variables.py`. +4. Delete `semantics.py` legacy bits (`is_v1`, the warn branches in + `check_user_nan_*`, `dim_coords_differ`). +5. Delete `config.py` symbols (`LEGACY_SEMANTICS`, the warning class, + the option key). +6. Update `arithmetics-design/goals.md` and `convention.md`. +7. Delete this file. diff --git a/arithmetics-design/open-items.md b/arithmetics-design/open-items.md new file mode 100644 index 00000000..2c94d407 --- /dev/null +++ b/arithmetics-design/open-items.md @@ -0,0 +1,70 @@ +# Open items for the v1 convention + +The convention itself is **fully specified and implemented**: +[`convention.md`](convention.md) Β§1–§13 plus the object-scope and +coordinate-alignment intro rules all have a v1 implementation. The high-level +schedule lives in [`goals.md`](goals.md) (opt-in β†’ default β†’ 1.0); this file +tracks the concrete items per stage. + +One open **design** decision remains β€” [#744] MultiIndex storage. No open +arithmetic-rule questions. + +## Stage 1 β€” release v1 (opt-in) + +Legacy stays the default. The transition surface must already be **complete** +here ([`goals.md`] step 1: *warn on legacy, raise on v1*). + +- [ ] **Resolve [#744] β€” MultiIndex storage β€” before #717 merges.** First-class + `pd.MultiIndex` vs. a flat dim + auxiliary level coords. Β§11's + stacked-MultiIndex rule and the storage of a `snapshot`-style + `(period, timestep)` dim both depend on it β€” and Β§11 **ships in #717**, so + changing the model after v1 is released would change observable Β§11 behaviour. + Must be locked (and its implementation in the same release cut) before #717 + ships v1, not deferred to the default flip. +- [ ] Land [#717] (v1 semantics) β†’ `master` +- [ ] **Transition surface complete** β€” every behaviour-change site raises under + v1 *and* warns under legacy (`warn_legacy`, naming the fix). This is the + "no silent change" guarantee ([`goals.md`] transitioning goal #3): shipping v1 + with a gap would silently change any model that opts in. +- [ ] Changelog note β€” v1 available via `options['semantics'] = 'v1'`; legacy + remains the default; link [`convention.md`]. + +## Stage 2 β€” make v1 the default (legacy opt-out) + +- [ ] Write the **migration guide** ([`docs-plan.md`] is only an outline). +- [ ] Flip the `options['semantics']` default to v1. +- [ ] Close [#714] once v1 ships as default. + +## Stage 3 β€” linopy 1.0, remove legacy + +- [ ] The **strip**: the concentrated MI machinery (the `alignment.py` + level-projection subsystem, netcdf MI (de)serialize) *and* the scattered + surface (`assign_multiindex_safe` Γ—39, `isinstance(MultiIndex)` guards, MI + branches) β€” plus every other `grep "LEGACY: remove at 1.0"` marker. + Dependency-ordered checklist in [`legacy-removal.md`]. +- [ ] Reframe [`convention.md`] / [`goals.md`] β€” drop the "v1"/legacy framing + once there is only one convention. + +## Not blocking (follow-ups) + +- [ ] Test cleanup β€” public-API assertions (`.indexes` / `.coords` / `.sizes` / + `.coord_dims`) instead of internal `.data` / `.coeffs.coords`; assert the + **full** error/warning text in v1-raise and legacy-warn tests; de-dup repeated + index/setup fixtures. +- [x] Pin the **legacy** side of every v1/legacy result divergence. + +## Decisions + +- [ ] **[#744] β€” MultiIndex storage** β€” the one open design decision; gates the + **v1 release** (Β§11 ships in #717 and depends on it), not just the default flip. +- [ ] **When v1 becomes the default** β€” pick the release; gated on the migration + guide. + + +[#744]: https://github.com/PyPSA/linopy/issues/744 +[#714]: https://github.com/PyPSA/linopy/issues/714 +[#717]: https://github.com/PyPSA/linopy/pull/717 +[`goals.md`]: goals.md +[`docs-plan.md`]: docs-plan.md +[`legacy-removal.md`]: legacy-removal.md +[`convention.md`]: convention.md diff --git a/benchmarks/conftest.py b/benchmarks/conftest.py index b9ef6014..f873279d 100644 --- a/benchmarks/conftest.py +++ b/benchmarks/conftest.py @@ -2,6 +2,8 @@ from __future__ import annotations +import os +from collections.abc import Iterator from typing import TYPE_CHECKING import pytest @@ -12,6 +14,33 @@ import linopy from benchmarks.registry import BenchSpec + +@pytest.fixture(autouse=True) +def _bench_semantics() -> Iterator[None]: + """ + Force the arithmetic semantics for the whole run from + ``LINOPY_BENCH_SEMANTICS`` (``legacy`` / ``v1``), restoring after. + + The memory-report job runs the build twice β€” once at the default (legacy), + once with this set to ``v1`` β€” to A/B the v1 cost. Leaving it unset keeps the + default, so the node ids the legacy run and the CodSpeed run produce are + identical (``test_build[…]``, no suffix): CodSpeed keeps its per-id history + against master, and the two report JSONs line up under ``benchmem compare``. + """ + import linopy + + mode = os.environ.get("LINOPY_BENCH_SEMANTICS") + if not mode: + yield + return + old = linopy.options["semantics"] + linopy.options["semantics"] = mode + try: + yield + finally: + linopy.options["semantics"] = old + + # Test modules the CodSpeed instruments measure (edit to change coverage). # build + the two export paths: to_lp (LP text) and to_solver (direct handoff, # which also exercises matrix-gen). matrices is dropped β€” a subset of to_solver; diff --git a/benchmarks/drivers/test_pypsa_carbon_management.py b/benchmarks/drivers/test_pypsa_carbon_management.py index 209416ba..77185337 100644 --- a/benchmarks/drivers/test_pypsa_carbon_management.py +++ b/benchmarks/drivers/test_pypsa_carbon_management.py @@ -12,6 +12,12 @@ pypsa = pytest.importorskip("pypsa") +@pytest.fixture(autouse=True) +def _skip_under_v1(_bench_semantics: None) -> None: + if lp.options["semantics"] == "v1": + pytest.skip("PyPSA emits NaN-valued constants that v1 rejects by design.") + + @pytest.fixture(scope="module") def network() -> Any: try: diff --git a/benchmarks/models/expression_arithmetic.py b/benchmarks/models/expression_arithmetic.py index 0d5af581..706af425 100644 --- a/benchmarks/models/expression_arithmetic.py +++ b/benchmarks/models/expression_arithmetic.py @@ -3,6 +3,7 @@ from __future__ import annotations import numpy as np +import xarray as xr import linopy from benchmarks.registry import BenchSpec, register @@ -20,7 +21,11 @@ def build_expression_arithmetic(n: int) -> linopy.Model: z = m.add_variables(coords=[range(n)], dims=["j"], name="z") # Expression arithmetic: broadcasting y (dim i) and z (dim j) against x (dim i,j) - coeffs = np.linspace(-1, 1, n * n).reshape(n, n) + coeffs = xr.DataArray( + np.linspace(-1, 1, n * n).reshape(n, n), + dims=["i", "j"], + coords={"i": range(n), "j": range(n)}, + ) expr1 = x * coeffs + y - z expr2 = 2 * x - 3 * y + z combined = expr1 + expr2 diff --git a/benchmarks/models/milp.py b/benchmarks/models/milp.py index f6058cc8..f3782ffb 100644 --- a/benchmarks/models/milp.py +++ b/benchmarks/models/milp.py @@ -20,6 +20,7 @@ from __future__ import annotations import numpy as np +import xarray as xr import linopy from benchmarks.registry import ( @@ -38,7 +39,11 @@ def build_milp(n: int) -> linopy.Model: cap = 100.0 # capacity per module Y_MAX = 5 # max modules per facility - transport = rng.uniform(1, 20, size=(n, n)) # per-unit shipping cost + transport = xr.DataArray( + rng.uniform(1, 20, size=(n, n)), + dims=["facility", "customer"], + coords={"facility": facilities, "customer": customers}, + ) fixed = rng.uniform(50, 200, size=n) # cost per facility module demand = rng.uniform(20, 80, size=n) # demand at each customer diff --git a/linopy/__init__.py b/linopy/__init__.py index b813f71d..10251de2 100644 --- a/linopy/__init__.py +++ b/linopy/__init__.py @@ -13,7 +13,7 @@ # we need to extend their __mul__ functions with a quick special case import linopy.monkey_patch_xarray # noqa: F401 from linopy.alignment import align -from linopy.config import options +from linopy.config import LinopySemanticsWarning, options from linopy.constants import ( EQUAL, GREATER_EQUAL, @@ -57,6 +57,7 @@ "GREATER_EQUAL", "LESS_EQUAL", "LinearExpression", + "LinopySemanticsWarning", "Model", "Objective", "OetcHandler", diff --git a/linopy/alignment.py b/linopy/alignment.py index 13126b85..1beb718e 100644 --- a/linopy/alignment.py +++ b/linopy/alignment.py @@ -24,7 +24,6 @@ from collections.abc import Callable, Hashable, Iterable, Mapping, Sequence from functools import partial from typing import TYPE_CHECKING, Any, Literal, NamedTuple, overload -from warnings import warn import numpy as np import pandas as pd @@ -42,11 +41,32 @@ # Added in xarray 2025.6.0; it subclasses ValueError on newer versions. CoordinateValidationError = ValueError # type: ignore[assignment, misc] -from linopy.constants import ( - HELPER_DIMS, - EvolvingAPIWarning, -) -from linopy.types import CoordsLike, DimsLike +from linopy.constants import HELPER_DIMS +from linopy.types import UNLABELED_TYPES, CoordsLike, DimsLike + + +def as_constant(other: Any) -> Any: + """ + Normalize a degenerate operand for arithmetic on entry. + + Two normalizations let the operators treat every numeric operand the + same way downstream: + + - a Python ``list`` carries array data but no numeric operators + (``-[1, 2]`` is a ``TypeError``, ``[1, 2] * x`` repeats the list), so + it becomes a numpy array; + - a 0-d numpy array (``np.array(1)``) is unwrapped to a Python scalar so + it takes the scalar fast-path instead of size-pairing. + + Everything else passes through unchanged β€” typed constants, DataArrays, + Variables, and Expressions all already behave. + """ + if isinstance(other, list): + other = np.asarray(other) + if isinstance(other, np.ndarray) and other.ndim == 0: + return other.item() + return other + if TYPE_CHECKING: from linopy.expressions import LinearExpression, QuadraticExpression @@ -390,10 +410,13 @@ def as_dataarray( if isinstance(arr, np.number): arr = float(arr) if dims is None: + # A scalar broadcasts over the coords' dims, but never over a + # helper dim (e.g. ``_term``) β€” those are storage book-keeping, + # not user axes. if isinstance(coords, Coordinates): - dims = coords.dims + dims = [d for d in coords.dims if d not in HELPER_DIMS] elif is_dict_like(coords) and np.ndim(arr) == 0: - dims = list(coords.keys()) + dims = [d for d in coords.keys() if d not in HELPER_DIMS] arr = DataArray(arr, coords=coords, dims=dims, **kwargs) elif not isinstance(arr, DataArray): @@ -512,18 +535,24 @@ def _project_onto_multiindex_levels( return arr, projections -def _warn_implicit_projections(projections: list[_LevelProjection]) -> None: +def _enforce_implicit_projections(projections: list[_LevelProjection]) -> None: """ - Deprecation warnings for implicit MultiIndex-level projections. + Semantics policy for implicit MultiIndex-level projections. - The same check in every mode (scenario B of the #732 / #737 discussion): - implicit projection is deprecated and raises under the v1 convention. The - strict path raises on coverage gaps before reaching here, so only partial - levels warn there; the non-strict path warns for both. + Implicit projection is legacy-only behavior (scenario B of the #732 / + #737 discussion): under legacy semantics it emits a deprecation warning + (#738: via ``warn_legacy``); under the v1 convention it raises β€” a dim + naming a level of a stacked MultiIndex dim is a shared-dim / aux-coord + concern (sections 8 and 11), and the projection must be written + explicitly by the caller. - TODO(#738): migrate to ``warn_legacy()`` / ``LinopySemanticsWarning`` - once the v1 semantics infrastructure (#717) lands. + The strict path raises on coverage gaps before reaching here, so only + partial levels arrive there; the non-strict path sees both. """ + # Deferred import: linopy.semantics imports xarray/pandas machinery that + # in turn may import this module's consumers; keep the seam lazy. + from linopy.semantics import is_v1, warn_legacy + for p in projections: if p.is_partial or p.has_gap: kind = ( @@ -532,116 +561,240 @@ def _warn_implicit_projections(projections: list[_LevelProjection]) -> None: else f"filling uncovered level combinations with NaN " f"(from level(s) {p.levels})" ) - warn( + if is_v1(): + raise ValueError( + f"multiindex-projection: implicitly {kind} onto MultiIndex " + f"dimension {p.dim!r} is not supported under the v1 " + f"convention (sections 8 and 11). Project the input onto " + f"the dimension explicitly, e.g. select with the " + f"dimension's level values." + ) + warn_legacy( f"multiindex-projection: implicitly {kind} onto MultiIndex " f"dimension {p.dim!r}. This is deprecated and will raise under " f"the v1 convention; project the input onto the dimension " f"explicitly (select with the dimension's level values) to " - f"keep current behavior.", - EvolvingAPIWarning, - stacklevel=3, + f"keep current behavior." ) -def _broadcast_to_coords( - arr: Any, - coords: CoordsLike | None = None, - dims: DimsLike | None = None, - **kwargs: Any, -) -> tuple[DataArray, list[_LevelProjection]]: +def _pair_axes_by_size( + shape: tuple[int, ...], sizes: dict[Hashable, int] +) -> tuple[list[Hashable] | None, str | None]: """ - Convert ``arr`` and broadcast it against ``coords`` (shared mechanics). + Pair each axis of an unlabeled array with the operand dim of matching size. - Returns the broadcast DataArray together with the MultiIndex-level - projections performed along the way, so the public entry points can - apply their own policy (warn or raise) to partial projections and - coverage gaps. + The pairing must be determined by the sizes alone (v1 convention, + coordinate-alignment intro): every axis size must match exactly one + operand dim, and no two axes may share a size. Returns + ``(dims, None)`` on success or ``(None, problem)`` where ``problem`` + describes why the pairing is impossible or ambiguous. """ - if coords is None: - return as_dataarray(arr, coords, dims, **kwargs), [] + by_size: dict[int, list[Hashable]] = {} + for d, n in sizes.items(): + by_size.setdefault(n, []).append(d) + + axes_per_size: dict[int, int] = {} + for s in shape: + axes_per_size[s] = axes_per_size.get(s, 0) + 1 + + for s, n_axes in axes_per_size.items(): + candidates = by_size.get(s, []) + if len(candidates) < n_axes: + return None, ( + f"no unambiguous dimension match for an axis of length {s}: " + f"the operand has dimensions {dict(sizes)}." + ) + if len(candidates) > 1 or n_axes > 1: + return None, ( + f"axis of length {s} could pair with any of " + f"{sorted(candidates, key=str)} β€” sizes alone cannot decide." + ) - expected = _coords_to_dict(coords, dims=dims) - if not expected: - return as_dataarray(arr, coords, dims, **kwargs), [] + return [by_size[s][0] for s in shape], None - if isinstance(arr, pd.Series | pd.DataFrame): - converted = _named_pandas_to_dataarray(arr) - if converted is not None: - arr = converted - if not isinstance(arr, DataArray): - # numpy/polars/unnamed-pandas inputs are positional β€” their only - # meaningful information is the values; any axis labels are - # auto-generated. Default dims to coords' keys so the conversion - # labels axes correctly (instead of dim_0/dim_1), then re-assign - # coords from expected so positional inputs align to coords by - # position. A shape mismatch surfaces here as a clear xarray - # "conflicting sizes" error rather than a confusing - # "coordinates do not match" further down. - if dims is None: - dims = list(expected) - arr = as_dataarray(arr, coords, dims=dims, **kwargs) - # Skip MultiIndex dims β€” re-assigning a PandasMultiIndex coord emits - # a FutureWarning and isn't needed (the conversion already used it). - arr = arr.assign_coords( - { - d: expected[d] - for d in arr.dims - if d in expected and not isinstance(arr.indexes.get(d), pd.MultiIndex) - } +def _dims_for_unlabeled_operand( + shape: tuple[int, ...], expected: dict[Hashable, Any] +) -> list[Hashable]: + """ + Choose dim names for an unlabeled (numpy / list / polars) input. + + Used everywhere an unlabeled array meets a known set of dims β€” bounds + and masks in ``add_variables`` / ``add_constraints``, and arithmetic + operands (#736). + + v1 (convention, coordinate-alignment intro): axes pair with the dims by + size; ambiguity or a missing match raises, with wrap-in-a-DataArray as + the documented resolution. Legacy: axes pair with the leading dims + positionally; a deprecation warning fires whenever the v1 pairing would + differ from or reject the positional one. + """ + from linopy.semantics import is_v1, warn_legacy + + # A 0-d operand has no axes to pair β€” it broadcasts over every dim, so it + # carries no dim names (matching a bare scalar). + if len(shape) == 0: + return [] + + # Helper dims (e.g. ``_term``) are storage book-keeping, never user axes, + # so they are not pairing candidates. + candidates = {d: v for d, v in expected.items() if d not in HELPER_DIMS} + sizes = {d: len(_as_index(v)) for d, v in candidates.items()} + paired, problem = _pair_axes_by_size(shape, sizes) + positional = list(candidates)[: len(shape)] + + if is_v1(): + if problem is not None: + raise ValueError( + f"Cannot pair an unlabeled array of shape {tuple(shape)} with " + f"the operand's dimensions: {problem} Wrap the array in an " + f"xarray.DataArray with explicit dims to name its axes." + ) + assert paired is not None + return paired + + # LEGACY: remove at 1.0 β€” positional pairing plus the transition warning. + if problem is not None: + warn_legacy( + f"An unlabeled array of shape {tuple(shape)} was paired with the " + f"operand's leading dimension(s) {positional} by position. Under " + f"the v1 convention this raises: {problem} Wrap the array in an " + f"xarray.DataArray with explicit dims to keep it working." + ) + elif paired != positional: + warn_legacy( + f"An unlabeled array of shape {tuple(shape)} was paired with the " + f"operand's leading dimension(s) {positional} by position. Under " + f"the v1 convention it pairs by size instead β€” with {paired} β€” " + f"which gives a different result. Wrap the array in an " + f"xarray.DataArray with explicit dims to make the pairing explicit." ) + return positional - arr, projections = _project_onto_multiindex_levels(arr, expected) +def _matmul_operand_to_dataarray( + other: Any, coords: Coordinates, coord_dims: tuple[Hashable, ...] +) -> DataArray: + """ + Convert a non-expression ``@`` operand, pairing unlabeled axes by size. + + Shared by ``LinearExpression.__matmul__`` and + ``QuadraticExpression.__matmul__``: :func:`_dims_for_positional_input` + decides which dims the contraction collapses (#736) β€” by size under v1, + positionally with a warning under legacy. Unlike the broadcast pipeline + this only converts (no reindex / expand / transpose). + """ + expected = {d: coords[d] for d in coord_dims} + dims = _dims_for_positional_input(other, expected, None) + return as_dataarray(other, coords=coords, dims=dims) + + +def _dims_for_positional_input( + arr: Any, expected: dict[Hashable, Any], dims: DimsLike | None +) -> DimsLike | None: + """ + Resolve the dim names a non-DataArray input's axes adopt. + + An explicit ``dims`` is honored as given. Otherwise an unlabeled + array (numpy / list / polars) pairs its axes with ``expected`` by + size (#736); any other input falls back to the coords dims, minus + helper dims like ``_term`` which are never user axes. + """ + if dims is not None: + return dims + if isinstance(arr, UNLABELED_TYPES) and np.ndim(arr) >= 1: + return _dims_for_unlabeled_operand(np.shape(arr), expected) + return [d for d in expected if d not in HELPER_DIMS] + + +def _label_input( + arr: Any, expected: dict[Hashable, Any], dims: DimsLike | None, **kwargs: Any +) -> DataArray: + """ + Convert a non-DataArray input to a DataArray labelled against ``expected``. + + The converter is handed ``expected`` (the normalized nameβ†’values dict), + not the raw sequence-form coords, so it selects coords by name β€” a + sequence would zip dims to coords by position, which is wrong once + size-pairing has chosen a non-leading dim. + """ + dims = _dims_for_positional_input(arr, expected, dims) + arr = as_dataarray(arr, expected, dims=dims, **kwargs) + # Re-assign non-MultiIndex coords from ``expected`` (a MultiIndex coord + # re-assignment emits a FutureWarning and the conversion already used it). + return arr.assign_coords( + { + d: expected[d] + for d in arr.dims + if d in expected and not isinstance(arr.indexes.get(d), pd.MultiIndex) + } + ) + + +def _reindex_reordered_dims(arr: DataArray, expected: dict[Hashable, Any]) -> DataArray: + """ + Reindex shared dims whose values match ``expected`` in a different order. + + Disagreeing value *sets* are left for downstream xarray alignment; + only a pure reordering of the same values is conformed here. + """ for dim, coord_values in expected.items(): - if dim not in arr.dims: - continue - if isinstance(arr.indexes.get(dim), pd.MultiIndex): + if dim not in arr.dims or isinstance(arr.indexes.get(dim), pd.MultiIndex): continue expected_idx = _as_index(coord_values) actual_idx = arr.coords[dim].to_index() if actual_idx.equals(expected_idx): continue - # Same values, different order β†’ reindex to match expected order. - # Different value sets are left alone for downstream xarray alignment. if len(actual_idx) == len(expected_idx) and set(actual_idx) == set( expected_idx ): arr = arr.reindex({dim: expected_idx}) + return arr + - # expand_dims prepends new dimensions and their coordinate variables; - # the subsequent transpose restores coords order. Both are no-ops when - # the array already matches. Reconstruct so the DataArray's coords - # iteration order also follows coords (a Dataset built from this picks - # up its dim order from coord insertion). +def _expand_missing_dims(arr: DataArray, expected: dict[Hashable, Any]) -> DataArray: + """ + Broadcast ``arr`` over ``expected`` dims it does not yet carry. + + A MultiIndex-backed dim is broadcast against a proper ``Coordinates`` + template: plain ``expand_dims`` would drop its level coords and leave a + degenerate flat index that fails to align downstream. The exception is + when ``arr`` already carries one of the MultiIndex's level names β€” + broadcasting would then raise on the conflicting index, so fall back to + ``expand_dims``. + """ expand = {k: v for k, v in expected.items() if k not in arr.dims} - if expand: - # expand_dims drops the level coords of a MultiIndex-backed dim, - # leaving a degenerate flat index that fails to align downstream. - # Broadcast against a proper Coordinates template instead. - plain = {} - for dim, coord_values in expand.items(): - mi = _as_multiindex(coord_values) - # Fall back to expand_dims when arr already carries one of the - # MultiIndex's level names as its own coord: broadcasting against - # the level coords would raise on the conflicting index. - if mi is None or set(mi.names) & (set(arr.coords) | set(arr.dims)): - plain[dim] = coord_values - continue - template = DataArray( - np.zeros(len(mi)), - coords=Coordinates.from_pandas_multiindex(mi, dim), - dims=[dim], - ) - arr, _ = broadcast(arr, template) - if plain: - arr = arr.expand_dims(plain) + if not expand: + return arr + plain = {} + for dim, coord_values in expand.items(): + mi = _as_multiindex(coord_values) + if mi is None or set(mi.names) & (set(arr.coords) | set(arr.dims)): + plain[dim] = coord_values + continue + template = DataArray( + np.zeros(len(mi)), + coords=Coordinates.from_pandas_multiindex(mi, dim), + dims=[dim], + ) + arr, _ = broadcast(arr, template) + if plain: + arr = arr.expand_dims(plain) + return arr + +def _order_like_coords(arr: DataArray, expected: dict[Hashable, Any]) -> DataArray: + """ + Transpose ``arr`` to ``coords`` dim order, then match coord iteration order. + + The reconstruction makes a Dataset built from ``arr`` pick up its dim + order from coord insertion, not just the transpose. + """ target_dims = tuple(d for d in expected if d in arr.dims) + tuple( d for d in arr.dims if d not in expected ) arr = arr.transpose(*target_dims) - coord_order = [c for c in target_dims if c in arr.coords] + [ c for c in arr.coords if c not in target_dims ] @@ -651,7 +804,43 @@ def _broadcast_to_coords( coords={c: arr.coords[c] for c in coord_order}, name=arr.name, ) + return arr + +def _broadcast_to_coords( + arr: Any, + coords: CoordsLike | None = None, + dims: DimsLike | None = None, + **kwargs: Any, +) -> tuple[DataArray, list[_LevelProjection]]: + """ + Convert ``arr`` and broadcast it against ``coords`` (shared mechanics). + + Returns the broadcast DataArray together with the MultiIndex-level + projections performed along the way, so the public entry points can + apply their own policy (warn or raise) to partial projections and + coverage gaps. Unlabeled inputs pair their axes with the coords dims by + size (#736); see :func:`_label_input`. + """ + if coords is None: + return as_dataarray(arr, coords, dims, **kwargs), [] + + expected = _coords_to_dict(coords, dims=dims) + if not expected: + return as_dataarray(arr, coords, dims, **kwargs), [] + + if isinstance(arr, pd.Series | pd.DataFrame): + converted = _named_pandas_to_dataarray(arr) + if converted is not None: + arr = converted + + if not isinstance(arr, DataArray): + arr = _label_input(arr, expected, dims, **kwargs) + + arr, projections = _project_onto_multiindex_levels(arr, expected) + arr = _reindex_reordered_dims(arr, expected) + arr = _expand_missing_dims(arr, expected) + arr = _order_like_coords(arr, expected) return arr, projections @@ -709,9 +898,9 @@ def broadcast_to_coords( index names, e.g. ``period`` / ``timestep``) and *level combinations* (its elements β€” one tuple per position, e.g. ``(2030, 't1')``). Inputs indexed by levels instead of the dim itself are implicitly projected - onto the dim's level combinations. These projections are deprecated in - both modes and emit an :class:`~linopy.EvolvingAPIWarning`; the v1 - convention will require them to be explicit. Two cases: + onto the dim's level combinations. These projections are legacy-only: + under legacy semantics they emit a :class:`~linopy.LinopySemanticsWarning` + deprecation; under the v1 convention they raise. Two cases: - input misses a whole level β†’ broadcasts across it; warns in both modes. - input gives some level combinations no value (a *coverage gap*) β†’ @@ -743,7 +932,7 @@ def broadcast_to_coords( """ if not strict: da, projections = _broadcast_to_coords(arr, coords, dims, **kwargs) - _warn_implicit_projections(projections) + _enforce_implicit_projections(projections) return da if label is None: @@ -771,7 +960,7 @@ def broadcast_to_coords( f"{p.dim!r}: {preview}. The input is indexed by level(s) " f"{p.levels} and must cover every combination." ) - _warn_implicit_projections(projections) + _enforce_implicit_projections(projections) validate_alignment(da, coords, dims=dims, label=label) return da diff --git a/linopy/config.py b/linopy/config.py index 5d269c4e..13731e58 100644 --- a/linopy/config.py +++ b/linopy/config.py @@ -9,6 +9,28 @@ from typing import Any +LEGACY_SEMANTICS = "legacy" +V1_SEMANTICS = "v1" +VALID_SEMANTICS = {LEGACY_SEMANTICS, V1_SEMANTICS} + +LEGACY_SEMANTICS_MESSAGE = ( + "The 'legacy' semantics are deprecated and will be removed in " + "linopy 1.0. Set linopy.options['semantics'] = 'v1' to opt in " + "to the new behaviour, or silence this warning with:\n" + " import warnings; warnings.filterwarnings(" + "'ignore', category=LinopySemanticsWarning)" +) + + +class LinopySemanticsWarning(FutureWarning): + """ + Emitted when code runs under the legacy arithmetic semantics. + + Subclasses ``FutureWarning`` rather than ``DeprecationWarning`` so it is + shown to end users by default; the legacy-to-v1 transition changes + results, not just an API surface. + """ + class OptionSettings: """Runtime configuration knobs (e.g. display widths). Use as a context manager or set values directly via ``options(key=value)``.""" @@ -30,6 +52,11 @@ def set_value(self, **kwargs: Any) -> None: for k, v in kwargs.items(): if k not in self._defaults: raise KeyError(f"{k} is not a valid setting.") + if k == "semantics" and v not in VALID_SEMANTICS: + raise ValueError( + f"Invalid semantics: {v!r}. " + f"Must be one of {sorted(VALID_SEMANTICS)}." + ) self._current_values[k] = v def get_value(self, name: str) -> Any: @@ -62,4 +89,5 @@ def __repr__(self) -> str: options = OptionSettings( display_max_rows=14, display_max_terms=6, + semantics=LEGACY_SEMANTICS, ) diff --git a/linopy/examples.py b/linopy/examples.py index 6e1cfb15..549574bf 100644 --- a/linopy/examples.py +++ b/linopy/examples.py @@ -2,6 +2,7 @@ This module contains examples of linear programming models using the linopy library. """ +import xarray as xr from numpy import arange from linopy import Model @@ -73,7 +74,7 @@ def benchmark_model(n: int = 10, integerlabels: bool = False) -> Model: naxis, maxis = [arange(n, dtype=float), arange(n).astype(str)] x = m.add_variables(coords=[naxis, maxis]) y = m.add_variables(coords=[naxis, maxis]) - m.add_constraints(x - y >= naxis) + m.add_constraints(x - y >= xr.DataArray(naxis, dims=["dim_0"])) m.add_constraints(x + y >= 0) m.add_objective((2 * x).sum() + y.sum()) return m diff --git a/linopy/expressions.py b/linopy/expressions.py index ea8588d2..9205ef08 100644 --- a/linopy/expressions.py +++ b/linopy/expressions.py @@ -44,7 +44,13 @@ from types import EllipsisType, NotImplementedType from linopy import constraints, variables -from linopy.alignment import as_dataarray, broadcast_to_coords, fill_missing_coords +from linopy.alignment import ( + _matmul_operand_to_dataarray, + as_constant, + as_dataarray, + broadcast_to_coords, + fill_missing_coords, +) from linopy.common import ( EmptyDeprecationWrapper, LocIndexer, @@ -67,7 +73,9 @@ to_dataframe, to_polars, ) -from linopy.config import options +from linopy.config import ( + options, +) from linopy.constants import ( CV_DIM, EQUAL, @@ -80,6 +88,19 @@ STACKED_TERM_DIM, TERM_DIM, ) +from linopy.semantics import ( + _legacy_coord_mismatch_message, + _legacy_coord_reorder_message, + _legacy_nan_rhs_constraint_message, + _shared_dim_mismatch_message, + absorb_absence, + check_user_nan, + conform_merge_dims, + enforce_aux_conflict, + first_mismatched_dim, + is_v1, + warn_legacy, +) from linopy.types import ( CONSTANT_TYPES, ConstantLike, @@ -606,11 +627,18 @@ def _multiply_by_linear_expression( ds = other.data[["coeffs", "vars"]].sel(_term=0).broadcast_like(self.data) ds = assign_multiindex_safe(ds, const=other.const) res = merge([self, ds], dim=FACTOR_DIM, cls=QuadraticExpression) - # deal with cross terms c1 * v2 + c2 * v1 + # deal with cross terms c1 * v2 + c2 * v1. The ``const`` factors + # are internal Β§6-propagated fields β€” NaN at absent slots, not + # user-supplied data. ``fillna(0)`` makes them safe to pass back + # through the public-API ``*`` (which Β§5-checks user NaN); the + # zeroed cross-term contribution at an absent slot adds nothing, + # and ``res`` already carries the absence marker from the + # FACTOR_DIM merge above (NaN + 0 = NaN), so absence survives + # and ``absorb_absence`` enforces the storage invariant. if self.has_constant: - res = res + self.const * other.reset_const() + res = res + other.reset_const() * self.const.fillna(0) if other.has_constant: - res = res + self.reset_const() * other.const + res = res + self.reset_const() * other.const.fillna(0) return cast(QuadraticExpression, res) def _align_constant( @@ -629,7 +657,9 @@ def _align_constant( fill_value : float, default: 0 Fill value for missing coordinates. join : str, optional - Alignment method. If None, uses size-aware default behavior. + Alignment method. If None, the default is determined by + ``options["semantics"]`` β€” ``"exact"`` under ``v1``, the + legacy size-aware behavior under ``legacy``. Returns ------- @@ -638,37 +668,143 @@ def _align_constant( aligned : DataArray The aligned constant. needs_data_reindex : bool - Whether the expression's data needs reindexing. - """ + Whether the expression's data needs reindexing. ``False`` whenever the + shared-dim indexes already agree β€” the ``override`` / ``left`` paths + and ``exact`` join once ``first_mismatched_dim`` confirms the match β€” + so the caller skips a full-dataset ``reindex_like`` deepcopy on the + hot multiply/add path. Only the label-changing joins (inner / outer / + right) return ``True``. + """ + # Β§11: aux-coord conflict is independent of dim alignment β€” fires + # on every join path. Gating it behind ``join is None`` (alongside + # the Β§8 dim check) would leave ``join="override"`` etc. silently + # dropping the conflicting coord, which is the #295 bug v1 is + # meant to close. v1 raises; legacy warns. + enforce_aux_conflict([self.const, other], stacklevel=4) if join is None: - if other.sizes == self.const.sizes: - return self.const, other.assign_coords(coords=self.coords), False + if is_v1(): + join = "exact" + else: + # LEGACY: remove at 1.0 β€” see arithmetics-design/legacy-removal.md. + # Legacy default: positional when sizes match, else left-join. + mismatch = first_mismatched_dim(self.const, other) + if other.sizes == self.const.sizes: + if mismatch is not None: + warn_legacy( + _legacy_coord_mismatch_message( + "this operator's constant operand", *mismatch + ), + stacklevel=4, + ) + return self.const, other.assign_coords(coords=self.coords), False + warn_legacy( + _legacy_coord_mismatch_message( + "this operator's constant operand", + *(mismatch or (None, None, None)), + ), + stacklevel=4, + ) + return ( + self.const, + other.reindex_like(self.const, fill_value=fill_value), + False, + ) + + if join == "override": + # Β§10: ``override`` is the *explicit* form of the legacy + # positional alignment. Positional pairing requires the + # shared dims to have matching sizes on both sides β€” without + # that the labels we're about to ``assign_coords`` on would + # be a coord-rename of mismatched data, which is exactly the + # silent-wrong-answer hazard ``override`` was tightened up + # to prevent. Raise with a clear message instead of letting + # xarray broadcast / error opaquely downstream. + shared = set(self.const.dims) & set(other.dims) + bad = sorted( + (d for d in shared if self.const.sizes[d] != other.sizes[d]), + key=str, + ) + if bad: + sizes = ", ".join( + f"{d!r}: left={self.const.sizes[d]}, right={other.sizes[d]}" + for d in bad + ) + raise ValueError( + f"join='override' requires matching sizes on shared " + f"dimensions, but sizes differ on {sizes}. Use " + f"join='inner' / 'outer' / 'left' / 'right' to combine " + f"by label, or reshape one side first." + ) + return self.const, other.assign_coords(coords=self.coords), False + if join == "left": return ( self.const, other.reindex_like(self.const, fill_value=fill_value), False, ) - elif join == "override": - return self.const, other.assign_coords(coords=self.coords), False - else: - self_const, aligned = xr.align( - self.const, - other, - join=join, - fill_value=fill_value, - ) - return self_const, aligned, True + # ``xr.align(..., join="exact")`` raises with a wording that's not + # API-stable across xarray releases; matching on ``"exact" in str(e)`` + # would silently degrade if upstream rephrases. Do the Β§8 check + # ourselves and raise the canonical ``_shared_dim_mismatch_message`` + # (same text as the v1-default and merge paths). Other joins + # (inner / outer / right) handle coord mismatches via the join + # mode and don't error here. + if join == "exact": + mismatch = first_mismatched_dim(self.const, other) + if mismatch is not None: + raise ValueError(_shared_dim_mismatch_message(*mismatch)) + return self.const, other, False + self_const, aligned = xr.align( + self.const, + other, + join=join, + fill_value=fill_value, + ) + return self_const, aligned, True def _add_constant( self, other: ConstantLike, join: JoinOptions | None = None ) -> Self: - # NaN values in self.const or other are filled with 0 (additive identity) - # so that missing data does not silently propagate through arithmetic. + if is_v1(): + return self._add_constant_v1(other, join) + return self._add_constant_legacy(other, join) + + def _add_constant_v1(self, other: ConstantLike, join: JoinOptions | None) -> Self: + # Β§6: absence propagates β€” self.const NaN stays NaN, no fillna(0). + # Β§5: user NaN raised in check_user_nan; never reaches the math here. + if np.isscalar(other) and join is None: + if isinstance(other, float) and np.isnan(other): + check_user_nan() + return self.assign(const=self.const + other) + da = broadcast_to_coords(other, coords=self.coords, strict=False) + if da.isnull().any(): + check_user_nan() + self_const, da, needs_data_reindex = self._align_constant( + da, fill_value=0, join=join + ) + if needs_data_reindex: + return self.__class__( + self.data.reindex_like(self_const, fill_value=self._fill_value).assign( + const=self_const + da + ), + self.model, + ) + return self.assign(const=self_const + da) + + # LEGACY: remove at 1.0 β€” see arithmetics-design/legacy-removal.md. + def _add_constant_legacy( + self, other: ConstantLike, join: JoinOptions | None + ) -> Self: + # NaN values in self.const or other are silently filled with 0 + # (additive identity) so missing data does not propagate through + # arithmetic. ``check_user_nan`` only warns under legacy. if np.isscalar(other) and join is None: + if isinstance(other, float) and np.isnan(other): + check_user_nan() return self.assign(const=self.const.fillna(0) + other) - da = broadcast_to_coords( - other, coords=self.coords, dims=self.coord_dims, strict=False - ) + da = broadcast_to_coords(other, coords=self.coords, strict=False) + if da.isnull().any(): + check_user_nan() self_const, da, needs_data_reindex = self._align_constant( da, fill_value=0, join=join ) @@ -688,18 +824,60 @@ def _apply_constant_op( other: ConstantLike, op: Callable[[DataArray, DataArray], DataArray], fill_value: float, + op_kind: str, join: JoinOptions | None = None, ) -> Self: - """ - Apply a constant operation (mul, div, etc.) to this expression with a scalar or array. + """Apply a constant operation (mul, div) to this expression.""" + if is_v1(): + return self._apply_constant_op_v1(other, op, fill_value, op_kind, join) + return self._apply_constant_op_legacy(other, op, fill_value, op_kind, join) - NaN values are filled with neutral elements before the operation: - - factor (other) is filled with fill_value (0 for mul, 1 for div) - - coeffs and const are filled with 0 (additive identity) - """ - factor = broadcast_to_coords( - other, coords=self.coords, dims=self.coord_dims, strict=False + def _apply_constant_op_v1( + self, + other: ConstantLike, + op: Callable[[DataArray, DataArray], DataArray], + fill_value: float, + op_kind: str, + join: JoinOptions | None, + ) -> Self: + # Β§6: NaN in coeffs/const propagates through op (NaN * x = NaN). + # Β§5: user NaN raised before we get here. + if isinstance(other, float) and np.isnan(other): + check_user_nan(op_kind=op_kind) + factor = broadcast_to_coords(other, coords=self.coords, strict=False) + if factor.isnull().any(): + check_user_nan(op_kind=op_kind) + self_const, factor, needs_data_reindex = self._align_constant( + factor, fill_value=fill_value, join=join ) + if needs_data_reindex: + data = self.data.reindex_like(self_const, fill_value=self._fill_value) + return self.__class__( + assign_multiindex_safe( + data, + coeffs=op(data.coeffs, factor), + const=op(self_const, factor), + ), + self.model, + ) + return self.assign(coeffs=op(self.coeffs, factor), const=op(self_const, factor)) + + # LEGACY: remove at 1.0 β€” see arithmetics-design/legacy-removal.md. + def _apply_constant_op_legacy( + self, + other: ConstantLike, + op: Callable[[DataArray, DataArray], DataArray], + fill_value: float, + op_kind: str, + join: JoinOptions | None, + ) -> Self: + # NaN values are silently filled with neutral elements before the op: + # factor β†’ fill_value (0 for mul, 1 for div), coeffs/const β†’ 0. + if isinstance(other, float) and np.isnan(other): + check_user_nan(op_kind=op_kind) + factor = broadcast_to_coords(other, coords=self.coords, strict=False) + if factor.isnull().any(): + check_user_nan(op_kind=op_kind) self_const, factor, needs_data_reindex = self._align_constant( factor, fill_value=fill_value, join=join ) @@ -720,14 +898,19 @@ def _apply_constant_op( def _multiply_by_constant( self, other: ConstantLike, join: JoinOptions | None = None ) -> Self: - return self._apply_constant_op(other, operator.mul, fill_value=0, join=join) + return self._apply_constant_op( + other, operator.mul, fill_value=0, op_kind="mul", join=join + ) def _divide_by_constant( self, other: ConstantLike, join: JoinOptions | None = None ) -> Self: - return self._apply_constant_op(other, operator.truediv, fill_value=1, join=join) + return self._apply_constant_op( + other, operator.truediv, fill_value=1, op_kind="div", join=join + ) def __div__(self, other: SideLike) -> Self: + other = as_constant(other) try: if isinstance(other, SUPPORTED_EXPRESSION_TYPES): raise TypeError( @@ -1182,15 +1365,46 @@ def to_constraint( which are moved to the left-hand-side and constant values which are moved to the right-hand side. """ + rhs = as_constant(rhs) if self.is_constant and is_constant(rhs): raise ValueError( f"Both sides of the constraint are constant. At least one side must contain variables. {self} {rhs}" ) - if isinstance(rhs, CONSTANT_TYPES): - rhs = broadcast_to_coords( - rhs, coords=self.coords, dims=self.coord_dims, strict=False + if is_v1(): + # Β§5 + Β§12: the RHS is a user-supplied constant just like any + # operand in arithmetic. Validate NaN here (raise) and let + # ``self.sub(rhs)`` do the Β§8 alignment β€” no silent + # reindex_like-pad that would mask a coordinate mismatch. + # An absent slot in ``self.const`` (propagated from Β§6) flows + # through ``sub`` into the RHS and reaches downstream + # auto-mask handling as "no constraint at this row" (Β§12). + if isinstance(rhs, CONSTANT_TYPES): + rhs = broadcast_to_coords(rhs, coords=self.coords, strict=False) + extra_dims = set(rhs.dims) - set(self.coord_dims) + if extra_dims: + logger.warning( + f"Constant RHS contains dimensions {extra_dims} not present " + f"in the expression, which might lead to inefficiencies. " + f"Consider collapsing the dimensions by taking min/max." + ) + if rhs.isnull().any(): + check_user_nan() + all_to_lhs = self.sub(rhs, join=join).data + computed_rhs = -all_to_lhs.const + data = assign_multiindex_safe( + all_to_lhs[["coeffs", "vars"]], sign=sign, rhs=computed_rhs ) + return constraints.Constraint(data, model=self.model) + + # LEGACY: remove at 1.0 β€” see arithmetics-design/legacy-removal.md. + # Legacy auto-mask path: NaN RHS is silently preserved as "no + # constraint at this row" (the legacy reindex_like-pad fills + # subset coords with NaN, then `sub` would fill them with 0 as + # part of normal arithmetic, so we restore the original NaN mask + # afterward). + if isinstance(rhs, CONSTANT_TYPES): + rhs = broadcast_to_coords(rhs, coords=self.coords, strict=False) extra_dims = set(rhs.dims) - set(self.coord_dims) if extra_dims: @@ -1199,11 +1413,17 @@ def to_constraint( f"in the expression, which might lead to inefficiencies. " f"Consider collapsing the dimensions by taking min/max." ) + # Two legacy paths can introduce NaN into ``rhs`` and the + # ``rhs_nan_mask`` below treats both the same (auto-mask drops + # the row). The fix-hints differ though, so check each cause + # before the reindex obscures them, and warn the right one. + mismatch = first_mismatched_dim(self.const, rhs) + if mismatch is not None: + warn_legacy(_legacy_coord_mismatch_message("constraint RHS", *mismatch)) + if bool(rhs.isnull().any()): + warn_legacy(_legacy_nan_rhs_constraint_message()) rhs = rhs.reindex_like(self.const, fill_value=np.nan) - # Remember where RHS is NaN (meaning "no constraint") before the - # subtraction, which may fill NaN with 0 as part of normal - # expression arithmetic. if isinstance(rhs, DataArray): rhs_nan_mask = rhs.isnull() else: @@ -1212,8 +1432,6 @@ def to_constraint( all_to_lhs = self.sub(rhs, join=join).data computed_rhs = -all_to_lhs.const - # Restore NaN at positions where the original constant RHS had no - # value so that downstream code still treats them as unconstrained. if rhs_nan_mask is not None and rhs_nan_mask.any(): computed_rhs = xr.where(rhs_nan_mask, np.nan, computed_rhs) @@ -1230,12 +1448,21 @@ def reset_const(self) -> Self: def isnull(self) -> DataArray: """ - Get a boolean mask with true values where there is only missing values in an expression. - - Returns - ------- - xr.DataArray - """ + Get a boolean mask reporting which slots are absent. + + Under v1 (Β§3), a slot is absent iff ``const`` is NaN β€” every + operator propagates const NaN whenever an operand is absent, so + this is the universal signal. A term with ``vars == -1`` is a + *dead term* (contributes nothing), not a slot-level absence + signal; ``fillna(value)`` can revive an absent slot to a present + constant while leaving the dead sentinel term in place. The + legacy convention has no real absence concept (``const`` is + always filled with 0); the historical AND of "all vars + sentinel" and "const NaN" is preserved verbatim under legacy. + """ + if is_v1(): + return self.const.isnull() + # LEGACY: remove at 1.0 β€” see arithmetics-design/legacy-removal.md. helper_dims = set(self.vars.dims).intersection(HELPER_DIMS) return (self.vars == -1).all(helper_dims) & self.const.isnull() @@ -1699,6 +1926,7 @@ def __add__( Note: If other is a numpy array or pandas object without axes names, dimension names of self will be filled in other """ + other = as_constant(other) if isinstance(other, QuadraticExpression): return other.__add__(self) @@ -1734,6 +1962,7 @@ def __sub__( | LinearExpression | QuadraticExpression, ) -> LinearExpression | QuadraticExpression: + other = as_constant(other) try: return self.__add__(-other) except TypeError: @@ -1758,6 +1987,7 @@ def __mul__( """ Multiply the expr by a factor. """ + other = as_constant(other) if isinstance(other, QuadraticExpression): return other.__rmul__(self) @@ -1803,8 +2033,9 @@ def __matmul__( """ Matrix multiplication with other, similar to xarray dot. """ + other = as_constant(other) if not isinstance(other, LinearExpression | variables.Variable): - other = as_dataarray(other, coords=self.coords, dims=self.coord_dims) + other = _matmul_operand_to_dataarray(other, self.coords, self.coord_dims) common_dims = list(set(self.coord_dims).intersection(other.dims)) return (self * other).sum(dim=common_dims) @@ -2215,6 +2446,7 @@ def __mul__(self, other: SideLike) -> QuadraticExpression: """ Multiply the expr by a factor. """ + other = as_constant(other) if isinstance(other, SUPPORTED_EXPRESSION_TYPES): raise TypeError( "unsupported operand type(s) for *: " @@ -2236,6 +2468,7 @@ def __add__(self, other: SideLike) -> QuadraticExpression: Note: If other is a numpy array or pandas object without axes names, dimension names of self will be filled in other """ + other = as_constant(other) try: if isinstance(other, CONSTANT_TYPES): return self._add_constant(other) @@ -2262,6 +2495,7 @@ def __sub__(self, other: SideLike) -> QuadraticExpression: Note: If other is a numpy array or pandas object without axes names, dimension names of self will be filled in other """ + other = as_constant(other) try: return self.__add__(-other) except TypeError: @@ -2285,12 +2519,13 @@ def __matmul__( """ Matrix multiplication with other, similar to xarray dot. """ + other = as_constant(other) if isinstance(other, SUPPORTED_EXPRESSION_TYPES): raise TypeError( "Higher order non-linear expressions are not yet supported." ) - other = as_dataarray(other, coords=self.coords, dims=self.coord_dims) + other = _matmul_operand_to_dataarray(other, self.coords, self.coord_dims) common_dims = list(set(self.coord_dims).intersection(other.dims)) return (self * other).sum(dim=common_dims) @@ -2427,7 +2662,7 @@ def as_expression( try: obj = broadcast_to_coords(obj, strict=False, **kwargs) except ValueError as e: - raise ValueError("Cannot convert to LinearExpression") from e + raise ValueError(f"Cannot convert to LinearExpression: {e}") from e return LinearExpression(obj, model) @@ -2523,6 +2758,28 @@ def merge( model = exprs[0].model + data = [e.data if isinstance(e, linopy_types) else e for e in exprs] + data = [fill_missing_coords(ds, fill_helper_dims=True) for ds in data] + + # Β§8: v1 aligns a reorder by label and raises a dim mismatch (before the + # Β§11 aux check, else a MultiIndex mismatch reads as a level-coord + # conflict). Legacy keeps positional alignment and only warns. + if join is None: + data, mismatch, reorder = conform_merge_dims(data, concat_dim=dim) + if is_v1(): + if mismatch is not None: + raise ValueError(_shared_dim_mismatch_message(*mismatch)) + elif mismatch is not None: # LEGACY: remove at 1.0 + warn_legacy( + _legacy_coord_mismatch_message(f"merge along dim {dim!r}", *mismatch) + ) + elif reorder is not None: # LEGACY: remove at 1.0 + warn_legacy( + _legacy_coord_reorder_message(f"merge along dim {dim!r}", *reorder) + ) + + enforce_aux_conflict(data) # Β§11 + if join is not None: override = join == "override" elif issubclass(cls, linopy_types) and dim in HELPER_DIMS: @@ -2533,9 +2790,6 @@ def merge( else: override = False - data = [e.data if isinstance(e, linopy_types) else e for e in exprs] - data = [fill_missing_coords(ds, fill_helper_dims=True) for ds in data] - if not kwargs: kwargs = { "coords": "minimal", @@ -2556,12 +2810,26 @@ def merge( if dim == TERM_DIM: ds = xr.concat([d[["coeffs", "vars"]] for d in data], dim, **kwargs) subkwargs = {**kwargs, "fill_value": 0} - const = xr.concat([d["const"] for d in data], dim, **subkwargs).sum(TERM_DIM) + # Under v1, Β§6 requires that an absent slot in any operand stays + # absent in the result; ``sum(skipna=False)`` propagates NaN + # rather than collapsing it to 0. + skipna = not is_v1() + const = xr.concat([d["const"] for d in data], dim, **subkwargs).sum( + TERM_DIM, skipna=skipna + ) ds = assign_multiindex_safe(ds, const=const) elif dim == FACTOR_DIM: ds = xr.concat([d[["vars"]] for d in data], dim, **kwargs) - coeffs = xr.concat([d["coeffs"] for d in data], dim, **kwargs).prod(FACTOR_DIM) - const = xr.concat([d["const"] for d in data], dim, **kwargs).prod(FACTOR_DIM) + # Β§6 also applies to the quadratic build: an absent factor must + # stay absent (``prod(skipna=False)`` β†’ NaN) rather than collapse + # to multiplicative identity 1. Matches the TERM_DIM branch above. + skipna = not is_v1() + coeffs = xr.concat([d["coeffs"] for d in data], dim, **kwargs).prod( + FACTOR_DIM, skipna=skipna + ) + const = xr.concat([d["const"] for d in data], dim, **kwargs).prod( + FACTOR_DIM, skipna=skipna + ) ds = assign_multiindex_safe(ds, coeffs=coeffs, const=const) else: ds = xr.concat(data, dim, **kwargs) @@ -2569,6 +2837,9 @@ def merge( for d in set(HELPER_DIMS) & set(ds.coords): ds = ds.reset_index(d, drop=True) + if is_v1(): + ds = absorb_absence(ds) + return cls(ds, model) diff --git a/linopy/piecewise.py b/linopy/piecewise.py index 25a0ce17..2f5afde9 100644 --- a/linopy/piecewise.py +++ b/linopy/piecewise.py @@ -1714,12 +1714,22 @@ def _add_incremental( ) if n_pieces >= 2: + # ``piece_dim`` is a *positional* relation here: the constraint pairs + # each lower piece with the next-higher one. The two ``isel`` slices + # share the dim but carry different labels (first n-1 vs last n-1 of + # piece_index), which v1 Β§8 rejects. Relabel the high slice onto the + # low slice's labels so alignment-by-label gives the intended + # positional pairing β€” convention Β§10's explicit-positional path. delta_lo = delta_var.isel({piece_dim: slice(None, -1)}, drop=True) - delta_hi = delta_var.isel({piece_dim: slice(1, None)}, drop=True) + delta_hi = delta_var.isel({piece_dim: slice(1, None)}, drop=True).assign_coords( + {piece_dim: delta_lo.coords[piece_dim]} + ) model.add_constraints( delta_hi <= delta_lo, name=f"{name}{PWL_FILL_ORDER_SUFFIX}" ) - binary_hi = binary_var.isel({piece_dim: slice(1, None)}, drop=True) + binary_hi = binary_var.isel( + {piece_dim: slice(1, None)}, drop=True + ).assign_coords({piece_dim: delta_lo.coords[piece_dim]}) model.add_constraints( binary_hi <= delta_lo, name=f"{name}{PWL_BINARY_ORDER_SUFFIX}" ) @@ -1727,7 +1737,10 @@ def _add_incremental( def _incremental_weighted(bp: DataArray) -> LinearExpression: steps = bp.diff(dim).rename({dim: piece_dim}) steps[piece_dim] = piece_index - bp0 = bp.isel({dim: 0}) + # ``drop=True`` keeps the breakpoint coord from sticking around as a + # scalar on ``bp0_term`` β€” otherwise Β§11 rejects it as an aux-coord + # conflict against the constraint LHS. + bp0 = bp.isel({dim: 0}, drop=True) bp0_term: DataArray | LinearExpression = bp0 if active is not None: bp0_term = bp0 * active diff --git a/linopy/semantics.py b/linopy/semantics.py new file mode 100644 index 00000000..38e06695 --- /dev/null +++ b/linopy/semantics.py @@ -0,0 +1,445 @@ +""" +v1 semantics helpers. + +Single home for the predicates, validators, and storage-invariant +enforcement that the v1 arithmetic convention requires. Importing from +here keeps ``expressions.py`` focused on the operator dispatch and lets a +future legacy removal be a single-file delete. + +See ``arithmetics-design/convention.md`` for the rules these helpers +implement and ``arithmetics-design/goals.md`` for the design intent. +""" + +from __future__ import annotations + +import os +import sys +from collections.abc import Sequence +from typing import Any +from warnings import warn + +import numpy as np +import pandas as pd +from xarray import DataArray, Dataset + +from linopy.config import ( + V1_SEMANTICS, + LinopySemanticsWarning, + options, +) +from linopy.constants import HELPER_DIMS + + +def _user_nan_message() -> str: + """User-NaN error text β€” distinguishes the two intents a user might have.""" + return ( + "NaN found in a user-supplied constant. linopy treats this as " + "ambiguous: if you meant a *data error*, fix it with .fillna(value); " + "if you meant *absent at this slot*, mark it on the variable " + "instead (mask=, .where(cond), .reindex(...), .shift(...))." + ) + + +def _shared_dim_mismatch_message(dim: str, left: Any, right: Any) -> str: + """Shared-dim error text β€” names the dim and shows the disagreeing labels.""" + return ( + f"Coordinate mismatch on shared dimension {dim!r}: " + f"left={_short_repr(left)}, right={_short_repr(right)}. " + "Resolve with `.sel(...)` / `.reindex(...)` to align before " + "combining, with `.assign_coords(...)` to relabel one side " + "(positional alignment, made explicit), with `linopy.align(...)` " + "to pre-align several operands at once, or by passing an explicit " + "`join=` argument to `.add` / `.sub` / `.mul` / `.div` / `.le` / " + "`.ge` / `.eq` (accepts inner / outer / left / right / override)." + ) + + +def _aux_conflict_message(name: str, left: Any, right: Any, kind: str) -> str: + """ + Aux-coord error text β€” names the coord, the failure mode (shape vs + value), and shows the disagreeing values. + """ + if kind == "shape": + problem = ( + f"Auxiliary coordinate {name!r} has differing shapes across " + f"operands: left.shape={np.shape(left)}, " + f"right.shape={np.shape(right)}. " + ) + else: + problem = ( + f"Auxiliary coordinate {name!r} has conflicting values across " + f"operands: left={_short_repr(left)}, right={_short_repr(right)}. " + ) + return ( + problem + "xarray would silently drop the conflict; linopy raises so the " + f"caller resolves it. Use `.drop_vars({name!r})` to remove the " + f"coord, `.assign_coords({name}=...)` to relabel one side, or " + "`.isel(..., drop=True)` if the coord was introduced by a " + "scalar isel." + ) + + +# --------------------------------------------------------------------------- +# Legacy-deprecation warnings β€” actionable, per-site (goal #2 in goals.md: +# tell the user *what* will change for the op they just ran). +# Each helper returns the message; ``warn_legacy(msg)`` issues it. +# --------------------------------------------------------------------------- + + +_OPT_IN_HINT = ( + "\n Opt in: linopy.options['semantics'] = 'v1'" + "\n Silence: warnings.filterwarnings('ignore', " + "category=LinopySemanticsWarning)" +) + + +# Per-op opening clause for ``_legacy_nan_constant_message`` β€” operand +# noun and the historical fill value (`+`/`*` filled with 0; `/` filled +# with 1, a different fill that's worth calling out at the warn site). +_LEGACY_NAN_FILL_CLAUSE = { + "add": ( + "NaN in the constant operand was silently treated as 0 by legacy" + " (additive identity)." + ), + "mul": ( + "NaN in the multiplicative factor was silently treated as 0 by" + " legacy (so the variable was zeroed out at that slot)." + ), + "div": ( + "NaN in the divisor was silently treated as 1 by legacy (a" + " different fill from `+`/`*` which use 0)." + ), +} + + +def _legacy_nan_constant_message(op_kind: str) -> str: + """Legacy NaN-fill warning for `+`/`*`/`/`, keyed by ``op_kind``.""" + return ( + _LEGACY_NAN_FILL_CLAUSE[op_kind] + " Under v1 this raises ValueError." + "\n Resolve: `.fillna(value)` (data error)" + "\n or `mask=` / `.where(cond)` / `.reindex(...)` " + "on the variable (intended absence)." + _OPT_IN_HINT + ) + + +def _legacy_coord_mismatch_message( + context: str, + dim: str | None = None, + left: Any = None, + right: Any = None, +) -> str: + """ + Mismatched dim coords silently aligned (positional or left-join). + + When ``dim`` / ``left`` / ``right`` are given, the message names the + offending dim and shows the diff β€” same shape as the v1-raise text + so the user sees the same information at warn time as at raise time. + """ + diff = ( + f"\n Dim: {dim!r}: left={_short_repr(left)}, right={_short_repr(right)}" + if dim is not None + else "" + ) + return ( + f"Coordinate mismatch in {context} silently aligned by legacy" + " (positional when sizes match, otherwise left-join)." + " Under v1 this raises ValueError." + + diff + + "\n Resolve: `.sel(...)` / `.reindex(...)` to align" + "\n `.assign_coords(...)` to relabel one side" + "\n `linopy.align(...)` to pre-align several operands" + "\n or pass an explicit `join=` argument." + _OPT_IN_HINT + ) + + +def _legacy_coord_reorder_message(context: str, dim: str, left: Any, right: Any) -> str: + """Same labels, different order β€” aligned positionally by legacy; v1 reindexes.""" + return ( + f"Coordinate order mismatch in {context} aligned positionally by legacy." + " Under v1 the same labels in a different order align by label (a" + " reindex), giving a different result." + f"\n Dim: {dim!r}: left={_short_repr(left)}, right={_short_repr(right)}" + "\n Resolve: `.sel(...)` / `.reindex(...)` to align" + "\n `.assign_coords(...)` to relabel one side" + "\n or pass an explicit `join=` argument." + _OPT_IN_HINT + ) + + +def _legacy_aux_conflict_message(name: str, left: Any, right: Any, kind: str) -> str: + """ + Conflicting aux coord silently dropped by xarray under legacy. + + The diff line names the failure mode (shape vs value) β€” same shape + as the v1-raise text so the user sees the same information at warn + time as at raise time. + """ + if kind == "shape": + diff = f"\n Shapes: left={np.shape(left)}, right={np.shape(right)}" + else: + diff = f"\n Values: left={_short_repr(left)}, right={_short_repr(right)}" + return ( + f"Auxiliary coordinate {name!r} was conflicting across operands" + " and silently dropped by legacy (xarray's default)." + " Under v1 this raises ValueError." + + diff + + f"\n Resolve: `.drop_vars({name!r})`" + f"\n `.assign_coords({name}=...)` to relabel one side" + "\n or `.isel(..., drop=True)` if a scalar isel " + "introduced it." + _OPT_IN_HINT + ) + + +def _legacy_nan_rhs_constraint_message() -> str: + """Constraint RHS NaN silently kept as 'no constraint at this row'.""" + return ( + "NaN in the constraint RHS was silently kept as 'no constraint" + " at this row' by legacy auto-mask. Under v1 this raises" + " ValueError." + "\n Resolve: `mask=` on the variable for explicit per-row " + "masking" + "\n or `.fillna(value)` if the NaN was a data error." + _OPT_IN_HINT + ) + + +def _legacy_masked_variable_message(name: str) -> str: + """A masked/shifted/reindexed variable used in arithmetic under legacy.""" + return ( + f"Variable {name!r} has absent slots (from `mask=` / `.where()`" + " / `.shift()` / `.reindex()`). Under legacy each absent slot" + " contributes 0 to the resulting expression's terms (so `x + y" + " >= 10` reduces to `x >= 10` there). Under v1 the absence" + " propagates through arithmetic instead (`x + y` becomes absent" + " at the slot and the constraint drops)." + f"\n Resolve: wrap with `{name}.fillna(0)` for the legacy" + " behaviour under v1" + "\n (no fix needed if you only use the variable in a" + " constraint LHS alone β€” `y >= 0` drops the same way in both)." + _OPT_IN_HINT + ) + + +_LINOPY_ROOT = os.path.dirname(os.path.abspath(__file__)) + + +def warn_legacy(message: str, *, stacklevel: int | None = None) -> None: + """ + Emit a `LinopySemanticsWarning` whose source-frame points at the + first call-stack frame *outside* the linopy package. + + Static ``stacklevel`` doesn't fit here β€” the call-chain depth from + ``warn_legacy`` to the user's code varies per site (e.g. masked-var + via ``__add__`` is 5 frames deep, via ``Variable.fillna`` is 4). On + Python 3.12+ we use the stdlib ``skip_file_prefixes`` argument + (implemented and tested in CPython); on 3.11 we fall back to a + static ``stacklevel=5``, good enough for the common merge chain. + Pass an explicit ``stacklevel`` to override (e.g. for tests). + """ + if stacklevel is not None: + warn(message, LinopySemanticsWarning, stacklevel=stacklevel) + elif sys.version_info >= (3, 12): + warn( + message, + LinopySemanticsWarning, + skip_file_prefixes=(_LINOPY_ROOT,), + ) + else: + warn(message, LinopySemanticsWarning, stacklevel=5) + + +def _short_repr(values: Any, limit: int = 6) -> str: + """Render an array-like as a short, readable string for error messages.""" + arr = np.asarray(values) + if arr.ndim == 0: + return repr(arr.item()) + flat = arr.ravel() + if flat.size <= limit: + return repr(flat.tolist()) + head = ", ".join(repr(v) for v in flat[:limit].tolist()) + return f"[{head}, ... ({flat.size} total)]" + + +def is_v1() -> bool: + """True iff the current semantics is v1.""" + return options["semantics"] == V1_SEMANTICS + + +def check_user_nan(*, op_kind: str = "add") -> None: + """ + Enforce Β§5 for a user-supplied constant (scalar or array). + + v1 raises ``ValueError`` with the generic user-NaN message; legacy + warns with operator-specific text (``"add"`` covers +/-, ``"mul"``, + ``"div"`` β€” they differ in which fill value legacy applied). + """ + if is_v1(): + raise ValueError(_user_nan_message()) + warn_legacy(_legacy_nan_constant_message(op_kind), stacklevel=5) + + +def enforce_aux_conflict(datasets: Sequence[Any], *, stacklevel: int = 5) -> None: + """ + Enforce Β§11 across the given operands: v1 raises on aux-coord + conflict, legacy warns (xarray would silently drop it). + """ + conflict = conflicting_aux_coord(datasets) + if conflict is None: + return + if is_v1(): + raise ValueError(_aux_conflict_message(*conflict)) + warn_legacy(_legacy_aux_conflict_message(*conflict), stacklevel=stacklevel) + + +def dim_coords_differ(a: DataArray, b: DataArray) -> bool: + """True if a and b share a dimension whose coordinate labels disagree.""" + return first_mismatched_dim(a, b) is not None + + +def first_mismatched_dim(a: DataArray, b: DataArray) -> tuple[str, Any, Any] | None: + """ + Return ``(dim, a_labels, b_labels)`` for the first shared dim that + disagrees on coordinate labels OR size, or ``None`` if all agree. + + Uses ``indexes[dim]`` (the bare pandas Index) rather than + ``coords[dim]`` β€” a coord DataArray's ``equals`` compares attached + aux coords too, which gives a false positive when only one operand + carries an aux coord on the shared dim (Β§11's territory, not Β§8's). + """ + for dim in set(a.dims) & set(b.dims): + if dim in a.indexes and dim in b.indexes: + if not a.indexes[dim].equals(b.indexes[dim]): + return str(dim), a.indexes[dim].values, b.indexes[dim].values + elif a.sizes[dim] != b.sizes[dim]: + return str(dim), None, None + return None + + +def conform_merge_dims( + datasets: Sequence[Dataset], concat_dim: str +) -> tuple[list[Dataset], tuple[str, Any, Any] | None, tuple[str, Any, Any] | None]: + """ + Inspect shared user dims for a merge, in a single pass over the operands. + + Returns ``(data, mismatch, reorder)``. A shared user dim whose labels are + the first operand's in a different order (same set, including a stacked + MultiIndex's tuples) is a *reorder*; one whose label set differs is a + *mismatch* β€” each reported as ``(dim, first_labels, other_labels)`` (first + found). Under v1, reorders are aligned to the first operand's order in the + returned data (via positional ``isel`` β€” ``reindex`` cannot reorder a + MultiIndex by tuple) and ``reorder`` is ``None``; the caller raises on + ``mismatch``. Under legacy, nothing is aligned and the caller warns: + ``reorder`` (v1 would align by label) or ``mismatch`` (v1 would raise). + Helper dims (``_term``, ``_factor``) and the concat dim are excluded; bare + dimension indexes are compared, so auxiliary coords stay Β§11's job. + """ + datasets = list(datasets) + if len(datasets) < 2: + return datasets, None, None + skip = set(HELPER_DIMS) | {concat_dim} + indexed = [ + {k: d.indexes[k] for k in d.dims if k not in skip and k in d.indexes} + for d in datasets + ] + shared = set.intersection(*(set(p) for p in indexed)) + if not shared: + return datasets, None, None + + permute = is_v1() + out = [datasets[0]] + mismatch: tuple[str, Any, Any] | None = None + reorder: tuple[str, Any, Any] | None = None + for i in range(1, len(datasets)): + plan: dict[Any, Any] = {} + for d in shared: + ref, idx = indexed[0][d], indexed[i][d] + if ref.equals(idx): + continue + positions = idx.get_indexer(ref) if len(idx) == len(ref) else None + if positions is not None and (positions >= 0).all(): + if permute: + plan[d] = positions + elif reorder is None: + reorder = (str(d), ref.values, idx.values) + elif mismatch is None: + mismatch = (str(d), ref.values, idx.values) + out.append(datasets[i].isel(plan) if plan else datasets[i]) + return out, mismatch, reorder + + +def conflicting_aux_coord( + datasets: Sequence[Any], +) -> tuple[str, Any, Any, str] | None: + """ + Find an auxiliary (non-dim) coord that two or more operands carry with + disagreeing values. + + Returns ``(name, left_values, right_values, kind)`` for the first + conflict found, or ``None`` if every shared aux coord agrees. ``kind`` + is ``"shape"`` if the two operands carry differently-shaped values for + the coord (e.g. one is a vector, the other a scalar), or ``"value"`` + if shapes agree but the values themselves disagree. The two failure + modes get different error text downstream. + + Per Β§11, an aux coord either propagates (values agree across operands) + or surfaces as an error; xarray's default silently drops the conflict + and is what this check intercepts under v1. When only one operand + carries the coord (``len(present) < 2``), it propagates from that + operand unchanged. + """ + if not datasets: + return None + all_names: set[str] = set() + for d in datasets: + all_names.update(d.coords) + for name in all_names: + present = [ + d.coords[name].values + for d in datasets + if name in d.coords and name not in d.dims + ] + # Β§11 asymmetric-presence: when only one operand carries the coord, + # it propagates unchanged β€” no conflict to surface. + if len(present) < 2: + continue + ref = present[0] + for vals in present[1:]: + if ref.shape != vals.shape: + return str(name), ref, vals, "shape" + if not _aux_values_equal(ref, vals): + return str(name), ref, vals, "value" + return None + + +def _aux_values_equal(a: np.ndarray, b: np.ndarray) -> bool: + """ + Equality for aux-coord value arrays with NaN-equal-NaN semantics + on every dtype. + + ``np.array_equal(..., equal_nan=True)`` only works on float dtypes + (it calls ``isnan`` which crashes on object/string). Aux coords on + object dtype can embed ``np.nan`` placeholders (e.g. ragged category + labels), and we want two operands with identical NaN placement to + compare equal β€” pandas' element-equality already treats NaN as + self-equal for object dtypes, so route through ``pd.Series.equals``. + """ + if np.issubdtype(a.dtype, np.floating): + return bool(np.array_equal(a, b, equal_nan=True)) + return bool(pd.Series(a.ravel()).equals(pd.Series(b.ravel()))) + + +def absorb_absence(ds: Dataset) -> Dataset: + """ + Enforce the v1 dead-term invariant on a merged dataset. + + ``const.isnull()`` at a slot β‡’ every term at that slot must have + ``coeffs = NaN`` and ``vars = -1``. After ``merge`` concatenates two + expressions along ``_term``, a slot that's absent in one operand + still carries the *other* operand's valid term in its row; this + helper masks those away so the Β§1/Β§2 storage invariant holds. + """ + if "const" not in ds or "coeffs" not in ds or "vars" not in ds: + return ds + mask = ds["const"].isnull() + if not bool(mask.any()): + return ds + coeffs = ds["coeffs"].where(~mask, np.nan) + vars_ = ds["vars"].where(~mask, -1) + return ds.assign(coeffs=coeffs, vars=vars_) diff --git a/linopy/sos_reformulation.py b/linopy/sos_reformulation.py index 4abfb755..eb2b982d 100644 --- a/linopy/sos_reformulation.py +++ b/linopy/sos_reformulation.py @@ -184,8 +184,13 @@ def reformulate_sos2( added_constraints = [first_name] + # Scalar isel keeps ``sos_dim`` as a leftover non-dim coord whose value + # differs between ``x``/``M`` (indexed at ``n-1``) and ``z`` (indexed at + # ``n-2``). v1 Β§11 rejects that aux-coord conflict, so we ``drop=True`` + # to remove ``sos_dim`` from the comparison entirely. model.add_constraints( - x_expr.isel({sos_dim: 0}) <= M.isel({sos_dim: 0}) * z_expr.isel({sos_dim: 0}), + x_expr.isel({sos_dim: 0}, drop=True) + <= M.isel({sos_dim: 0}, drop=True) * z_expr.isel({sos_dim: 0}, drop=True), name=first_name, ) @@ -211,8 +216,9 @@ def reformulate_sos2( added_constraints.append(mid_name) model.add_constraints( - x_expr.isel({sos_dim: n - 1}) - <= M.isel({sos_dim: n - 1}) * z_expr.isel({sos_dim: n - 2}), + x_expr.isel({sos_dim: n - 1}, drop=True) + <= M.isel({sos_dim: n - 1}, drop=True) + * z_expr.isel({sos_dim: n - 2}, drop=True), name=last_name, ) added_constraints.extend([last_name, card_name]) diff --git a/linopy/types.py b/linopy/types.py index 6b4cf712..21a8dddd 100644 --- a/linopy/types.py +++ b/linopy/types.py @@ -39,6 +39,10 @@ | pl.Series ) CONSTANT_TYPES: tuple[type, ...] = get_args(ConstantLike) + +UnlabeledLike: TypeAlias = numpy.ndarray | list | pl.Series +UNLABELED_TYPES = get_args(UnlabeledLike) + SignLike: TypeAlias = str | numpy.ndarray | DataArray | Series | DataFrame MaskLike: TypeAlias = numpy.ndarray | DataArray | Series | DataFrame PathLike: TypeAlias = str | Path diff --git a/linopy/variables.py b/linopy/variables.py index 0eed704d..ff01fd74 100644 --- a/linopy/variables.py +++ b/linopy/variables.py @@ -31,7 +31,7 @@ from xarray.core.utils import Frozen import linopy.expressions as expressions -from linopy.alignment import broadcast_to_coords +from linopy.alignment import as_dataarray, broadcast_to_coords from linopy.common import ( LabelPositionIndex, LocIndexer, @@ -63,6 +63,16 @@ STASHED_UPPER, TERM_DIM, ) +from linopy.semantics import ( + _legacy_coord_mismatch_message, + _legacy_masked_variable_message, + _shared_dim_mismatch_message, + check_user_nan, + enforce_aux_conflict, + first_mismatched_dim, + is_v1, + warn_legacy, +) from linopy.types import ( ConstantLike, DimsLike, @@ -328,14 +338,51 @@ def to_linexpr( linopy.LinearExpression Linear expression with the variables and coefficients. """ - coefficient = broadcast_to_coords( - coefficient, coords=self.coords, dims=self.dims, strict=False - ) - coefficient = coefficient.reindex_like(self.labels, fill_value=0) - coefficient = coefficient.fillna(0) + # Β§8: check on the raw coefficient, before broadcast aligns it away. + if not np.isscalar(coefficient): + coeff_da = as_dataarray(coefficient) + enforce_aux_conflict([self.labels, coeff_da], stacklevel=4) + mismatch = first_mismatched_dim(self.labels, coeff_da) + if mismatch is not None: + if is_v1(): + raise ValueError(_shared_dim_mismatch_message(*mismatch)) + warn_legacy( + _legacy_coord_mismatch_message( + "this operator's constant operand", *mismatch + ), + stacklevel=4, + ) + coefficient = broadcast_to_coords(coefficient, coords=self.coords, strict=False) + # Β§5: user-supplied NaN in the coefficient must raise (v1) / warn + # (legacy) β€” it's the multiplicative analogue of ``x + nan_data`` + # and otherwise enters the expression silently. The default + # coefficient ``1`` carries no NaN, so the check is a no-op there. + if coefficient.isnull().any(): + check_user_nan(op_kind="mul") + if is_v1(): + # Under v1 the LinearExpression must carry absence (NaN at + # `labels == -1`) so Β§6 propagation through downstream + # arithmetic works. + coefficient = coefficient.reindex_like(self.labels, fill_value=np.nan) + absent = self.labels == -1 + coefficient = coefficient.where(~absent) + else: + # LEGACY: warn if the variable carries absent slots β€” those + # silently contribute 0 here, but v1 will propagate the + # absence and produce a different result downstream. This is + # the origin of the most common legacy↔v1 divergence (masked + # variables in arithmetic) that no other warn-site catches. + has_absence = bool((self.labels == -1).any()) + if has_absence: + warn_legacy(_legacy_masked_variable_message(self.name)) + coefficient = coefficient.reindex_like(self.labels, fill_value=0) + coefficient = coefficient.fillna(0) ds = Dataset({"coeffs": coefficient, "vars": self.labels}).expand_dims( TERM_DIM, -1 ) + if is_v1(): + const = DataArray(np.where(absent, np.nan, 0.0), coords=self.labels.coords) + ds = ds.assign(const=const) return expressions.LinearExpression(ds, self.model) def __repr__(self) -> str: @@ -417,7 +464,6 @@ def __mul__(self, other: SideLike) -> ExpressionLike: try: if isinstance(other, Variable | ScalarVariable): return self.to_linexpr() * other - return self.to_linexpr(other) except TypeError: return NotImplemented @@ -1232,19 +1278,28 @@ def where( def fillna( self, - fill_value: ScalarVariable | dict[str, str | float | int] | Variable | Dataset, - ) -> Variable: + fill_value: int + | float + | ScalarVariable + | dict[str, str | float | int] + | Variable + | Dataset, + ) -> Variable | expressions.LinearExpression: """ - Fill missing values with a variable. + Fill missing (absent) slots. - This operation call ``xarray.DataArray.fillna`` but ensures preserving - the linopy.Variable type. + A numeric ``fill_value`` substitutes a *constant* for the absent + variable slots, so the result is a :class:`LinearExpression` (a + constant is not a variable). A Variable / ScalarVariable + ``fill_value`` keeps the result a Variable. Parameters ---------- - fill_value : Variable/ScalarVariable - Variable to use for filling. + fill_value : numeric, Variable, or ScalarVariable + Value to fill the absent slots with. """ + if isinstance(fill_value, int | float | np.integer | np.floating): + return self.to_linexpr().fillna(fill_value) return self.where(~self.isnull(), fill_value) def ffill(self, dim: str, limit: None = None) -> Variable: @@ -1356,6 +1411,48 @@ def equals(self, other: Variable) -> bool: shift = varwrap(Dataset.shift, fill_value=_fill_value) + def reindex( + self, + indexers: Mapping[Any, Any] | None = None, + **indexers_kwargs: Any, + ) -> Variable: + """ + Reindex the variable to a new set of coordinates. + + New positions are marked absent (``labels = -1``, + ``lower = upper = NaN``); existing positions are preserved. This + is one of the named mechanisms in convention.md Β§4 for creating + absence. + """ + return self.__class__( + self.data.reindex(indexers, fill_value=self._fill_value, **indexers_kwargs), + self.model, + self.name, + ) + + def reindex_like( + self, + other: Any, + **kwargs: Any, + ) -> Variable: + """ + Reindex the variable to another object's coordinates. + + New positions are marked absent with the sentinel fill values + (see :meth:`reindex`). + """ + if isinstance(other, DataArray): + ref = other.to_dataset(name="__tmp__") + elif isinstance(other, Dataset): + ref = other + else: + ref = other.data + return self.__class__( + self.data.reindex_like(ref, fill_value=self._fill_value, **kwargs), + self.model, + self.name, + ) + swap_dims = varwrap(Dataset.swap_dims) set_index = varwrap(Dataset.set_index) @@ -1366,7 +1463,10 @@ def equals(self, other: Variable) -> bool: stack = varwrap(Dataset.stack) - unstack = varwrap(Dataset.unstack) + # ``fill_value=_fill_value`` so missing (region, year) combinations end up + # as the absent-slot sentinel (labels=-1, lower=upper=NaN) instead of as + # NaN labels β€” Β§2 storage invariant + Β§4 absence-creation guarantee. + unstack = varwrap(Dataset.unstack, fill_value=_fill_value) iterate_slices = iterate_slices diff --git a/pyproject.toml b/pyproject.toml index 8ddbfa07..108d8c18 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -95,6 +95,7 @@ benchmarks = [ "dask==2026.6.0", "pytest==9.1.1", "pytest-benchmark==5.2.3", + "pytest-benchmem[plot-static]==0.4.7", "pytest-memray==1.8.0", "pytest-codspeed==5.0.3", ] diff --git a/test/conftest.py b/test/conftest.py index 067452d2..6439168f 100644 --- a/test/conftest.py +++ b/test/conftest.py @@ -3,6 +3,8 @@ from __future__ import annotations import os +import warnings +from collections.abc import Generator from typing import TYPE_CHECKING import pandas as pd @@ -11,6 +13,16 @@ if TYPE_CHECKING: from linopy import Model, Variable +# ``linopy`` is intentionally NOT imported at module level β€” doing so +# loads it from site-packages before pytest's ``--doctest-modules`` +# collection walks the ``linopy/`` source tree, and the resulting +# __file__ mismatch breaks the whole run on Windows CI (and elsewhere). +# Same reasoning as the ``filterwarnings`` comment in ``pyproject.toml``. +# Values mirror ``linopy.config.LEGACY_SEMANTICS`` / ``V1_SEMANTICS``. +_LEGACY_SEMANTICS = "legacy" +_V1_SEMANTICS = "v1" +_VALID_SEMANTICS = {_LEGACY_SEMANTICS, _V1_SEMANTICS} + def pytest_addoption(parser: pytest.Parser) -> None: """Add custom command line options.""" @@ -25,6 +37,10 @@ def pytest_addoption(parser: pytest.Parser) -> None: def pytest_configure(config: pytest.Config) -> None: """Configure pytest with custom markers and behavior.""" config.addinivalue_line("markers", "gpu: marks tests as requiring GPU hardware") + for sem in sorted(_VALID_SEMANTICS): + config.addinivalue_line( + "markers", f"{sem}: run this test only under the {sem} semantics" + ) # Set environment variable so test modules can check if GPU tests are enabled # This is needed because parametrize happens at import time @@ -63,6 +79,35 @@ def pytest_collection_modifyitems( item.add_marker(pytest.mark.gpu) +@pytest.fixture(autouse=True, params=[_LEGACY_SEMANTICS, _V1_SEMANTICS]) +def semantics(request: pytest.FixtureRequest) -> Generator[str, None, None]: + """ + Run every test under both arithmetic semantics by default. + + A test marked with a semantics name (``@pytest.mark.legacy`` or + ``@pytest.mark.v1``) runs only under that semantics. Under ``legacy``, + ``LinopySemanticsWarning`` is suppressed so test output stays clean; + ``test_convention.py`` verifies the warnings are actually emitted. + """ + # Deferred import (see top-of-file comment). + from linopy.config import LinopySemanticsWarning, options + + item = request.node + for sem in _VALID_SEMANTICS: + if item.get_closest_marker(sem) and request.param != sem: + pytest.skip(f"{sem}-only test") + + old = options["semantics"] + options["semantics"] = request.param + if request.param == _LEGACY_SEMANTICS: + with warnings.catch_warnings(): + warnings.simplefilter("ignore", LinopySemanticsWarning) + yield request.param + else: + yield request.param + options["semantics"] = old + + @pytest.fixture def m() -> Model: from linopy import Model diff --git a/test/test_alignment.py b/test/test_alignment.py index 73d2cfdb..a5a850af 100644 --- a/test/test_alignment.py +++ b/test/test_alignment.py @@ -28,7 +28,7 @@ from xarray import DataArray from xarray.testing.assertions import assert_equal -from linopy import EvolvingAPIWarning, LinearExpression, Model, Variable +from linopy import LinearExpression, LinopySemanticsWarning, Model, Variable from linopy.alignment import ( _coords_to_dict, align, @@ -649,6 +649,25 @@ def test_extra_coord_entries_broadcast_in(self) -> None: assert list(da.coords["dim_0"].values) == ["a", "b"] assert list(da.coords["dim_2"].values) == ["A", "B"] + @pytest.mark.v1 + def test_explicit_dims_bypass_size_pairing(self) -> None: + """ + An explicit ``dims`` is honored positionally, like xarray β€” even when + size-pairing would otherwise be ambiguous (both coords dims size 4). + Pins the "infer order only when the user didn't name it" rule. + """ + coords = {"a": [0, 1, 2, 3], "b": [4, 5, 6, 7]} # both size 4 + + # No dims β†’ size-pairing can't decide β†’ raises. + with pytest.raises(ValueError, match=r"sizes alone cannot decide"): + broadcast_to_coords(np.arange(4), coords=coords, strict=False) + + # Explicit dims=['a'] β†’ the axis is labeled 'a' positionally, no + # pairing, no raise (matches xarray's positional dims assignment). + da = broadcast_to_coords(np.arange(4), coords=coords, dims=["a"], strict=False) + assert set(da.dims) == {"a", "b"} + assert (da.sel(b=4).values == np.arange(4)).all() + # --------------------------------------------------------------------------- # Implicit MultiIndex-level projection β€” the legacy/v1 fork point @@ -659,12 +678,13 @@ class TestBroadcastToCoordsMultiIndexProjection: """ Inputs indexed by levels of a stacked MultiIndex dim are projected onto it. - Implicit projection is deprecated (scenario B, #732/#737): it warns under - both modes today and will raise under the v1 convention. Coverage gaps - raise under strict mode. When #717 lands, the deprecation tests here fork - into legacy (warn) and v1 (raise) variants. + Implicit projection is legacy-only (scenario B, #732/#737): under legacy + semantics it warns (LinopySemanticsWarning); under the v1 convention it + raises (sections 8 and 11). Coverage gaps raise under strict mode in both. + Tests fork accordingly: @pytest.mark.legacy (warn) / @pytest.mark.v1 (raise). """ + @pytest.mark.legacy def test_broadcasts_single_level( self, mi_coords: xr.Coordinates, by_level1: DataArray ) -> None: @@ -675,7 +695,7 @@ def test_broadcasts_single_level( 'snapshot' MultiIndex by a weighting indexed only by 'period'. Each level combination of the MultiIndex must pick up its level's value. """ - with pytest.warns(EvolvingAPIWarning, match=r"broadcasting level subset"): + with pytest.warns(LinopySemanticsWarning, match=r"broadcasting level subset"): da = broadcast_to_coords( by_level1, coords=mi_coords, dims=["dim_3"], strict=False ) @@ -687,6 +707,17 @@ def test_broadcasts_single_level( assert da.sel(dim_3=(2, "a")).item() == 20.0 assert da.sel(dim_3=(2, "b")).item() == 20.0 + @pytest.mark.v1 + def test_broadcasts_single_level_raises_v1( + self, mi_coords: xr.Coordinates, by_level1: DataArray + ) -> None: + """v1: the implicit level projection must be written explicitly.""" + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + broadcast_to_coords( + by_level1, coords=mi_coords, dims=["dim_3"], strict=False + ) + + @pytest.mark.legacy def test_stacks_full_levels(self, mi_coords: xr.Coordinates) -> None: """ A constant indexed by all MI level names stacks element-wise into the MI dim. @@ -702,7 +733,7 @@ def test_stacks_full_levels(self, mi_coords: xr.Coordinates) -> None: weights = pd.Series([10.0, 20.0], index=subset) with pytest.warns( - EvolvingAPIWarning, match=r"filling uncovered level combinations" + LinopySemanticsWarning, match=r"filling uncovered level combinations" ): da = broadcast_to_coords( weights, coords=mi_coords, dims=["dim_3"], strict=False @@ -715,6 +746,17 @@ def test_stacks_full_levels(self, mi_coords: xr.Coordinates) -> None: assert np.isnan(da.sel(dim_3=(1, "b")).item()) assert np.isnan(da.sel(dim_3=(2, "a")).item()) + @pytest.mark.v1 + def test_stacks_full_levels_raises_v1(self, mi_coords: xr.Coordinates) -> None: + """v1: implicit NaN-fill of uncovered combinations must be explicit.""" + subset = pd.MultiIndex.from_tuples( + [(1, "a"), (2, "b")], names=["level1", "level2"] + ) + weights = pd.Series([10.0, 20.0], index=subset) + + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + broadcast_to_coords(weights, coords=mi_coords, dims=["dim_3"], strict=False) + def test_full_coverage_is_silent( self, mi_coords: xr.Coordinates, mi_index: pd.MultiIndex ) -> None: @@ -723,12 +765,12 @@ def test_full_coverage_is_silent( Aligning an input that reconstructs the whole MultiIndex onto its dim is equivalent to the input already carrying that dim (future Β§11), so it must - not emit the EvolvingAPIWarning the partial/gap projections do. + not warn or raise the way partial/gap projections do. """ full = pd.Series([1.0, 2.0, 3.0, 4.0], index=mi_index) with warnings.catch_warnings(): - warnings.simplefilter("error", EvolvingAPIWarning) + warnings.simplefilter("error", LinopySemanticsWarning) da = broadcast_to_coords( full, coords=mi_coords, dims=["dim_3"], strict=False ) @@ -794,18 +836,30 @@ def test_unrelated_mi_series_still_unstacks(self) -> None: assert set(da.dims) == {"time", "foo", "bar"} + @pytest.mark.legacy def test_partially_named_mi_levels(self) -> None: """A None level name in the MultiIndex is skipped during projection.""" mi = pd.MultiIndex.from_product([[1, 2], ["a", "b"]], names=("level1", None)) mi.name = "dim_3" by_level1 = DataArray([10.0, 20.0], coords={"level1": [1, 2]}, dims=["level1"]) - with pytest.warns(EvolvingAPIWarning, match=r"broadcasting level subset"): + with pytest.warns(LinopySemanticsWarning, match=r"broadcasting level subset"): da = broadcast_to_coords(by_level1, coords={"dim_3": mi}, strict=False) assert da.dims == ("dim_3",) assert da.values.tolist() == [10.0, 10.0, 20.0, 20.0] + @pytest.mark.v1 + def test_partially_named_mi_levels_raises_v1(self) -> None: + """v1: projection onto a partially named MultiIndex raises like any other.""" + mi = pd.MultiIndex.from_product([[1, 2], ["a", "b"]], names=("level1", None)) + mi.name = "dim_3" + by_level1 = DataArray([10.0, 20.0], coords={"level1": [1, 2]}, dims=["level1"]) + + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + broadcast_to_coords(by_level1, coords={"dim_3": mi}, strict=False) + + @pytest.mark.legacy def test_gap_detection_with_extra_dims(self, mi_coords: xr.Coordinates) -> None: """Gaps are detected per level combination even when the input has extra dims.""" arr = DataArray( @@ -815,7 +869,7 @@ def test_gap_detection_with_extra_dims(self, mi_coords: xr.Coordinates) -> None: ) with pytest.warns( - EvolvingAPIWarning, match=r"filling uncovered level combinations" + LinopySemanticsWarning, match=r"filling uncovered level combinations" ): da = broadcast_to_coords( arr, coords=mi_coords, dims=["dim_3"], strict=False @@ -823,6 +877,20 @@ def test_gap_detection_with_extra_dims(self, mi_coords: xr.Coordinates) -> None: assert set(da.dims) == {"dim_3", "extra"} + @pytest.mark.v1 + def test_gap_detection_with_extra_dims_raises_v1( + self, mi_coords: xr.Coordinates + ) -> None: + """v1: gap-producing projections raise regardless of extra dims.""" + arr = DataArray( + [[[1.0, np.nan], [2.0, 2.0]], [[3.0, 3.0], [4.0, 4.0]]], + dims=["level1", "level2", "extra"], + coords={"level1": [1, 2], "level2": ["a", "b"], "extra": [0, 1]}, + ) + + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + broadcast_to_coords(arr, coords=mi_coords, dims=["dim_3"], strict=False) + def test_strict_gap_error_truncates_long_missing_list(self) -> None: """More than 5 missing combinations are truncated in the error message.""" idx = pd.MultiIndex.from_product( @@ -843,6 +911,7 @@ def test_strict_gap_error_truncates_long_missing_list(self) -> None: # --- strict-mode policy on MI projections (deprecation / gaps) --- + @pytest.mark.legacy def test_strict_partial_level_warns( self, mi_coords: xr.Coordinates, by_level1: DataArray ) -> None: @@ -850,10 +919,10 @@ def test_strict_partial_level_warns( Per-level bounds broadcast across the MI dim, with the deprecation warning. Scenario B (#732 / #737 discussion): implicit MI-level projection is - deprecated everywhere, including the strict (bounds/mask) path, and will - raise under the v1 convention. + deprecated everywhere, including the strict (bounds/mask) path, and + raises under the v1 convention. """ - with pytest.warns(EvolvingAPIWarning, match=r"broadcasting level subset"): + with pytest.warns(LinopySemanticsWarning, match=r"broadcasting level subset"): da = broadcast_to_coords( by_level1, mi_coords, dims=["dim_3"], label="lower bound" ) @@ -861,6 +930,17 @@ def test_strict_partial_level_warns( assert da.sel(dim_3=(1, "b")).item() == 10.0 assert da.sel(dim_3=(2, "a")).item() == 20.0 + @pytest.mark.v1 + def test_strict_partial_level_raises_v1( + self, mi_coords: xr.Coordinates, by_level1: DataArray + ) -> None: + """v1: per-level bounds must be projected onto the dim explicitly.""" + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + broadcast_to_coords( + by_level1, mi_coords, dims=["dim_3"], label="lower bound" + ) + + @pytest.mark.legacy def test_strict_rejects_coverage_gap(self, mi_coords: xr.Coordinates) -> None: """A coverage gap warns on the broadcast rung but raises on the strict rung.""" subset = pd.MultiIndex.from_tuples( @@ -869,13 +949,29 @@ def test_strict_rejects_coverage_gap(self, mi_coords: xr.Coordinates) -> None: weights = pd.Series([10.0, 20.0], index=subset) with pytest.warns( - EvolvingAPIWarning, match=r"filling uncovered level combinations" + LinopySemanticsWarning, match=r"filling uncovered level combinations" ): broadcast_to_coords(weights, coords=mi_coords, dims=["dim_3"], strict=False) with pytest.raises(ValueError, match=r"no value for .* level combination"): broadcast_to_coords(weights, mi_coords, dims=["dim_3"], label="lower bound") + @pytest.mark.v1 + def test_strict_rejects_coverage_gap_v1(self, mi_coords: xr.Coordinates) -> None: + """v1: a coverage gap raises on both rungs (different errors, same outcome).""" + subset = pd.MultiIndex.from_tuples( + [(1, "a"), (2, "b")], names=["level1", "level2"] + ) + weights = pd.Series([10.0, 20.0], index=subset) + + # Broadcast rung: the implicit projection itself is rejected. + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + broadcast_to_coords(weights, coords=mi_coords, dims=["dim_3"], strict=False) + + # Strict rung: the coverage gap is rejected before the projection policy. + with pytest.raises(ValueError, match=r"no value for .* level combination"): + broadcast_to_coords(weights, mi_coords, dims=["dim_3"], label="lower bound") + def test_strict_rejects_unnamed_mi_mismatch(self) -> None: """ A MultiIndex input with unnamed levels cannot be projected by level name, diff --git a/test/test_constraint.py b/test/test_constraint.py index 10febf15..f74bddc7 100644 --- a/test/test_constraint.py +++ b/test/test_constraint.py @@ -599,6 +599,7 @@ def test_constraint_rhs_setter_broadcasts_missing_dim() -> None: assert (con.rhs.sel(i=1) == 2.0).all() +@pytest.mark.legacy def test_constraint_rhs_setter_projects_multiindex_level() -> None: """ Rhs indexed by one MultiIndex level is projected onto the stacked dim. @@ -617,13 +618,30 @@ def test_constraint_rhs_setter_projects_multiindex_level() -> None: rhs_by_level = xr.DataArray( [10.0, 20.0], coords={"level1": [1, 2]}, dims=["level1"] ) - with pytest.warns(linopy.EvolvingAPIWarning, match="broadcasting level subset"): + with pytest.warns(linopy.LinopySemanticsWarning, match="broadcasting level subset"): con.rhs = rhs_by_level assert con.rhs.sel(dim_3=(1, "b")).item() == 10.0 assert con.rhs.sel(dim_3=(2, "a")).item() == 20.0 +@pytest.mark.v1 +def test_constraint_rhs_setter_mi_level_raises_v1() -> None: + """v1: an rhs indexed by an MI level must be projected explicitly.""" + idx = pd.MultiIndex.from_product([[1, 2], ["a", "b"]], names=("level1", "level2")) + idx.name = "dim_3" + coords = xr.Coordinates.from_pandas_multiindex(idx, "dim_3") + m = Model() + x = m.add_variables(coords=coords, name="x") + con = m.add_constraints(1 * x >= 0, name="con") + + rhs_by_level = xr.DataArray( + [10.0, 20.0], coords={"level1": [1, 2]}, dims=["level1"] + ) + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + con.rhs = rhs_by_level + + def test_constraint_labels_setter_invalid(c: linopy.constraints.CSRConstraint) -> None: # Test that assigning labels raises AttributeError (Constraint is frozen) with pytest.raises(AttributeError): diff --git a/test/test_constraints.py b/test/test_constraints.py index acc41b2e..2de429d1 100644 --- a/test/test_constraints.py +++ b/test/test_constraints.py @@ -144,7 +144,6 @@ def test_constraint_assignment_with_reindex() -> None: @pytest.mark.parametrize( "rhs_factory", [ - pytest.param(lambda m, v: v, id="numpy"), pytest.param(lambda m, v: xr.DataArray(v, dims=["dim_0"]), id="dataarray"), pytest.param(lambda m, v: pd.Series(v, index=v), id="series"), pytest.param( @@ -168,10 +167,37 @@ def test_constraint_rhs_lower_dim(rhs_factory: Any) -> None: assert c.shape == (10, 10) +@pytest.mark.legacy +def test_constraint_rhs_unlabeled_lower_dim_legacy() -> None: + # An unlabeled array rhs pairs positionally with the leading dim under + # legacy; both dims are size 10 so the pairing is a size-coincidence. + m = Model() + naxis = np.arange(10, dtype=float) + maxis = np.arange(10).astype(str) + x = m.add_variables(coords=[naxis, maxis]) + y = m.add_variables(coords=[naxis, maxis]) + + c = m.add_constraints(x - y >= naxis) + assert c.shape == (10, 10) + + +@pytest.mark.v1 +def test_constraint_rhs_unlabeled_lower_dim_ambiguous_raises_v1() -> None: + # v1: both dims are size 10, so an unlabeled length-10 rhs cannot be + # paired by size β€” it raises (wrap in a DataArray to disambiguate). + m = Model() + naxis = np.arange(10, dtype=float) + maxis = np.arange(10).astype(str) + x = m.add_variables(coords=[naxis, maxis]) + y = m.add_variables(coords=[naxis, maxis]) + + with pytest.raises(ValueError, match=r"sizes alone cannot decide"): + m.add_constraints(x - y >= naxis) + + @pytest.mark.parametrize( "rhs_factory", [ - pytest.param(lambda m: np.ones((5, 3)), id="numpy"), pytest.param(lambda m: pd.DataFrame(np.ones((5, 3))), id="dataframe"), ], ) @@ -186,6 +212,29 @@ def test_constraint_rhs_higher_dim_constant_warns( assert "dimensions" in caplog.text +@pytest.mark.legacy +def test_constraint_rhs_unlabeled_higher_dim_warns_legacy(caplog: Any) -> None: + # Legacy: an unlabeled (5, 3) rhs pairs axis 0 with dim_0 positionally and + # broadcasts the extra axis, warning about the extra dimension. + m = Model() + x = m.add_variables(coords=[range(5)], name="x") + + with caplog.at_level("WARNING", logger="linopy.expressions"): + m.add_constraints(x >= np.ones((5, 3))) + assert "dimensions" in caplog.text + + +@pytest.mark.v1 +def test_constraint_rhs_unlabeled_higher_dim_raises_v1() -> None: + # v1: the (5, 3) rhs has an axis of length 3 matching no dim of x β€” it + # cannot be paired by size and raises. + m = Model() + x = m.add_variables(coords=[range(5)], name="x") + + with pytest.raises(ValueError, match=r"no unambiguous dimension match"): + m.add_constraints(x >= np.ones((5, 3))) + + def test_constraint_rhs_higher_dim_dataarray_reindexes() -> None: """DataArray RHS with extra dims reindexes to expression coords (no raise).""" m = Model() @@ -355,7 +404,13 @@ def test_sanitize_infinities() -> None: m.add_constraints(y <= -np.inf, name="con_wrong_neg_inf") +@pytest.mark.legacy class TestConstraintCoordinateAlignment: + """ + Constraint-side counterpart of TestCoordinateAlignment β€” legacy-only; + v1 raises on subset/superset/NaN-RHS (see convention.md Β§5/Β§8/Β§12). + """ + @pytest.fixture(params=["xarray", "pandas_series"], ids=["da", "series"]) def subset(self, request: Any) -> xr.DataArray | pd.Series: if request.param == "xarray": diff --git a/test/test_convention.py b/test/test_convention.py new file mode 100644 index 00000000..b4f2c3e5 --- /dev/null +++ b/test/test_convention.py @@ -0,0 +1,29 @@ +"""Tests for the v1 semantics option and the test harness.""" + +from __future__ import annotations + +import pytest + +import linopy + + +class TestSemanticsOption: + def test_default_is_legacy(self) -> None: + linopy.options.reset() + assert linopy.options["semantics"] == "legacy" + + def test_invalid_value_raises(self) -> None: + with pytest.raises(ValueError, match="Invalid semantics"): + linopy.options["semantics"] = "exact" + + +class TestHarness: + """The autouse ``semantics`` fixture and the legacy / v1 markers.""" + + @pytest.mark.legacy + def test_legacy_marker(self) -> None: + assert linopy.options["semantics"] == "legacy" + + @pytest.mark.v1 + def test_v1_marker(self) -> None: + assert linopy.options["semantics"] == "v1" diff --git a/test/test_dualization.py b/test/test_dualization.py index c819b3fe..f7bb1c43 100644 --- a/test/test_dualization.py +++ b/test/test_dualization.py @@ -491,13 +491,19 @@ def test_strong_duality_array_variable() -> None: m = Model() x = m.add_variables(lower=0, coords=[pd.RangeIndex(n)], name="x") - for i in range(4): - lhs: Any = sum(float(A[i, j]) * x[j] for j in range(n)) - m.add_constraints(lhs <= float(b[i]), name=f"c{i}") - + # Label the matrix: the column axis shares x's dim, the row axis is the + # constraint index. + A_da = xr.DataArray( + A, + dims=["con", "dim_0"], + coords={"con": pd.RangeIndex(4), "dim_0": x.indexes["dim_0"]}, + ) + lhs = (A_da * x).sum("dim_0") # expression indexed by `con` # Bounded because x >= 0 and A x <= b with A > 0. - obj: Any = sum(float(c_obj[j]) * x[j] for j in range(n)) - m.add_objective(obj, sense="max") + m.add_constraints(lhs <= b, name="c") # raw numpy `b` pairs with `con` by size + + c_da = xr.DataArray(c_obj, dims=["dim_0"], coords={"dim_0": x.indexes["dim_0"]}) + m.add_objective((c_da * x).sum(), sense="max") primal_obj = _solve(m) dual_obj = _solve(m.dualize()) diff --git a/test/test_legacy_violations.py b/test/test_legacy_violations.py new file mode 100644 index 00000000..0247571d --- /dev/null +++ b/test/test_legacy_violations.py @@ -0,0 +1,2555 @@ +""" +Legacy convention violations and v1 fixes. + +Pairs ``@pytest.mark.legacy`` tests that document the surprising legacy +behaviour against ``@pytest.mark.v1`` tests that pin the v1 fix. Each +class corresponds to a section of ``arithmetics-design/convention.md`` +and to one or more linopy bug reports. + +Slice A β€” constant operand path (Β§5, Β§8, Β§9): + Β§8 Shared dimensions must match exactly β†’ #708 / #586 / #550 + Β§5 User-supplied NaN raises β†’ #713 / PyPSA #1683 + Β§9 Non-shared dimensions broadcast β†’ (positive regression guard) + +Slice B β€” expression-OP-expression / variable-OP-variable (Β§8 via `merge`): + Β§8 Shared dimensions must match exactly β†’ #708 / #570 (expr+expr branch) + +Slice C β€” absence propagation (Β§3, Β§6, Β§7): + Β§6 Absence propagates through every operator β†’ #712 (absent-as-zero) + Β§3 isnull() is the unifying predicate β†’ #711 + Β§7 fillna()/.where() resolve absence β†’ (positive coverage) + +Slice D β€” Variable.reindex / .reindex_like (Β§4 absence creation): + Β§4 Reindexing extends coords and marks new slots absent + +Slice E β€” named-method join= + constraint RHS (Β§10, Β§12): + Β§10 .add/.sub/.mul/.div/.le/.ge/.eq accept explicit join= + Β§12 NaN in constraint RHS raises (v1) β†’ PyPSA #1683 + Β§12 Coord mismatch in RHS raises (v1) β†’ #707 + +Slice G β€” reductions skip absent slots (Β§13): + Β§13 sum / groupby.sum skip absent, sum of none is the zero expression + +Slice F β€” auxiliary-coordinate conflicts (Β§11): + Β§11 Non-dim coord conflict raises (v1) β†’ #295 + Β§11 Non-conflicting aux coords propagate through arithmetic + +Slice H β€” unlabeled-operand pairing (coordinate-alignment intro): + Unlabeled operands (numpy / list / polars) pair with the linopy + operand's dims by size; ambiguity or no-match raises (v1) β†’ #736 + +Slice H β€” object scope (convention preamble): + Non-linopy operands behave exactly like constant-only expressions +""" + +from __future__ import annotations + +import contextlib +import operator +import warnings +from collections.abc import Generator +from typing import Any, cast + +import numpy as np +import pandas as pd +import polars as pl +import pytest +import xarray as xr + +from linopy import LinearExpression, Model, QuadraticExpression +from linopy.config import LinopySemanticsWarning +from linopy.expressions import merge +from linopy.testing import assert_conequal, assert_linequal, assert_quadequal +from linopy.variables import Variable + + +@pytest.fixture +def m() -> Model: + return Model() + + +@pytest.fixture +def time() -> pd.RangeIndex: + return pd.RangeIndex(5, name="time") + + +@pytest.fixture +def x(m: Model, time: pd.RangeIndex) -> Variable: + return m.add_variables(lower=0, coords=[time], name="x") + + +@pytest.fixture +def unsilenced() -> Generator[None, None, None]: + """Drop the autouse fixture's LinopySemanticsWarning filter for one test.""" + with warnings.catch_warnings(): + warnings.simplefilter("always", LinopySemanticsWarning) + yield + + +# ===================================================================== +# Β§8 β€” Shared dimensions must match exactly (constant operand) +# ===================================================================== + + +_OPS = { + "add": operator.add, + "sub": operator.sub, + "mul": operator.mul, + "div": operator.truediv, +} + + +class TestExactAlignmentConstant: + @pytest.mark.v1 + @pytest.mark.parametrize("op", ["add", "sub", "mul", "div"]) + def test_same_size_different_labels_raises(self, x: Variable, op: str) -> None: + """ + #708 / #550 β€” same shape, different labels: legacy aligns by + position; v1 raises. Holds for every binary operator. + """ + other = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + with pytest.raises(ValueError, match="Coordinate mismatch on shared dimension"): + _OPS[op](x, other) + + @pytest.mark.v1 + @pytest.mark.parametrize("op", ["add", "sub", "mul", "div"]) + def test_subset_constant_raises(self, x: Variable, op: str) -> None: + """ + #711 / #708 β€” constant covers only some of the variable's + coords. Legacy left-joins (silently drops the gap); v1 raises. + """ + subset = xr.DataArray( + [10.0, 20.0], dims=["time"], coords={"time": pd.Index([1, 3], name="time")} + ) + with pytest.raises(ValueError, match="Coordinate mismatch on shared dimension"): + _OPS[op](x, subset) + + @pytest.mark.legacy + def test_add_same_size_different_labels_silent(self, x: Variable) -> None: + """Document the legacy behaviour: silent positional alignment.""" + other = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + # Legacy keeps left coords; the user's intended pairing by label is lost. + result = x + other + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + assert result.const.values.tolist() == [1.0, 2.0, 3.0, 4.0, 5.0] + + @pytest.mark.legacy + def test_add_subset_constant_silent(self, x: Variable) -> None: + """Document the legacy behaviour: silent left-join (gaps β†’ 0).""" + subset = xr.DataArray( + [10.0, 20.0], dims=["time"], coords={"time": pd.Index([1, 3], name="time")} + ) + result = x + subset + # Legacy reindex_like fills the missing positions with 0 (additive fill). + assert result.const.sel(time=0).item() == 0.0 + assert result.const.sel(time=1).item() == 10.0 + assert result.const.sel(time=3).item() == 20.0 + + +class TestBroadcastNonSharedDim: + """ + Β§9 β€” a dimension that exists only in one operand broadcasts freely. + Runs under both semantics: this is unchanged behaviour. + """ + + def test_add_broadcast_introduces_new_dim(self, x: Variable) -> None: + bcast = xr.DataArray( + [10.0, 20.0], dims=["scenario"], coords={"scenario": [0, 1]} + ) + result = x + bcast + assert set(result.const.dims) == {"time", "scenario"} + assert result.const.sizes == {"time": 5, "scenario": 2} + + def test_mul_broadcast_introduces_new_dim(self, x: Variable) -> None: + bcast = xr.DataArray([2.0, 3.0], dims=["scenario"], coords={"scenario": [0, 1]}) + result = x * bcast + assert set(result.coeffs.dims) == {"time", "scenario", "_term"} + + +# ===================================================================== +# Coordinate-alignment intro β€” unlabeled operands pair by size (#736) +# ===================================================================== + + +# The three unlabeled types must behave identically; each constructor takes a +# value sequence so the same cases parametrize over numpy / list / polars. +UNLABELED_1D = [ + pytest.param(lambda v: np.asarray(v, dtype=float), id="numpy"), + pytest.param(lambda v: [float(x) for x in v], id="list"), + pytest.param(lambda v: pl.Series([float(x) for x in v]), id="polars"), +] + + +class TestUnlabeledPairing: + """ + Unlabeled operands (numpy arrays, lists, polars Series) carry no labels, + so they pair with the linopy operand's dims *by size*. Ambiguity (two + dims share the size, or the array is square) or no size match raises + under v1; legacy pairs with the leading dims positionally and warns when + the v1 pairing would differ or reject. + """ + + @pytest.fixture + def xy(self) -> Variable: + # dims of distinct sizes so a 1-d operand pairs unambiguously + m = Model() + return m.add_variables( + coords=[pd.RangeIndex(3, name="a"), pd.RangeIndex(4, name="b")], name="xy" + ) + + @pytest.fixture + def square(self) -> Variable: + # both dims size 4 β†’ a 1-d length-4 operand is ambiguous + m = Model() + return m.add_variables( + coords=[pd.RangeIndex(4, name="p"), pd.RangeIndex(4, name="q")], name="sq" + ) + + @pytest.fixture + def wide(self) -> Variable: + # four dims of distinct sizes β€” a lower-rank operand pairs a subset + m = Model() + return m.add_variables( + coords=[ + pd.RangeIndex(3, name="a"), + pd.RangeIndex(4, name="b"), + pd.RangeIndex(5, name="c"), + pd.RangeIndex(6, name="d"), + ], + name="wide", + ) + + # -- 1-d operands ----------------------------------------------------- + + @pytest.mark.v1 + @pytest.mark.parametrize("make", UNLABELED_1D) + def test_v1_pairs_by_size(self, xy: Variable, make: Any) -> None: + # length-4 array pairs with dim "b" (size 4), not the leading "a" (3) + result = (1 * xy) + make(range(4)) + assert set(result.const.dims) == {"a", "b"} + assert result.const.sizes == {"a": 3, "b": 4} + + @pytest.mark.v1 + @pytest.mark.parametrize("make", UNLABELED_1D) + def test_v1_ambiguous_square_raises(self, square: Variable, make: Any) -> None: + with pytest.raises(ValueError, match=r"sizes alone cannot decide"): + (1 * square) + make(range(4)) + + @pytest.mark.v1 + @pytest.mark.parametrize("make", UNLABELED_1D) + def test_v1_no_size_match_raises(self, xy: Variable, make: Any) -> None: + with pytest.raises(ValueError, match=r"no unambiguous dimension match"): + (1 * xy) + make(range(7)) + + def test_dataarray_wrapping_resolves_ambiguity(self, square: Variable) -> None: + # the documented escape hatch: name the axis with a DataArray + result = (1 * square) + xr.DataArray(np.arange(4.0), dims=["p"]) + assert set(result.const.dims) == {"p", "q"} + + # -- multi-dim operands (numpy only β€” list / polars are 1-d) ---------- + + @pytest.mark.v1 + def test_v1_size_order_independent(self, xy: Variable) -> None: + # a 2-d (4, 3) operand pairs (b, a) by size regardless of axis order. + result = (1 * xy) + np.ones((4, 3)) + assert result.const.sizes == {"a": 3, "b": 4} + + @pytest.mark.v1 + def test_v1_multidim_square_ambiguous_raises(self, square: Variable) -> None: + # a 2-d (4, 4) operand against dims (p: 4, q: 4) β€” sizes cannot tell + # (p, q) from (q, p). Exercises the multi-axis-same-size branch the 1-d + # cases never reach. + with pytest.raises(ValueError, match=r"sizes alone cannot decide"): + (1 * square) + np.ones((4, 4)) + + @pytest.mark.v1 + def test_v1_multidim_no_size_match_raises(self, xy: Variable) -> None: + # a 2-d (5, 3) operand against (a: 3, b: 4) β€” the length-5 axis matches + # no dim. + with pytest.raises(ValueError, match=r"no unambiguous dimension match"): + (1 * xy) + np.ones((5, 3)) + + @pytest.mark.v1 + def test_v1_lower_rank_operand_pairs_subset_and_broadcasts( + self, wide: Variable + ) -> None: + # a 2-d (4, 5) operand against four dims pairs its axes with (b, c) by + # size and broadcasts over the unpaired (a, d). + result = (1 * wide) + np.ones((4, 5)) + assert result.const.sizes == {"a": 3, "b": 4, "c": 5, "d": 6} + + @pytest.mark.v1 + def test_v1_ambiguous_axis_within_higher_rank_raises(self) -> None: + # a 2-d (4, 5) operand where the length-4 axis matches two dims of the + # operand (a: 4, b: 4) β€” ambiguous even though the length-5 axis is + # unique. + m = Model() + y = m.add_variables( + coords=[ + pd.RangeIndex(4, name="a"), + pd.RangeIndex(4, name="b"), + pd.RangeIndex(5, name="c"), + ], + name="y", + ) + with pytest.raises(ValueError, match=r"sizes alone cannot decide"): + (1 * y) + np.ones((4, 5)) + + # -- matmul ----------------------------------------------------------- + + @pytest.mark.v1 + @pytest.mark.parametrize("make", UNLABELED_1D) + def test_v1_matmul_pairs_by_size(self, xy: Variable, make: Any) -> None: + # matmul contracts the paired dim: length-4 array pairs with "b" + result = (1 * xy) @ make(range(4)) + assert set(result.coord_dims) == {"a"} + + # -- construction (add_variables bounds) ------------------------------ + + @pytest.mark.v1 + def test_v1_add_variables_bound_pairs_by_size(self) -> None: + # The rule is the same for construction inputs: a bare-numpy bound + # pairs with the matching dim by size, not positionally (where the + # length-5 array would hit the leading dim "a" and conflict). + m = Model() + x = m.add_variables( + coords=[pd.RangeIndex(4, name="a"), pd.RangeIndex(5, name="time")], + lower=np.arange(5.0), + name="x", + ) + assert dict(x.lower.sizes) == {"a": 4, "time": 5} + assert (x.lower.isel(a=0).values == np.arange(5.0)).all() + + @pytest.mark.v1 + def test_v1_add_variables_ambiguous_bound_raises(self) -> None: + m = Model() + with pytest.raises(ValueError, match=r"sizes alone cannot decide"): + m.add_variables( + coords=[pd.RangeIndex(4, name="p"), pd.RangeIndex(4, name="q")], + lower=np.arange(4.0), + name="x", + ) + + @pytest.mark.legacy + def test_legacy_no_divergence_is_silent( + self, xy: Variable, unsilenced: None + ) -> None: + # length-3 array: legacy pairs positionally with the leading dim "a"; + # v1 would pair it with "a" too (only "a" is size 3), so there is no + # divergence and no warning. Pins the silent case. + with warnings.catch_warnings(): + warnings.simplefilter("error", LinopySemanticsWarning) + result = (1 * xy) + np.arange(3.0) + assert result.const.sizes == {"a": 3, "b": 4} + assert (result.const.isel(b=0).values == np.arange(3.0)).all() + + @pytest.mark.legacy + def test_legacy_warns_when_v1_would_differ( + self, xy: Variable, unsilenced: None + ) -> None: + # length-4 array: legacy pairs positionally with "a" (size 3) β†’ error, + # but the warning fires first explaining the v1 divergence. + # legacy positional pairing assigns the len-4 array to "a" (size 3), + # which then conflicts β€” assert it's that shape error, nothing incidental. + with pytest.warns(LinopySemanticsWarning, match=r"pairs by size instead"): + with contextlib.suppress(ValueError): + (1 * xy) + np.arange(4.0) + + @pytest.mark.legacy + def test_legacy_ambiguous_pairs_positionally_with_warning( + self, square: Variable, unsilenced: None + ) -> None: + # The square (p:4, q:4) case where v1 *raises* β€” legacy must instead + # pair positionally with the leading dim and warn, never raise. This is + # the biggest legacy/v1 divergence and the strongest no-regression guard. + with pytest.warns(LinopySemanticsWarning, match=r"this raises"): + result = (1 * square) + np.arange(4.0) + assert result.const.sizes == {"p": 4, "q": 4} + # paired with the leading dim "p": the array varies along p, broadcast over q + assert (result.const.isel(q=0).values == np.arange(4.0)).all() + + @pytest.mark.legacy + def test_legacy_add_variables_bound_positional(self, unsilenced: None) -> None: + # Regression guard for the unification: legacy add_variables still + # assigns an unlabeled bound positionally (here unambiguous β€” only "a" + # is size 5 β€” so it is also silent), producing the pre-#736 result. + m = Model() + with warnings.catch_warnings(): + warnings.simplefilter("error", LinopySemanticsWarning) + x = m.add_variables( + coords=[pd.RangeIndex(5, name="a"), pd.RangeIndex(3, name="b")], + lower=np.arange(5.0), + name="x", + ) + assert dict(x.lower.sizes) == {"a": 5, "b": 3} + assert (x.lower.isel(b=0).values == np.arange(5.0)).all() + + +# ===================================================================== +# Β§5 β€” User-supplied NaN raises (covers #713 and PyPSA #1683) +# ===================================================================== + + +class TestUserNaNRaises: + @pytest.mark.v1 + @pytest.mark.parametrize("op", ["add", "sub", "mul", "div"]) + def test_nan_dataarray_raises( + self, x: Variable, time: pd.RangeIndex, op: str + ) -> None: + # Use [2, NaN, 3, 4, 5] so div doesn't trip on a 0 divisor at slot 0. + nan_data = xr.DataArray( + [2.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + with pytest.raises(ValueError, match="NaN"): + _OPS[op](x, nan_data) + + @pytest.mark.v1 + @pytest.mark.parametrize("op", ["add", "sub", "mul", "div"]) + def test_nan_scalar_raises(self, x: Variable, op: str) -> None: + with pytest.raises(ValueError, match="NaN"): + _OPS[op](x, float("nan")) + + @pytest.mark.v1 + def test_pypsa_1683_inf_times_zero_raises( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """ + PyPSA #1683 β€” ``min_pu * nominal_fix`` with ``p_nom=inf`` and + ``p_min_pu=0`` yields a NaN bound. v1 surfaces this at construction, + not as a downstream solve failure. + """ + nominal_fix = xr.DataArray( + [np.inf, np.inf, np.inf, np.inf, np.inf], + dims=["time"], + coords={"time": time}, + ) + min_pu = xr.DataArray( + [1.0, 0.0, 1.0, 1.0, 1.0], dims=["time"], coords={"time": time} + ) + bound = min_pu * nominal_fix # 0 * inf = NaN at time=1 + assert np.isnan(bound.values[1]) + with pytest.raises(ValueError, match="NaN"): + x * bound + + @pytest.mark.v1 + def test_to_linexpr_coefficient_nan_raises( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """ + Β§5 on the direct ``Variable.to_linexpr(coefficient)`` entry β€” + callers can construct an expression bypassing the operator + overloads. NaN in the explicit coefficient is still user data + and must raise (otherwise the NaN flows into ``coeffs`` and + Β§6 silently propagates absence from what was actually a + data error). + """ + nan_coeff = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + with pytest.raises(ValueError, match="NaN"): + x.to_linexpr(nan_coeff) + + @pytest.mark.legacy + def test_add_nan_dataarray_silently_fills_with_zero( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """Document legacy: NaN in addend silently becomes 0 (#713).""" + nan_data = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + result = x + nan_data + assert result.const.sel(time=1).item() == 0.0 # NaN β†’ 0 + + @pytest.mark.legacy + def test_mul_nan_dataarray_silently_fills_with_zero( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """ + Document legacy: NaN in multiplier silently becomes 0 β€” variable + zeroed out at that slot (#713). + """ + nan_data = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + result = x * nan_data + assert result.coeffs.squeeze().sel(time=1).item() == 0.0 + + +# ===================================================================== +# Legacy emits LinopySemanticsWarning where v1 would diverge +# ===================================================================== + + +def _one_legacy_warning(*ops) -> str: # type: ignore[no-untyped-def] + """ + Run ``ops`` (a series of callables) under fresh warning capture + and return the first ``LinopySemanticsWarning``'s text. Test + helper β€” keeps each test focused on the message, not the plumbing. + """ + with warnings.catch_warnings(record=True) as ws: + warnings.simplefilter("always") + for op in ops: + op() + legacy = [w for w in ws if issubclass(w.category, LinopySemanticsWarning)] + assert legacy, "expected at least one LinopySemanticsWarning" + return str(legacy[0].message) + + +# Common tail shared by every legacy warning β€” separated so each +# expected-message test can focus on the part that's specific to the +# rule, without 4 lines of boilerplate per test. +_OPT_IN_HINT = ( + "\n Opt in: linopy.options['semantics'] = 'v1'" + "\n Silence: warnings.filterwarnings('ignore', " + "category=LinopySemanticsWarning)" +) + + +class TestLegacyWarning: + """ + Asserts the *full text* of each legacy warning. The point: the + test reads like a spec β€” a reviewer judges the message's helpfulness + by reading the test, and any change to a message surfaces as a diff + in the test. Goal #2 (actionable warnings) lives or dies here. + """ + + @pytest.mark.legacy + def test_coord_mismatch_const_operand_same_size( + self, x: Variable, unsilenced: None + ) -> None: + other = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + msg = _one_legacy_warning(lambda: x + other) + assert msg == ( + "Coordinate mismatch in this operator's constant operand " + "silently aligned by legacy (positional when sizes match, " + "otherwise left-join). Under v1 this raises ValueError." + "\n Dim: 'time': left=[0, 1, 2, 3, 4], " + "right=[10, 11, 12, 13, 14]" + "\n Resolve: `.sel(...)` / `.reindex(...)` to align" + "\n `.assign_coords(...)` to relabel one side" + "\n `linopy.align(...)` to pre-align several operands" + "\n or pass an explicit `join=` argument." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_coord_mismatch_const_operand_subset( + self, x: Variable, unsilenced: None + ) -> None: + subset = xr.DataArray( + [10.0, 20.0], dims=["time"], coords={"time": pd.Index([1, 3], name="time")} + ) + msg = _one_legacy_warning(lambda: x + subset) + assert msg == ( + "Coordinate mismatch in this operator's constant operand " + "silently aligned by legacy (positional when sizes match, " + "otherwise left-join). Under v1 this raises ValueError." + "\n Dim: 'time': left=[0, 1, 2, 3, 4], right=[1, 3]" + "\n Resolve: `.sel(...)` / `.reindex(...)` to align" + "\n `.assign_coords(...)` to relabel one side" + "\n `linopy.align(...)` to pre-align several operands" + "\n or pass an explicit `join=` argument." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_coord_mismatch_merge_path( + self, m: Model, time: pd.RangeIndex, unsilenced: None + ) -> None: + x_local = m.add_variables(lower=0, coords=[time], name="x_local") + other = m.add_variables( + lower=0, coords=[pd.Index([10, 11, 12, 13, 14], name="time")], name="other" + ) + msg = _one_legacy_warning(lambda: x_local + other) + assert msg == ( + "Coordinate mismatch in merge along dim '_term' silently " + "aligned by legacy (positional when sizes match, otherwise " + "left-join). Under v1 this raises ValueError." + "\n Dim: 'time': left=[0, 1, 2, 3, 4], " + "right=[10, 11, 12, 13, 14]" + "\n Resolve: `.sel(...)` / `.reindex(...)` to align" + "\n `.assign_coords(...)` to relabel one side" + "\n `linopy.align(...)` to pre-align several operands" + "\n or pass an explicit `join=` argument." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_nan_addend( + self, x: Variable, time: pd.RangeIndex, unsilenced: None + ) -> None: + nan_data = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + msg = _one_legacy_warning(lambda: x + nan_data) + assert msg == ( + "NaN in the constant operand was silently treated as 0 by " + "legacy (additive identity). Under v1 this raises ValueError." + "\n Resolve: `.fillna(value)` (data error)" + "\n or `mask=` / `.where(cond)` / `.reindex(...)` " + "on the variable (intended absence)." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_nan_multiplier( + self, x: Variable, time: pd.RangeIndex, unsilenced: None + ) -> None: + nan_factor = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + msg = _one_legacy_warning(lambda: x * nan_factor) + assert msg == ( + "NaN in the multiplicative factor was silently treated as 0 " + "by legacy (so the variable was zeroed out at that slot). " + "Under v1 this raises ValueError." + "\n Resolve: `.fillna(value)` (data error)" + "\n or `mask=` / `.where(cond)` / `.reindex(...)` " + "on the variable (intended absence)." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_nan_divisor( + self, x: Variable, time: pd.RangeIndex, unsilenced: None + ) -> None: + nan_divisor = xr.DataArray( + [2.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + msg = _one_legacy_warning(lambda: x / nan_divisor) + assert msg == ( + "NaN in the divisor was silently treated as 1 by legacy (a " + "different fill from `+`/`*` which use 0). Under v1 this " + "raises ValueError." + "\n Resolve: `.fillna(value)` (data error)" + "\n or `mask=` / `.where(cond)` / `.reindex(...)` " + "on the variable (intended absence)." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_aux_conflict(self, m: Model, unsilenced: None) -> None: + A = pd.Index([1, 2, 3], name="A") + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + msg = _one_legacy_warning(lambda: v + const) + assert msg == ( + "Auxiliary coordinate 'B' was conflicting across operands " + "and silently dropped by legacy (xarray's default). Under v1 " + "this raises ValueError." + "\n Values: left=[311, 311, 322], right=[400, 400, 500]" + "\n Resolve: `.drop_vars('B')`" + "\n `.assign_coords(B=...)` to relabel one side" + "\n or `.isel(..., drop=True)` if a scalar isel " + "introduced it." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_nan_constraint_rhs( + self, x: Variable, time: pd.RangeIndex, unsilenced: None + ) -> None: + nan_rhs = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + msg = _one_legacy_warning(lambda: x <= nan_rhs) + assert msg == ( + "NaN in the constraint RHS was silently kept as 'no " + "constraint at this row' by legacy auto-mask. Under v1 this " + "raises ValueError." + "\n Resolve: `mask=` on the variable for explicit per-row " + "masking" + "\n or `.fillna(value)` if the NaN was a data error." + + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_masked_variable_in_arithmetic( + self, time: pd.RangeIndex, unsilenced: None + ) -> None: + """ + The masked-variable warning is the one that catches the + ``2 * x + y`` (no fillna) divergence β€” no other site fires for + this case. Message names the variable and the fillna(0) fix. + """ + m = Model() + x = m.add_variables(lower=0, coords=[time], name="x") + mask = xr.DataArray( + [True, True, True, True, False], dims=["time"], coords={"time": time} + ) + m.add_variables(lower=0, coords=[time], name="y", mask=mask) + y = m.variables["y"] + msg = _one_legacy_warning(lambda: 2 * x + y) + assert msg == ( + "Variable 'y' has absent slots (from `mask=` / `.where()` / " + "`.shift()` / `.reindex()`). Under legacy each absent slot " + "contributes 0 to the resulting expression's terms (so `x + " + "y >= 10` reduces to `x >= 10` there). Under v1 the absence " + "propagates through arithmetic instead (`x + y` becomes " + "absent at the slot and the constraint drops)." + "\n Resolve: wrap with `y.fillna(0)` for the legacy " + "behaviour under v1" + "\n (no fix needed if you only use the variable " + "in a constraint LHS alone β€” `y >= 0` drops the same way in " + "both)." + _OPT_IN_HINT + ) + + @pytest.mark.legacy + def test_warning_stacklevel_points_to_user_call( + self, time: pd.RangeIndex, unsilenced: None + ) -> None: + """ + The warning's source frame must be the user's call site, not + a linopy internal β€” IDE jump-to-source on the warning depends + on it, and the rollout-warning is useless if it points at + ``linopy/expressions.py`` instead of the user's source. + """ + m = Model() + x = m.add_variables(lower=0, coords=[time], name="x") + mask = xr.DataArray( + [True, True, True, True, False], dims=["time"], coords={"time": time} + ) + y = m.add_variables(lower=0, coords=[time], name="y", mask=mask) + with warnings.catch_warnings(record=True) as ws: + warnings.simplefilter("always") + _ = 2 * x + y # this is the user's call site + relevant = [w for w in ws if issubclass(w.category, LinopySemanticsWarning)] + assert relevant, "expected the masked-variable warning to fire" + assert relevant[0].filename == __file__, ( + f"warning frame is {relevant[0].filename!r}, " + "should be the user's source file" + ) + + +# ===================================================================== +# Β§8 β€” Shared dimensions must match exactly (expr+expr / var+var, merge path) +# ===================================================================== + + +class TestExactAlignmentMerge: + @pytest.fixture + def x_other(self, m: Model) -> Variable: + # Same shape, different labels β€” legacy uses positional override. + return m.add_variables( + lower=0, + coords=[pd.Index([10, 11, 12, 13, 14], name="time")], + name="x_other", + ) + + @pytest.fixture + def x_subset(self, m: Model) -> Variable: + # Subset coords on the same dim β€” legacy outer-joins (and pads). + return m.add_variables( + lower=0, + coords=[pd.Index([1, 3], name="time")], + name="x_subset", + ) + + @pytest.mark.v1 + def test_var_plus_var_different_labels_raises( + self, x: Variable, x_other: Variable + ) -> None: + with pytest.raises(ValueError, match="Coordinate mismatch"): + x + x_other + + @pytest.mark.v1 + def test_expr_plus_expr_different_labels_raises( + self, x: Variable, x_other: Variable + ) -> None: + with pytest.raises(ValueError, match="Coordinate mismatch"): + (1 * x) + (1 * x_other) + + @pytest.mark.v1 + def test_var_plus_var_subset_raises(self, x: Variable, x_subset: Variable) -> None: + with pytest.raises(ValueError, match="Coordinate mismatch"): + x + x_subset + + @pytest.mark.v1 + def test_var_minus_var_different_labels_raises( + self, x: Variable, x_other: Variable + ) -> None: + with pytest.raises(ValueError, match="Coordinate mismatch"): + x - x_other + + def test_var_plus_var_same_coords_works( + self, m: Model, time: pd.RangeIndex + ) -> None: + """Same coords on a shared dim is fine β€” regression guard.""" + a = m.add_variables(lower=0, coords=[time], name="a") + b = m.add_variables(lower=0, coords=[time], name="b") + result = a + b + assert result.sizes["time"] == 5 + + def test_var_plus_var_broadcast_non_shared_dim_works( + self, m: Model, time: pd.RangeIndex + ) -> None: + """Β§9 regression guard for the merge path: non-shared dims broadcast.""" + a = m.add_variables(lower=0, coords=[time], name="a") + b = m.add_variables( + lower=0, coords=[pd.Index([0, 1], name="scenario")], name="b" + ) + result = a + b + assert set(result.coord_dims) == {"time", "scenario"} + + @pytest.mark.v1 + def test_var_plus_var_reordered_labels_align(self, m: Model) -> None: + a = m.add_variables(coords=[pd.Index(["costs", "penalty"], name="e")], name="a") + b = m.add_variables(coords=[pd.Index(["penalty", "costs"], name="e")], name="b") + result = (1 * a) + (1 * b) + assert list(result.indexes["e"]) == ["costs", "penalty"] + # The index alone is identical under legacy β€” assert the variable + # *pairing*: Β§8 sums a["costs"] with b["costs"], NOT b's leading + # (positional) entry b["penalty"]. + costs = {int(v) for v in result.vars.sel(e="costs").values} + assert costs == {int(a.labels.sel(e="costs")), int(b.labels.sel(e="costs"))} + + @pytest.mark.v1 + def test_reordered_constants_pair_by_label_not_position(self, m: Model) -> None: + ea = pd.Index(["costs", "penalty"], name="e") + eb = pd.Index(["penalty", "costs"], name="e") + a = m.add_variables(coords=[ea], name="a") + pd.Series([100.0, 200.0], index=ea) + b = m.add_variables(coords=[eb], name="b") + pd.Series([1.0, 2.0], index=eb) + result = a + b + assert float(result.const.sel(e="costs")) == 102.0 + assert float(result.const.sel(e="penalty")) == 201.0 + + @pytest.mark.v1 + def test_multi_operand_merge_reordered_pairs_by_label(self, m: Model) -> None: + ea = pd.Index(["x", "y", "z"], name="e") + er = pd.Index(["z", "y", "x"], name="e") + a = m.add_variables(coords=[ea], name="a") + pd.Series( + [1.0, 2.0, 3.0], index=ea + ) + b = m.add_variables(coords=[er], name="b") + pd.Series( + [10.0, 20.0, 30.0], index=er + ) + c = m.add_variables(coords=[ea], name="c") + pd.Series( + [100, 200, 300.0], index=ea + ) + result = merge([a, b, c], cls=type(a)) + assert float(result.const.sel(e="x")) == 131.0 + assert float(result.const.sel(e="z")) == 313.0 + + @pytest.mark.v1 + def test_quadratic_merge_reordered_aligns(self, m: Model) -> None: + ea = pd.Index(["x", "y", "z"], name="e") + er = pd.Index(["z", "y", "x"], name="e") + x = m.add_variables(coords=[ea], name="x") + y = m.add_variables(coords=[er], name="y") + result = (x * x) + (y * y) + assert list(result.indexes["e"]) == ["x", "y", "z"] + # by-label: the "x" slot holds x["x"]**2 and y["x"]**2, so the only + # variable labels there are x["x"] and y["x"] β€” not y's leading y["z"]. + labels = {int(v) for v in result.vars.sel(e="x").values.ravel() if v >= 0} + assert labels == {int(x.labels.sel(e="x")), int(y.labels.sel(e="x"))} + + @pytest.mark.v1 + def test_reordered_multiindex_aligns_by_tuple(self, m: Model) -> None: + mi1 = pd.MultiIndex.from_tuples( + [(1, "a"), (1, "b"), (2, "a")], names=["p", "s"] + ) + mi2 = pd.MultiIndex.from_tuples( + [(2, "a"), (1, "b"), (1, "a")], names=["p", "s"] + ) + mi1.name = "snap" + mi2.name = "snap" + x = m.add_variables(coords=[mi1], name="x") + pd.Series( + [1.0, 2.0, 3.0], index=mi1 + ) + y = m.add_variables(coords=[mi2], name="y") + pd.Series( + [10.0, 20.0, 30.0], index=mi2 + ) + result = x + y + got = dict(zip(map(tuple, result.indexes["snap"]), result.const.values)) + assert got == {(1, "a"): 31.0, (1, "b"): 22.0, (2, "a"): 13.0} + + @pytest.mark.legacy + def test_reordered_merge_positional_legacy(self, m: Model) -> None: + ea = pd.Index(["costs", "penalty"], name="e") + eb = pd.Index(["penalty", "costs"], name="e") + a = m.add_variables(coords=[ea], name="a") + pd.Series([100.0, 200.0], index=ea) + b = m.add_variables(coords=[eb], name="b") + pd.Series([1.0, 2.0], index=eb) + result = a + b + assert float(result.const.sel(e="costs")) == 101.0 + assert float(result.const.sel(e="penalty")) == 202.0 + + @pytest.mark.legacy + def test_reordered_merge_warns_legacy(self, m: Model, unsilenced: None) -> None: + ea = pd.Index(["costs", "penalty"], name="e") + eb = pd.Index(["penalty", "costs"], name="e") + a = m.add_variables(coords=[ea], name="a") + b = m.add_variables(coords=[eb], name="b") + with pytest.warns(LinopySemanticsWarning) as record: + (1 * a) + (1 * b) + msg = next( + str(w.message) for w in record if w.category is LinopySemanticsWarning + ) + assert msg == ( + "Coordinate order mismatch in merge along dim '_term' aligned " + "positionally by legacy. Under v1 the same labels in a different " + "order align by label (a reindex), giving a different result." + "\n Dim: 'e': left=['costs', 'penalty'], " + "right=['penalty', 'costs']" + "\n Resolve: `.sel(...)` / `.reindex(...)` to align" + "\n `.assign_coords(...)` to relabel one side" + "\n or pass an explicit `join=` argument." + "\n Opt in: linopy.options['semantics'] = 'v1'" + "\n Silence: warnings.filterwarnings('ignore', " + "category=LinopySemanticsWarning)" + ) + + @pytest.mark.legacy + def test_var_plus_var_reordered_pairs_positionally_legacy(self, m: Model) -> None: + """ + Legacy preserves master's #550 behaviour: var+var with reordered + labels pairs by position, so the "costs" slot wrongly sums + a["costs"] with b's leading entry b["penalty"]. The index still + reads as the left operand's order, which is what hides the + mispairing β€” so assert the variable pairing, not just the index. + """ + a = m.add_variables(coords=[pd.Index(["costs", "penalty"], name="e")], name="a") + b = m.add_variables(coords=[pd.Index(["penalty", "costs"], name="e")], name="b") + result = (1 * a) + (1 * b) + assert list(result.indexes["e"]) == ["costs", "penalty"] + costs = {int(v) for v in result.vars.sel(e="costs").values} + assert costs == {int(a.labels.sel(e="costs")), int(b.labels.sel(e="penalty"))} + + @pytest.mark.legacy + def test_quadratic_merge_reordered_pairs_positionally_legacy( + self, m: Model + ) -> None: + """ + Legacy pairs the quadratic merge by position too: the "x" slot + pairs x["x"]**2 with y's leading y["z"]**2, not y["x"]**2. + """ + ea = pd.Index(["x", "y", "z"], name="e") + er = pd.Index(["z", "y", "x"], name="e") + x = m.add_variables(coords=[ea], name="x") + y = m.add_variables(coords=[er], name="y") + result = (x * x) + (y * y) + assert list(result.indexes["e"]) == ["x", "y", "z"] + labels = {int(v) for v in result.vars.sel(e="x").values.ravel() if v >= 0} + assert labels == {int(x.labels.sel(e="x")), int(y.labels.sel(e="z"))} + + @pytest.mark.v1 + def test_different_multiindex_raises_dim_mismatch(self, m: Model) -> None: + mi1 = pd.MultiIndex.from_tuples( + [(1, "a"), (1, "b"), (2, "a")], names=["p", "s"] + ) + mi2 = pd.MultiIndex.from_tuples( + [(1, "a"), (1, "b"), (3, "c")], names=["p", "s"] + ) + mi1.name = "snap" + mi2.name = "snap" + x = m.add_variables(coords=[mi1], name="x") + y = m.add_variables(coords=[mi2], name="y") + with pytest.raises(ValueError, match="shared dimension 'snap'"): + (1 * x) + (1 * y) + + @pytest.mark.legacy + def test_var_plus_var_different_labels_silent( + self, x: Variable, x_other: Variable + ) -> None: + """ + Document legacy: same-shape var+var aligns by position via + override; the right-hand labels are silently dropped. + """ + result = x + x_other + # Left wins via override β†’ time coords are x's [0..4], even though + # x_other was time=[10..14]. The two terms are paired by position. + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + + @pytest.mark.legacy + def test_warn_on_var_plus_var_different_labels( + self, x: Variable, x_other: Variable, unsilenced: None + ) -> None: + with pytest.warns(LinopySemanticsWarning, match="merge along dim"): + x + x_other + + +# ===================================================================== +# Β§6 β€” Absence propagates through every operator +# Β§3 β€” isnull() reports absent slots (covers #712 absent-as-zero) +# ===================================================================== + + +class TestAbsencePropagation: + @pytest.fixture + def xs(self, x: Variable) -> Variable: + # x.shift(time=1) β†’ absent at time=0, present elsewhere. + return x.shift(time=1) + + @pytest.mark.v1 + def test_to_linexpr_marks_absent_with_nan_const(self, xs: Variable) -> None: + """ + Variable.to_linexpr() encodes absence as NaN const + NaN + coeff + vars=-1, so Β§6 has something to propagate. + """ + expr = xs.to_linexpr() + assert np.isnan(expr.const.values[0]) + assert np.isnan(expr.coeffs.values[0, 0]) + assert int(expr.vars.values[0, 0]) == -1 + assert not np.isnan(expr.const.values[1:]).any() + + @pytest.mark.v1 + def test_isnull_reports_absent_slot(self, xs: Variable) -> None: + """Β§3: isnull() reports the absent slot on a LinearExpression.""" + expr = xs.to_linexpr() + assert bool(expr.isnull().values[0]) + assert not bool(expr.isnull().values[1:].any()) + + @pytest.mark.v1 + @pytest.mark.parametrize("op", ["add", "sub", "mul", "div"]) + def test_scalar_op_preserves_absence(self, xs: Variable, op: str) -> None: + """ + #712 β€” `shifted OP scalar` stays absent at the shifted slot. + Holds for every binary operator: const and coeffs both NaN. + """ + result = _OPS[op](xs, 3) + assert np.isnan(result.const.values[0]) + assert np.isnan(result.coeffs.values[0, 0]) + assert bool(result.isnull().values[0]) + # And the present slots carry the expected per-op value. + expected_const = {"add": 3.0, "sub": -3.0, "mul": 0.0, "div": 0.0}[op] + assert (result.const.values[1:] == expected_const).all() + + @pytest.mark.v1 + def test_add_present_variable_propagates_absence( + self, xs: Variable, x: Variable + ) -> None: + """`x + xs` is absent wherever xs is, even though x is fine there.""" + result = xs + x + assert np.isnan(result.const.values[0]) + assert bool(result.isnull().values[0]) + assert not bool(result.isnull().values[1:].any()) + + @pytest.mark.v1 + def test_merge_absorbs_dead_terms_at_absent_slot( + self, xs: Variable, x: Variable + ) -> None: + """ + Β§1/Β§2 storage invariant β€” ``const.isnull()`` at a slot implies + every term at that slot has ``coeffs = NaN`` and ``vars = -1``. + ``xs + x`` merges xs's absent slot with x's live term; the live + term must be absorbed, not silently kept alongside a NaN const. + Regression guard for ``_absorb_absence`` (commit 4d87a05). + """ + result = xs + x + assert np.isnan(result.coeffs.values[0]).all() + assert (result.vars.values[0] == -1).all() + + @pytest.mark.v1 + def test_merge_absorbs_dead_terms_multi_operand( + self, m: Model, time: pd.RangeIndex + ) -> None: + """ + Same invariant on a 3-operand merge: a regression that absorbs + only on the binary path would still leave one live term at the + absent slot here. + """ + x = m.add_variables(lower=0, coords=[time], name="x") + y = m.add_variables(lower=0, coords=[time], name="y") + xs = x.shift(time=1) + result = (1 * x) + (1 * y) + xs + assert np.isnan(result.coeffs.values[0]).all() + assert (result.vars.values[0] == -1).all() + # And the present rows still carry all three live terms. + assert (~np.isnan(result.coeffs.values[1:])).all() + + @pytest.mark.v1 + def test_absent_distinguishable_from_zero(self, x: Variable, xs: Variable) -> None: + """ + #712 β€” under v1, ``x.shift(time=1) * 3`` and ``x * 0`` are + distinct: the first is absent, the second is a present zero. + """ + absent = xs * 3 + zero = x * 0 + assert bool(absent.isnull().values[0]) + assert not bool(zero.isnull().values[0]) + + @pytest.mark.v1 + @pytest.mark.parametrize( + "build", + [ + "var_mul_var", + "var_pow_2", + "expr_mul_var", + "expr_mul_expr", + "quad_plus_linexpr", + "quad_times_scalar", + ], + ) + def test_quadratic_absence_propagates( + self, xs: Variable, x: Variable, build: str + ) -> None: + """ + Β§6 on the quadratic build paths β€” every entry point that ends + in a QuadraticExpression must keep an absent factor absent. + + Regression for ``prod(skipna=True)`` on the FACTOR_DIM branch + and the cross-term ``self.const * other.reset_const()`` path, + plus the downstream operators on the resulting quadratic. + """ + builders = { + "var_mul_var": lambda: xs * x, + "var_pow_2": lambda: xs**2, + "expr_mul_var": lambda: (1 * xs) * x, + "expr_mul_expr": lambda: (1 * xs) * (1 * x), + "quad_plus_linexpr": lambda: (xs * x) + (2 * x), + "quad_times_scalar": lambda: (xs * x) * 3, + } + quad = builders[build]() + # absent slot stays absent in the resulting quadratic + assert bool(quad.isnull().values[0]) + # and Β§1/Β§2: every term at the absent slot has coeffs NaN and vars -1. + assert np.isnan(quad.coeffs.values[0]).all() + assert (quad.vars.values[0] == -1).all() + # And the present slots stay present (cross-term storage may carry + # vars=-1 as the "no second factor" sentinel inside a term β€” that's + # not absence, so check the slot-level isnull predicate, not vars). + assert not bool(quad.isnull().values[1:].any()) + + @pytest.mark.legacy + def test_legacy_collapses_absent_to_zero(self, xs: Variable) -> None: + """ + Document the #712 bug: legacy treats absent as 0 after `* 3`. + + The term ends up as ``coeffs=3 * vars=-1 + const=0`` β€” a + ``coeff*sentinel`` term that evaluates to 0 at the solver layer. + There is no NaN signal anywhere, so ``isnull()`` returns False and + downstream code can't tell ``xs * 3`` apart from ``x * 0``. + """ + result = xs * 3 + assert not np.isnan(result.const.values[0]) + assert not np.isnan(result.coeffs.values[0, 0]) + assert not bool(result.isnull().values[0]) + + @pytest.mark.legacy + def test_legacy_to_linexpr_fills_absent_with_zero(self, xs: Variable) -> None: + """ + Legacy counterpart of the NaN encoding: ``to_linexpr()`` stores + the absent slot as a present zero (``const=0``, ``coeff=1`` over + the ``vars=-1`` sentinel), so ``isnull()`` is blind to it. + """ + expr = xs.to_linexpr() + assert expr.const.values[0] == 0.0 + assert not np.isnan(expr.coeffs.values[0, 0]) + assert not bool(expr.isnull().values.any()) + + @pytest.mark.legacy + @pytest.mark.parametrize("op", ["add", "sub", "mul", "div"]) + def test_legacy_scalar_op_fills_absent(self, xs: Variable, op: str) -> None: + """ + Legacy fills the absent slot with 0 before the op, so the shifted + slot carries the same value as every present slot (vs v1's NaN). + """ + result = _OPS[op](xs, 3) + assert not bool(result.isnull().values[0]) + expected = {"add": 3.0, "sub": -3.0, "mul": 0.0, "div": 0.0}[op] + assert (result.const.values == expected).all() + + @pytest.mark.legacy + def test_legacy_add_present_variable_keeps_live_term( + self, xs: Variable, x: Variable + ) -> None: + """ + ``xs + x``: legacy keeps ``x[0]`` live at the absent slot (the + merge does not absorb it), so the slot is present, not NaN. + """ + result = xs + x + assert not bool(result.isnull().values[0]) + # x[0] survives as a live term alongside the xs=-1 sentinel. + assert int(x.labels.values[0]) in result.vars.values[0].tolist() + + @pytest.mark.legacy + def test_legacy_merge_keeps_live_terms_multi_operand( + self, m: Model, time: pd.RangeIndex + ) -> None: + """3-operand merge: legacy keeps x[0] and y[0] live at the absent slot.""" + x = m.add_variables(lower=0, coords=[time], name="x") + y = m.add_variables(lower=0, coords=[time], name="y") + xs = x.shift(time=1) + result = (1 * x) + (1 * y) + xs + assert not bool(result.isnull().values[0]) + live = result.vars.values[0][~np.isnan(result.coeffs.values[0])].tolist() + assert int(x.labels.values[0]) in live + assert int(y.labels.values[0]) in live + + @pytest.mark.legacy + @pytest.mark.parametrize( + "build", + [ + "var_mul_var", + "var_pow_2", + "expr_mul_var", + "expr_mul_expr", + "quad_plus_linexpr", + "quad_times_scalar", + ], + ) + def test_legacy_quadratic_collapses_absent( + self, xs: Variable, x: Variable, build: str + ) -> None: + """ + Every quadratic build path collapses the absent factor to a + present zero under legacy β€” no NaN signal anywhere (vs v1, which + keeps the slot absent). + """ + builders = { + "var_mul_var": lambda: xs * x, + "var_pow_2": lambda: xs**2, + "expr_mul_var": lambda: (1 * xs) * x, + "expr_mul_expr": lambda: (1 * xs) * (1 * x), + "quad_plus_linexpr": lambda: (xs * x) + (2 * x), + "quad_times_scalar": lambda: (xs * x) * 3, + } + quad = builders[build]() + assert not bool(quad.isnull().values[0]) + assert not np.isnan(quad.coeffs.values[0]).any() + + +class TestFillnaResolves: + """Β§7 β€” fillna()/.where() are how the caller resolves an absent slot.""" + + @pytest.fixture + def xs(self, x: Variable) -> Variable: + return x.shift(time=1) + + @pytest.mark.v1 + def test_expr_fillna_replaces_absent_const(self, xs: Variable) -> None: + result = xs.to_linexpr().fillna(42) + assert result.const.values[0] == 42.0 + assert result.const.values[1:].tolist() == [0.0, 0.0, 0.0, 0.0] + assert not bool(result.isnull().values.any()) + + @pytest.mark.v1 + def test_variable_fillna_numeric_returns_expression(self, xs: Variable) -> None: + """ + A constant fill is not a variable, so the return type is a + LinearExpression. + """ + from linopy import LinearExpression + + result = xs.fillna(42) + assert isinstance(result, LinearExpression) + assert result.const.values[0] == 42.0 + + @pytest.mark.v1 + def test_variable_fillna_zero_revives_slot_as_present_zero( + self, xs: Variable + ) -> None: + from linopy import LinearExpression + + result = xs.fillna(0) + assert isinstance(result, LinearExpression) # numeric fill β†’ expression + assert not bool(result.isnull().values[0]) + assert result.const.values[0] == 0.0 + + @pytest.mark.v1 + def test_outer_fillna_then_add_collapses_to_just_added( + self, m: Model, time: pd.RangeIndex + ) -> None: + """ + Interpretation A β€” once `(x + y.shift())` is absent at slot 0, + ``.fillna(0)`` revives the slot as the constant 0 (dead terms + stay dead), and a subsequent ``+ x`` re-introduces only ``x[0]``. + Compare ``x + y.shift().fillna(0) + x`` which would double-count + ``x`` at slot 0 β€” the placement of fillna is load-bearing. + """ + x = m.add_variables(lower=0, coords=[time], name="x") + y = m.add_variables(lower=0, coords=[time], name="y") + expr = (x + y.shift(time=1)).fillna(0) + x + + # At slot 0 the only live term is 1Β·x[0]; const is 0 β†’ result == x[0]. + coeffs0 = expr.coeffs.values[0] + vars0 = expr.vars.values[0] + live = ~np.isnan(coeffs0) + assert int(live.sum()) == 1 + assert float(coeffs0[live][0]) == 1.0 + assert int(vars0[live][0]) == int(x.labels.values[0]) + assert float(expr.const.values[0]) == 0.0 + + # At slots 1+ all three terms are live (x[i] + y[i-1] + x[i]) β€” the + # outer ``+ x`` is genuinely additive where y.shift was present. + assert int((~np.isnan(expr.coeffs.values[1])).sum()) == 3 + + @pytest.mark.v1 + def test_masked_variable_constraint_via_fillna(self) -> None: + """ + v1 counterpart of ``test_masked_variable_model`` β€” under Β§6 the + constraint ``x + y >= 10`` drops at the masked y slots, so the + caller must say ``y.fillna(0)`` to keep ``x >= 10`` there. + """ + m = Model() + lower = pd.Series(0, range(10)) + x = m.add_variables(lower, name="x") + mask = pd.Series([True] * 8 + [False, False]) + y = m.add_variables(lower, name="y", mask=mask) + m.add_constraints(x + y.fillna(0), ">=", 10) + m.add_constraints(y, ">=", 0) + m.add_objective(2 * x + y) + + # The constraint x + y.fillna(0) >= 10 binds at every slot. + rhs = m.constraints["con0"].rhs.values + assert not np.isnan(rhs).any() + + @pytest.mark.legacy + def test_legacy_expr_fillna_is_noop(self, xs: Variable) -> None: + """ + Legacy has already filled the absent slot with 0, so there is no + NaN for ``fillna(42)`` to replace β€” the 42 never lands (vs v1, + which puts 42 at the formerly-absent slot). + """ + result = xs.to_linexpr().fillna(42) + assert result.const.values.tolist() == [0.0, 0.0, 0.0, 0.0, 0.0] + assert not bool(result.isnull().values.any()) + + @pytest.mark.legacy + def test_legacy_variable_fillna_numeric_is_noop(self, xs: Variable) -> None: + """Same no-op on the Variable path: the fill value is ignored.""" + from linopy import LinearExpression + + result = xs.fillna(42) + assert isinstance(result, LinearExpression) + assert result.const.values[0] == 0.0 + + @pytest.mark.legacy + def test_legacy_outer_fillna_then_add_double_counts( + self, m: Model, time: pd.RangeIndex + ) -> None: + """ + Legacy never made slot 0 absent, so ``(x + y.shift()).fillna(0)`` + already carries x[0] (and the y sentinel); the outer ``+ x`` then + adds a *second* x[0] β€” three live terms at slot 0, double-counting + x (vs v1's single live x[0]). + """ + x = m.add_variables(lower=0, coords=[time], name="x") + y = m.add_variables(lower=0, coords=[time], name="y") + expr = (x + y.shift(time=1)).fillna(0) + x + + coeffs0 = expr.coeffs.values[0] + vars0 = expr.vars.values[0] + live = ~np.isnan(coeffs0) + assert int(live.sum()) == 3 + # x[0] appears twice β€” the double-count legacy can't avoid. + assert vars0[live].tolist().count(int(x.labels.values[0])) == 2 + assert float(expr.const.values[0]) == 0.0 + + +# ===================================================================== +# Β§4 β€” Variable.reindex / .reindex_like create absence +# ===================================================================== + + +class TestVariableReindex: + """ + Reindexing past the original coords marks the new positions + absent (labels=-1, lower/upper=NaN); Β§4 lists this as one of the + named mechanisms for creating absence. Runs under both semantics: + this is a new API that didn't exist on master. + """ + + def test_reindex_extends_with_absent( + self, x: Variable, time: pd.RangeIndex + ) -> None: + extended = pd.RangeIndex(8, name="time") + result = x.reindex(time=extended) + assert result.sizes["time"] == 8 + # Original slots 0..4 are preserved + assert int(result.labels.values[0]) == int(x.labels.values[0]) + # New slots 5..7 are absent + assert (result.labels.values[5:] == -1).all() + assert np.isnan(result.lower.values[5:]).all() + assert np.isnan(result.upper.values[5:]).all() + + def test_reindex_subset_drops_coords(self, x: Variable) -> None: + """ + Reindex to a strict subset shrinks the variable (no absence + introduced β€” those slots are just gone). + """ + result = x.reindex(time=pd.RangeIndex(3, name="time")) + assert result.sizes["time"] == 3 + assert not (result.labels.values == -1).any() + + def test_reindex_like_extends_with_absent(self, m: Model, x: Variable) -> None: + wider = m.add_variables( + lower=0, coords=[pd.RangeIndex(7, name="time")], name="wider" + ) + result = x.reindex_like(wider) + assert result.sizes["time"] == 7 + assert (result.labels.values[5:] == -1).all() + + @pytest.mark.v1 + def test_reindexed_variable_propagates_absence_in_arithmetic( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """ + Β§4 + Β§6 hand-off: a reindex-introduced absence flows through + the next operator and is visible via isnull(). + """ + wider = x.reindex(time=pd.RangeIndex(7, name="time")) + expr = wider * 3 + assert bool(expr.isnull().values[5:].all()) + assert not bool(expr.isnull().values[:5].any()) + + @pytest.mark.legacy + def test_legacy_reindexed_variable_fills_absent_in_arithmetic( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """ + Legacy collapses the reindex-introduced absence to 0, so the + extended slots are present zeros after ``* 3`` (vs v1, which + keeps them absent and visible via isnull()). + """ + wider = x.reindex(time=pd.RangeIndex(7, name="time")) + expr = wider * 3 + assert not bool(expr.isnull().values.any()) + assert expr.const.values.tolist() == [0.0] * 7 + + def test_where_creates_absence(self, x: Variable) -> None: + """Β§4 β€” ``.where(cond)`` marks slots absent in place.""" + cond = xr.DataArray( + [True, True, False, False, False], + dims=["time"], + coords={"time": x.coords["time"]}, + ) + masked = x.where(cond) + assert (masked.labels.values[2:] == -1).all() + assert not (masked.labels.values[:2] == -1).any() + + def test_unstack_creates_absence_at_missing_combinations(self, m: Model) -> None: + """ + Β§4 β€” ``.unstack`` of a non-rectangular MultiIndex leaves the + missing combinations as absent slots. + """ + # Three (region, year) observations that don't form a full grid: + # (DE, 2030) and (DE, 2040) exist but (FR, 2030) only β€” so + # unstacking (FR, 2040) becomes absent. + idx = pd.MultiIndex.from_tuples( + [("DE", 2030), ("DE", 2040), ("FR", 2030)], + names=("region", "year"), + ) + # Since #732, sequence-form MultiIndex coords entries must be named. + idx.name = "dim_0" + v = m.add_variables(coords=[idx], name="v") + unstacked = v.unstack("dim_0") + assert unstacked.sizes == {"region": 2, "year": 2} + # (FR, 2040) missing β†’ absent + assert int(unstacked.labels.sel(region="FR", year=2040).item()) == -1 + # The three present cells stay present + assert int(unstacked.labels.sel(region="DE", year=2030).item()) != -1 + + @pytest.mark.parametrize("method", ["roll", "sel", "isel"]) + def test_data_preserving_methods_do_not_create_absence( + self, x: Variable, method: str + ) -> None: + """ + Β§4 negative β€” operations that *move or select* existing data + never introduce absent slots. Pins the spec's contrast against + the absence-creating mechanisms. + """ + results = { + "roll": lambda: x.roll(time=2), + "sel": lambda: x.sel(time=[0, 2, 4]), + "isel": lambda: x.isel(time=[0, 2, 4]), + } + result = results[method]() + assert not (result.labels.values == -1).any() + + +# ===================================================================== +# Β§10 β€” named-method join= argument (opt-in alignment) +# ===================================================================== + + +class TestNamedMethodJoin: + """ + Under v1 the bare operators raise on coord mismatch (Β§8). The + named methods let the caller opt in to a specific join mode. + """ + + @pytest.fixture + def subset(self, time: pd.RangeIndex) -> xr.DataArray: + return xr.DataArray( + [10.0, 30.0], dims=["time"], coords={"time": pd.Index([1, 3], name="time")} + ) + + def test_add_join_inner_intersects(self, x: Variable, subset: xr.DataArray) -> None: + """`.add(other, join="inner")` picks the intersection of coords.""" + result = x.add(subset, join="inner") + assert list(result.coords["time"].values) == [1, 3] + + def test_add_join_outer_fills(self, x: Variable, subset: xr.DataArray) -> None: + """`.add(other, join="outer")` unions coords (gaps are filled).""" + result = x.add(subset, join="outer") + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + + def test_mul_join_inner(self, x: Variable, subset: xr.DataArray) -> None: + result = x.mul(subset, join="inner") + assert list(result.coords["time"].values) == [1, 3] + + @pytest.mark.v1 + def test_le_join_inner_on_subset_rhs( + self, x: Variable, subset: xr.DataArray + ) -> None: + """`.le(rhs, join="inner")` lets a subset RHS through cleanly.""" + result = x.le(subset, join="inner") + assert list(result.coords["time"].values) == [1, 3] + + @pytest.mark.legacy + def test_legacy_le_join_inner_keeps_all_coords( + self, x: Variable, subset: xr.DataArray + ) -> None: + """ + On the constraint path legacy ignores ``join="inner"`` and keeps + all left coords, leaving the unmatched RHS slots NaN (vs v1's + clean intersection to [1, 3]). + """ + result = x.le(subset, join="inner") + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + rhs = result.rhs.values + assert np.isnan(rhs[[0, 2, 4]]).all() + assert rhs[[1, 3]].tolist() == [10.0, 30.0] + + @pytest.mark.v1 + def test_bare_op_still_raises_on_mismatch( + self, x: Variable, subset: xr.DataArray + ) -> None: + """`x + subset` (no `join=`) still raises β€” opt-in is required.""" + with pytest.raises(ValueError, match="Coordinate mismatch on shared dimension"): + x + subset + + def test_add_join_override_aligns_positionally(self, x: Variable) -> None: + """ + ``join="override"`` is the explicit-positional mode β€” the right + operand's labels are dropped and the left's are reused. The mode + is opt-in precisely because it can silently mis-pair if the user + didn't mean it. + """ + relabelled = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + result = x.add(relabelled, join="override") + # Override keeps the left operand's labels β€” and silently re-uses + # the right's values at those positions. + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + assert result.const.values.tolist() == [1.0, 2.0, 3.0, 4.0, 5.0] + + def test_add_join_override_size_mismatch_raises(self, x: Variable) -> None: + """ + Β§10 / ``override`` documentation says "positional alignment, made + explicit". Positional pairing is only well-defined when the + shared-dim sizes match; with mismatched sizes ``override`` would + silently mis-pair (or raise opaquely from xarray) instead of + producing a clear error. Regression for the dropped legacy + ``other.sizes == self.const.sizes`` gate. + """ + shorter = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["time"], + coords={"time": pd.Index([0, 1, 2], name="time")}, + ) + with pytest.raises(ValueError, match="join='override' requires matching"): + x.add(shorter, join="override") + + def test_reindex_like_resolves_mismatch_before_bare_op(self, x: Variable) -> None: + """ + Β§10 names ``.reindex(...)`` / ``.reindex_like(...)`` as + canonical resolutions β€” pre-aligning lets the bare operator + accept the once-mismatched operand without ``join=``. + """ + other = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + aligned = other.reindex_like(x.labels, fill_value=0) + result = x + aligned # bare + succeeds because coords now match + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + + def test_assign_coords_resolves_mismatch_before_bare_op(self, x: Variable) -> None: + """ + ``.assign_coords(...)`` is the explicit-positional escape β€” + relabels one side outright so the bare operator's exact-join + check passes. + """ + other = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + relabelled = other.assign_coords(time=x.coords["time"]) + result = x + relabelled # bare + succeeds after relabel + assert list(result.coords["time"].values) == [0, 1, 2, 3, 4] + + +# ===================================================================== +# Β§12 β€” constraints follow the same rules +# ===================================================================== + + +_SIGNS = { + "le": operator.le, + "ge": operator.ge, + "eq": operator.eq, +} + + +class TestConstraintRHS: + @pytest.mark.v1 + @pytest.mark.parametrize("sign", ["le", "ge", "eq"]) + def test_subset_rhs_raises(self, x: Variable, sign: str) -> None: + """Β§12 β€” all three comparison signs align by Β§8 the same way.""" + subset = xr.DataArray( + [10.0, 20.0], + dims=["time"], + coords={"time": pd.Index([1, 3], name="time")}, + ) + with pytest.raises(ValueError, match="Coordinate mismatch on shared dimension"): + _SIGNS[sign](x, subset) + + @pytest.mark.v1 + @pytest.mark.parametrize("sign", ["le", "ge", "eq"]) + def test_nan_rhs_raises(self, x: Variable, time: pd.RangeIndex, sign: str) -> None: + """ + Β§5/Β§12 β€” a NaN in a user-supplied RHS raises for every sign, + never silently becomes "no constraint" the way legacy auto_mask + treats it. + """ + nan_rhs = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + with pytest.raises(ValueError, match="NaN"): + _SIGNS[sign](x, nan_rhs) + + @pytest.mark.v1 + @pytest.mark.parametrize("sign", ["le", "ge", "eq"]) + def test_absence_propagates_to_rhs_drops_constraint( + self, x: Variable, sign: str + ) -> None: + """ + Β§6 β†’ Β§12 for every sign: a constraint over an absent LHS slot + yields NaN RHS, which downstream auto-mask interprets as "no + constraint here". + """ + xs = x.shift(time=1) + # xs is absent at time=0; the constraint's RHS at that slot + # should be NaN (no constraint), not 10. + constraint = _SIGNS[sign](xs, 10) + rhs = constraint.rhs.values + assert np.isnan(rhs[0]) + assert (rhs[1:] == 10).all() + + @pytest.mark.legacy + @pytest.mark.parametrize("sign", ["le", "ge", "eq"]) + def test_legacy_absence_keeps_rhs_at_absent_slot( + self, x: Variable, sign: str + ) -> None: + """ + Legacy fills the absent LHS slot with 0, so the RHS stays 10 + everywhere and the constraint is emitted at that slot too (vs + v1, where the NaN RHS drops it). + """ + xs = x.shift(time=1) + constraint = _SIGNS[sign](xs, 10) + assert (constraint.rhs.values == 10).all() + + @pytest.mark.v1 + def test_pypsa_1683_nan_rhs_raises(self, x: Variable, time: pd.RangeIndex) -> None: + """ + PyPSA #1683 on the constraint side β€” ``min_pu * nominal_fix`` + with ``p_nom=inf`` and ``p_min_pu=0`` yields NaN at the bad slot; + v1 raises at construction instead of silently passing NaN to + the solver. + """ + nominal = xr.DataArray([np.inf] * 5, dims=["time"], coords={"time": time}) + min_pu = xr.DataArray( + [1.0, 0.0, 1.0, 1.0, 1.0], dims=["time"], coords={"time": time} + ) + bound = min_pu * nominal # 0*inf = NaN at time=1 + with pytest.raises(ValueError, match="NaN"): + x >= bound + + @pytest.mark.legacy + def test_nan_rhs_silently_treated_as_unconstrained( + self, x: Variable, time: pd.RangeIndex + ) -> None: + """ + Document the legacy auto_mask path: a NaN RHS is silently + kept as NaN and the constraint at that row is later dropped. + """ + nan_rhs = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + constraint = x <= nan_rhs + assert np.isnan(constraint.rhs.values[1]) + + @pytest.mark.legacy + def test_warn_on_nan_rhs( + self, x: Variable, time: pd.RangeIndex, unsilenced: None + ) -> None: + nan_rhs = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + with pytest.warns(LinopySemanticsWarning, match="no constraint at this row"): + x <= nan_rhs + + @pytest.mark.legacy + def test_warn_on_coord_mismatch_rhs_distinguishes_from_nan( + self, x: Variable, unsilenced: None + ) -> None: + """ + A subset RHS has no user NaN β€” legacy's ``reindex_like`` is what + introduces the NaN at the unmatched positions. The warning should + diagnose the *coord mismatch* (fix: ``.sel`` / ``.reindex``), not + the NaN-RHS auto-mask (fix: ``mask=`` / ``.fillna``). Regression + for the conflated warn text where both causes used the same + ``_legacy_nan_rhs_constraint_message``. + """ + subset = xr.DataArray( + [10.0, 20.0], + dims=["time"], + coords={"time": pd.Index([1, 3], name="time")}, + ) + with pytest.warns( + LinopySemanticsWarning, match="Coordinate mismatch in constraint RHS" + ): + x <= subset + + @pytest.mark.legacy + def test_both_warnings_fire_when_rhs_has_user_nan_and_mismatch( + self, x: Variable, unsilenced: None + ) -> None: + """ + Independent causes β€” when the RHS is both subset (mismatch) and + carries a user NaN, both fix-hints should surface so the caller + sees each problem with its own resolution. + """ + both = xr.DataArray( + [10.0, np.nan], + dims=["time"], + coords={"time": pd.Index([1, 3], name="time")}, + ) + with warnings.catch_warnings(record=True) as caught: + warnings.simplefilter("always", LinopySemanticsWarning) + x <= both + messages = [str(w.message) for w in caught] + assert any("Coordinate mismatch in constraint RHS" in m for m in messages) + assert any("no constraint at this row" in m for m in messages) + + +# ===================================================================== +# Β§13 β€” reductions skip absent slots (not propagate) +# ===================================================================== + + +class TestReductionsSkipAbsent: + """ + Per Β§13, ``sum`` / ``groupby.sum`` skip absent slots rather than + propagating them β€” the only asymmetry against Β§6's binary-operator + rule. The expected behaviour falls out of xarray's ``skipna=True`` + default; these tests pin it under v1 so future changes don't drift. + + Each v1 test is paired with its ``legacy`` counterpart so the + divergence is pinned on both sides: legacy fills the absent slot + with 0, so ``xs + 5`` stores ``5`` there too and that extra term is + counted by the sum (β†’ 25), whereas v1 skips it (β†’ 20). + + Scope: Β§13 also names ``mean``, ``resample``, and ``coarsen``, but + those are not yet exposed on ``LinearExpression`` (see #703). The + spec text is the rule they will follow when implemented; tests + belong with the implementation PR. + """ + + @pytest.fixture + def xs(self, x: Variable) -> Variable: + return x.shift(time=1) + + @pytest.mark.v1 + def test_sum_over_dim_skips_absent(self, xs: Variable) -> None: + """ + ``(xs + 5).sum('time')`` skips the absent slot at t=0 and + sums the four present 5s β†’ 20. + """ + result = (xs + 5).sum("time") + assert float(result.const) == 20.0 + + @pytest.mark.v1 + def test_sum_no_dim_skips_absent(self, xs: Variable) -> None: + result = (xs + 5).sum() + assert float(result.const) == 20.0 + + @pytest.mark.v1 + def test_sum_of_all_absent_is_zero(self, x: Variable) -> None: + """Β§13 β€” "the sum of none is the zero expression.""" "" + all_absent = x.shift(time=10).to_linexpr() + assert bool(all_absent.isnull().all().item()) + result = all_absent.sum("time") + assert float(result.const) == 0.0 + + @pytest.mark.v1 + def test_groupby_sum_skips_absent(self, xs: Variable) -> None: + """Each group's sum drops absent members, just like ``.sum``.""" + groups = xr.DataArray( + [0, 0, 1, 1, 1], dims=["time"], coords={"time": xs.coords["time"]} + ) + result = (xs + 5).groupby(groups).sum() + # group 0: [NaN, 5] β†’ 5; group 1: [5, 5, 5] β†’ 15 + assert result.const.values.tolist() == [5.0, 15.0] + + @pytest.mark.legacy + def test_sum_over_dim_fills_absent(self, xs: Variable) -> None: + """ + Legacy fills the absent slot at t=0 with 0, so ``xs + 5`` stores + ``5`` there and all five 5s are summed β†’ 25 (vs v1's 20). + """ + result = (xs + 5).sum("time") + assert float(result.const) == 25.0 + + @pytest.mark.legacy + def test_sum_no_dim_fills_absent(self, xs: Variable) -> None: + result = (xs + 5).sum() + assert float(result.const) == 25.0 + + @pytest.mark.legacy + def test_groupby_sum_fills_absent(self, xs: Variable) -> None: + """Legacy fills the absent member with 0, so group 0 sums [5, 5] β†’ 10.""" + groups = xr.DataArray( + [0, 0, 1, 1, 1], dims=["time"], coords={"time": xs.coords["time"]} + ) + result = (xs + 5).groupby(groups).sum() + # group 0: [5, 5] β†’ 10 (vs v1's 5); group 1: [5, 5, 5] β†’ 15 + assert result.const.values.tolist() == [10.0, 15.0] + + +# ===================================================================== +# Β§11 β€” auxiliary (non-dim) coordinate conflicts raise (covers #295) +# ===================================================================== + + +class TestAuxCoordConflict: + """ + Per Β§11, an auxiliary (non-dim) coord that two operands carry + with disagreeing values must raise β€” xarray silently drops the + conflict in arithmetic, which is the #295 bug. + """ + + @pytest.fixture + def A(self) -> pd.Index: + return pd.Index([1, 2, 3], name="A") + + @pytest.mark.v1 + def test_expr_plus_dataarray_aux_conflict_raises( + self, m: Model, A: pd.Index + ) -> None: + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v + const + + @pytest.mark.v1 + def test_var_plus_var_aux_conflict_raises(self, m: Model, A: pd.Index) -> None: + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + w = m.add_variables(lower=0, coords=[A], name="w").assign_coords( + B=("A", [400, 400, 500]) + ) + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v + w + + @pytest.mark.v1 + def test_aux_conflict_survives_reordered_dim(self, m: Model) -> None: + v = m.add_variables( + lower=0, coords=[pd.Index(["x", "y", "z"], name="A")], name="v" + ).assign_coords(B=("A", [1, 2, 3])) + w = m.add_variables( + lower=0, coords=[pd.Index(["z", "y", "x"], name="A")], name="w" + ).assign_coords(B=("A", [1, 2, 3])) + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v + w + + @pytest.mark.v1 + def test_mul_constant_aux_conflict_raises(self, m: Model, A: pd.Index) -> None: + """Same rule on the multiplication path β€” not just ``+``.""" + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [2.0, 3.0, 4.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v * const + + @pytest.mark.v1 + def test_constraint_aux_conflict_raises(self, m: Model, A: pd.Index) -> None: + """Β§11 reaches constraint construction via the same machinery.""" + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v == const + + @pytest.mark.v1 + def test_scalar_isel_aux_conflict_raises(self, m: Model, A: pd.Index) -> None: + """ + Scalar isels leave the indexed dim as a non-dim coord whose + value differs between operands picked at different positions. + """ + v = m.add_variables(lower=0, coords=[A], name="v") + a0 = (1 * v).isel({"A": 0}) # scalar A=1 + a1 = (1 * v).isel({"A": 1}) # scalar A=2 + with pytest.raises(ValueError, match="Auxiliary coordinate"): + a0 + a1 + + def test_isel_with_drop_true_avoids_conflict(self, m: Model, A: pd.Index) -> None: + """ + The Β§11 escape hatch the convention recommends: drop the + leftover scalar coord with ``isel(..., drop=True)``. + """ + v = m.add_variables(lower=0, coords=[A], name="v") + a0 = (1 * v).isel({"A": 0}, drop=True) + a1 = (1 * v).isel({"A": 1}, drop=True) + result = a0 + a1 # no aux coord β†’ no conflict + assert "A" not in result.coords + + def test_assign_coords_resolves_conflict(self, m: Model, A: pd.Index) -> None: + """ + Β§11's third escape hatch: relabel one side with + ``.assign_coords`` so the coord values agree across operands. + """ + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + relabelled = const.assign_coords(B=v.coords["B"]) + result = v + relabelled + assert result.coords["B"].values.tolist() == [311, 311, 322] + + @pytest.mark.v1 + def test_multi_operand_merge_aux_conflict_raises( + self, m: Model, A: pd.Index + ) -> None: + """ + The merge-path check inspects all operands, not just two β€” + a 3-way ``sum(...)`` where the third disagrees still raises. + """ + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + w = m.add_variables(lower=0, coords=[A], name="w").assign_coords( + B=("A", [311, 311, 322]) + ) + u = m.add_variables(lower=0, coords=[A], name="u").assign_coords( + B=("A", [999, 999, 999]) + ) + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v + w + u + + @pytest.mark.v1 + def test_aux_conflict_raises_under_explicit_join_constant( + self, m: Model, A: pd.Index + ) -> None: + """ + Β§11 is independent of Β§8 β€” an explicit ``join=`` must not + silence the aux-coord raise. Regression for the + ``if join is None:`` gating bug where ``.add(const, join="override")`` + would silently drop the conflicting coord. + """ + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + for join in ("override", "inner", "outer", "left", "right", "exact"): + with pytest.raises(ValueError, match="Auxiliary coordinate"): + v.add(const, join=join) + + @pytest.mark.v1 + def test_aux_conflict_raises_under_explicit_join_merge( + self, m: Model, A: pd.Index + ) -> None: + """ + Same rule on the merge path: ``linopy.merge([v, w], join="override")`` + with a conflicting aux coord must raise. + """ + import linopy + + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + w = m.add_variables(lower=0, coords=[A], name="w").assign_coords( + B=("A", [400, 400, 500]) + ) + for join in ("override", "inner", "outer", "left", "right", "exact"): + with pytest.raises(ValueError, match="Auxiliary coordinate"): + linopy.merge([1 * v, 1 * w], join=join) + + @pytest.mark.legacy + def test_aux_conflict_silently_keeps_left(self, m: Model, A: pd.Index) -> None: + """ + Document legacy: a conflict is silently resolved by keeping + the left operand's aux coord β€” the right operand's [400,400,500] + disappears with no signal to the caller. + """ + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + result = v + const + assert result.coords["B"].values.tolist() == [311, 311, 322] + + @pytest.mark.legacy + def test_warn_on_aux_conflict( + self, m: Model, A: pd.Index, unsilenced: None + ) -> None: + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + with pytest.warns(LinopySemanticsWarning, match=r"(?s)'B'.*silently dropped"): + v + const + + +class TestAuxCoordPropagation: + """ + Non-conflicting aux coords must propagate through arithmetic and + into constraints β€” the positive half of Β§11. + """ + + @pytest.fixture + def A(self) -> pd.Index: + return pd.Index([1, 2, 3], name="A") + + def test_aux_coord_survives_scalar_mul(self, m: Model, A: pd.Index) -> None: + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + assert "B" in (3 * v).coords + + def test_aux_coord_survives_scalar_add(self, m: Model, A: pd.Index) -> None: + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + assert "B" in (v + 5).coords + + def test_aux_coord_propagates_through_var_plus_var( + self, m: Model, A: pd.Index + ) -> None: + B = ("A", [311, 311, 322]) + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords(B=B) + w = m.add_variables(lower=0, coords=[A], name="w").assign_coords(B=B) + result = v + w + assert "B" in result.coords + assert result.coords["B"].values.tolist() == [311, 311, 322] + + def test_aux_coord_propagates_into_constraint(self, m: Model, A: pd.Index) -> None: + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + c = v <= 10 + assert "B" in c.coords + + def test_aux_coord_only_on_dataarray_propagates( + self, m: Model, A: pd.Index + ) -> None: + """ + ``x * a`` where ``a`` carries an aux coord and ``x`` doesn't β€” + the coord propagates through every binary operator and into the + constraint. Hits the `_align_constant` path (var-OP-DataArray) + distinct from the `merge` path tested below. + """ + x = m.add_variables(lower=0, coords=[A], name="x") + a = xr.DataArray( + [2.0, 3.0, 4.0], dims=["A"], coords={"A": A, "B": ("A", [10, 20, 30])} + ) + for expr in (x * a, x + a, x / a): + assert "B" in expr.coords + assert expr.coords["B"].values.tolist() == [10, 20, 30] + # And into the constraint + c = x <= a + assert "B" in c.coords + + def test_aux_coord_only_on_one_side_propagates(self, m: Model, A: pd.Index) -> None: + """Var+var counterpart of the above β€” hits the `merge` path.""" + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + w = m.add_variables(lower=0, coords=[A], name="w") # no B + result = v + w + assert "B" in result.coords + + def test_aux_coord_object_dtype_with_nan_compares_equal( + self, m: Model, A: pd.Index + ) -> None: + """ + Aux coords with object dtype can embed NaN placeholders (e.g. + ragged category labels). Two operands with identical NaN + placement must compare equal β€” `np.array_equal` alone treats + NaN as self-unequal on object dtype, so the Β§11 raise would + false-positive without the pandas-equals fallback. + """ + B = np.array([311, np.nan, 322], dtype=object) + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords(B=("A", B)) + w = m.add_variables(lower=0, coords=[A], name="w").assign_coords(B=("A", B)) + # Same B on both sides, NaN at the same slot β€” should propagate, not raise. + result = v + w + assert "B" in result.coords + + +# ===================================================================== +# Error-message content (raise self-description) +# ===================================================================== + + +class TestErrorMessageContent: + """ + The three v1 raises must be self-describing: name the dim or + coord and show the disagreeing values so the user can act on the + message without re-running with extra prints. Substring assertions + elsewhere don't cover this β€” these tests pin the rich content. + """ + + @pytest.fixture + def A(self) -> pd.Index: + return pd.Index([1, 2, 3], name="A") + + @pytest.mark.v1 + def test_user_nan_message_separates_intents(self, x: Variable) -> None: + """ + The Β§5 raise must not collapse `data error` and `absence` into a + single suggestion β€” they need different fixes. + """ + with pytest.raises(ValueError) as exc: + x + float("nan") + msg = str(exc.value) + assert "data error" in msg and ".fillna(value)" in msg + assert "absent" in msg and "mask=" in msg + assert ".reindex" in msg or ".where(cond)" in msg + + @pytest.mark.v1 + def test_shared_dim_message_names_dim_and_values( + self, m: Model, time: pd.RangeIndex + ) -> None: + """ + The merge-path Β§8 raise must name the offending dim and show + both sides' labels β€” otherwise the user can't tell which dim + out of many. + """ + other = m.add_variables( + lower=0, coords=[pd.Index([10, 11, 12, 13, 14], name="time")], name="other" + ) + x_local = m.add_variables(lower=0, coords=[time], name="x_local") + with pytest.raises(ValueError) as exc: + x_local + other + msg = str(exc.value) + assert "'time'" in msg + assert "[0, 1, 2, 3, 4]" in msg + assert "[10, 11, 12, 13, 14]" in msg + + @pytest.mark.v1 + def test_aux_conflict_message_names_coord_and_values( + self, m: Model, A: pd.Index + ) -> None: + """ + The Β§11 raise must name the conflicting coord and show both + sides' values β€” and mention `.assign_coords` as a fix, not only + `.drop_vars` and `isel(drop=True)`. + """ + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], + dims=["A"], + coords={"A": A, "B": ("A", [400, 400, 500])}, + ) + with pytest.raises(ValueError) as exc: + v + const + msg = str(exc.value) + assert "'B'" in msg + assert "conflicting values" in msg # value-mismatch failure mode + assert "[311, 311, 322]" in msg + assert "[400, 400, 500]" in msg + # All three resolution paths from Β§11 should be listed. + assert ".drop_vars" in msg + assert ".assign_coords" in msg + assert ".isel" in msg + + @pytest.mark.v1 + def test_aux_conflict_message_distinguishes_shape_vs_value( + self, m: Model, A: pd.Index + ) -> None: + """ + Shape mismatch and value disagreement are different failure + modes β€” surface that in the message text so the caller can + diagnose without re-reading both arrays. + """ + # scalar-isel leaves a 0-d aux coord on one side; the full vector + # on the other has a different shape, not a different value. + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + scalar_side = (1 * v).isel({"A": 0}) # B becomes a 0-d scalar coord + full_side = 1 * v + with pytest.raises(ValueError, match="differing shapes") as exc: + scalar_side + full_side + msg = str(exc.value) + assert "'B'" in msg + assert "shape" in msg + + +# ===================================================================== +# Rough edges β€” catches NaN that slips past the operator-level check +# ===================================================================== + + +class TestUserNaNEdgeCases: + """ + Regression guards for three NaN-entry routes that were untested. + The first two are already caught upstream (at the operator that + constructs the expression); the third needed an ``add_objective`` + boundary check because a hand-built expression with NaN const + skips the operator path entirely. + """ + + @pytest.mark.v1 + def test_nan_in_expression_used_in_objective_raises( + self, m: Model, time: pd.RangeIndex + ) -> None: + """ + ``add_objective((x * nan_costs).sum())`` raises at the ``*`` + before the objective even sees the expression β€” guard against + a regression that lets NaN-cost objectives slip through. + """ + x = m.add_variables(lower=0, coords=[time], name="x") + nan_costs = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + with pytest.raises(ValueError, match="NaN"): + m.add_objective((x * nan_costs).sum()) + + @pytest.mark.v1 + def test_nan_in_constraint_lhs_raises(self, m: Model, time: pd.RangeIndex) -> None: + """ + ``(x + nan_da) <= 5`` raises at the ``+`` on the LHS β€” the + RHS path is tested elsewhere; this guards the symmetric LHS + case. + """ + x = m.add_variables(lower=0, coords=[time], name="x") + nan_da = xr.DataArray( + [1.0, np.nan, 3.0, 4.0, 5.0], dims=["time"], coords={"time": time} + ) + with pytest.raises(ValueError, match="NaN"): + (x + nan_da) <= 5 + + +# ===================================================================== +# Object scope β€” non-linopy operands behave like constant expressions +# ===================================================================== + + +class TestObjectScope: + """ + Per the convention's object-scope statement, behaviour is + object-agnostic: ``x OP arr`` builds exactly what ``x OP arr_expr`` + builds, where ``arr_expr`` is the constant-only LinearExpression + holding ``arr``'s values and coordinates β€” whatever type ``arr`` + enters as, in either operand position. + """ + + _VALUES = [1.0, 2.0, 3.0, 4.0, 5.0] + + @pytest.fixture + def da(self, time: pd.RangeIndex) -> xr.DataArray: + return xr.DataArray(self._VALUES, dims=["time"], coords={"time": time}) + + def raw_and_wrapped( + self, kind: str, m: Model, time: pd.RangeIndex, da: xr.DataArray + ) -> tuple[Any, Any]: + """Return a raw constant of the given kind and its constant-expression twin.""" + from linopy import LinearExpression + + if kind == "dataarray": + return da, LinearExpression(da, m) + if kind == "series": + return pd.Series(self._VALUES, index=time), LinearExpression(da, m) + if kind == "scalar": + return 7.5, LinearExpression(7.5, m) + raise AssertionError(kind) + + @pytest.mark.parametrize("op", ["add", "sub", "mul", "radd", "rsub", "rmul"]) + @pytest.mark.parametrize("kind", ["dataarray", "series", "scalar"]) + def test_op_matches_const_expr_op( + self, + m: Model, + x: Variable, + time: pd.RangeIndex, + da: xr.DataArray, + kind: str, + op: str, + ) -> None: + raw, wrapped = self.raw_and_wrapped(kind, m, time, da) + forward = op in ("add", "sub", "mul") + opfunc = _OPS[op.removeprefix("r")] + if forward: + assert_linequal(opfunc(x, raw), opfunc(x, wrapped)) + else: + assert_linequal(opfunc(raw, x), opfunc(wrapped, x)) + + @pytest.mark.parametrize("kind", ["dataarray", "series", "scalar"]) + def test_distributive_law_mixed_types( + self, + m: Model, + x: Variable, + time: pd.RangeIndex, + da: xr.DataArray, + kind: str, + ) -> None: + """``(x + y) * arr`` distributes the same whether ``arr`` is raw or wrapped.""" + raw, wrapped = self.raw_and_wrapped(kind, m, time, da) + y = m.add_variables(lower=0, coords=[time], name="y") + assert_linequal((x + y) * raw, x * raw + y * raw) + assert_linequal((x + y) * raw, (x + y) * wrapped) + + @pytest.mark.parametrize("kind", ["dataarray", "series", "scalar"]) + def test_associative_law_mixed_types( + self, + m: Model, + x: Variable, + time: pd.RangeIndex, + da: xr.DataArray, + kind: str, + ) -> None: + """``(x + arr) + y`` and ``x + (arr + y)`` agree for raw constants.""" + raw, _ = self.raw_and_wrapped(kind, m, time, da) + y = m.add_variables(lower=0, coords=[time], name="y") + assert_linequal((x + raw) + y, x + (raw + y)) + + @pytest.mark.v1 + def test_coord_mismatch_raises_on_either_route(self, m: Model, x: Variable) -> None: + """Β§8 fires identically whether the mismatched constant is raw or wrapped.""" + from linopy import LinearExpression + + mismatched = xr.DataArray( + self._VALUES, + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + with pytest.raises(ValueError, match="Coordinate mismatch on shared dimension"): + x + mismatched + with pytest.raises(ValueError, match="Coordinate mismatch on shared dimension"): + x + LinearExpression(mismatched, m) + + def test_division_by_const_expr_is_type_error( + self, m: Model, x: Variable, da: xr.DataArray + ) -> None: + """ + The one type-decided footnote: a constant can be a divisor, an + expression cannot β€” even one holding only constants. + """ + from linopy import LinearExpression + + x / da # works: dividing by a constant + with pytest.raises(TypeError): + x / LinearExpression(da, m) + + +# ===================================================================== +# Cross-cutting guard: operations that MUST agree under both semantics +# ===================================================================== +# +# The autouse ``semantics`` fixture runs each test under a single mode, so +# a per-op test that under-asserts (e.g. checks only ``.indexes``) can pass +# under both modes while the actual result silently diverges β€” that is how +# the reordered-merge mispairing slipped through review. This guard builds +# each mode-invariant operation under *both* semantics and compares them +# with linopy.testing's strict structural helpers (which align by coords +# and compare vars/coeffs/const), so a regression that makes one of these +# paths semantics-dependent fails loudly. Genuinely divergent operations +# belong in the per-section classes above, not here. + + +def _build_under_both(builder: Any) -> tuple[Any, Any]: + """Build ``builder()`` under legacy then v1; return ``(legacy, v1)``.""" + from linopy.config import options + + saved = options["semantics"] + out = {} + try: + for sem in ("legacy", "v1"): + options["semantics"] = sem + with warnings.catch_warnings(): + warnings.simplefilter("ignore", LinopySemanticsWarning) + out[sem] = builder() + finally: + options["semantics"] = saved + return out["legacy"], out["v1"] + + +def _eq_time() -> pd.RangeIndex: + return pd.RangeIndex(5, name="time") + + +def _eq_da() -> xr.DataArray: + return xr.DataArray(np.arange(1.0, 6.0), dims=["time"], coords={"time": _eq_time()}) + + +def _eq_subset() -> xr.DataArray: + return xr.DataArray( + [10.0, 30.0], dims=["time"], coords={"time": pd.Index([1, 3], name="time")} + ) + + +def _op_merge_same_coords() -> LinearExpression: + m = Model() + return m.add_variables(lower=0, coords=[_eq_time()], name="a") + m.add_variables( + lower=0, coords=[_eq_time()], name="b" + ) + + +def _op_merge_broadcast() -> LinearExpression: + m = Model() + a = m.add_variables(lower=0, coords=[_eq_time()], name="a") + b = m.add_variables(lower=0, coords=[pd.Index([0, 1], name="scenario")], name="b") + return a + b + + +def _op_quadratic_same_coords() -> QuadraticExpression: + m = Model() + x = m.add_variables(coords=[_eq_time()], name="x") + y = m.add_variables(coords=[_eq_time()], name="y") + return cast(QuadraticExpression, (x * x) + (y * y)) + + +def _op_add_join_inner() -> LinearExpression: + return cast( + LinearExpression, + Model() + .add_variables(coords=[_eq_time()], name="x") + .add(_eq_subset(), join="inner"), + ) + + +def _op_add_join_outer() -> LinearExpression: + return cast( + LinearExpression, + Model() + .add_variables(coords=[_eq_time()], name="x") + .add(_eq_subset(), join="outer"), + ) + + +def _op_add_join_override() -> LinearExpression: + relabelled = xr.DataArray( + [1.0, 2.0, 3.0, 4.0, 5.0], + dims=["time"], + coords={"time": pd.Index([10, 11, 12, 13, 14], name="time")}, + ) + return cast( + LinearExpression, + Model() + .add_variables(coords=[_eq_time()], name="x") + .add(relabelled, join="override"), + ) + + +def _op_associative() -> LinearExpression: + m = Model() + x = m.add_variables(coords=[_eq_time()], name="x") + y = m.add_variables(lower=0, coords=[_eq_time()], name="y") + return (x + _eq_da()) + y + + +def _op_distributive() -> LinearExpression: + m = Model() + x = m.add_variables(coords=[_eq_time()], name="x") + y = m.add_variables(lower=0, coords=[_eq_time()], name="y") + return (x + y) * _eq_da() + + +def _op_raw_matches_wrapped() -> LinearExpression: + m = Model() + x = m.add_variables(coords=[_eq_time()], name="x") + return x + LinearExpression(_eq_da(), m) + + +def _op_aux_coord_resolved() -> LinearExpression: + m = Model() + A = pd.Index([1, 2, 3], name="A") + v = m.add_variables(lower=0, coords=[A], name="v").assign_coords( + B=("A", [311, 311, 322]) + ) + const = xr.DataArray( + [10.0, 20.0, 30.0], dims=["A"], coords={"A": A, "B": ("A", [400, 400, 500])} + ) + return v + const.assign_coords(B=v.coords["B"]) + + +def _op_clean_constraint() -> Any: + m = Model() + x = m.add_variables(coords=[_eq_time()], name="x") + return x <= _eq_da() + + +_EQUIVALENT_OPS = [ + pytest.param(_op_merge_same_coords, assert_linequal, id="merge_same_coords"), + pytest.param(_op_merge_broadcast, assert_linequal, id="merge_broadcast"), + pytest.param(_op_quadratic_same_coords, assert_quadequal, id="quadratic_same"), + pytest.param(_op_add_join_inner, assert_linequal, id="add_join_inner"), + pytest.param(_op_add_join_outer, assert_linequal, id="add_join_outer"), + pytest.param(_op_add_join_override, assert_linequal, id="add_join_override"), + pytest.param(_op_associative, assert_linequal, id="associative_law"), + pytest.param(_op_distributive, assert_linequal, id="distributive_law"), + pytest.param(_op_raw_matches_wrapped, assert_linequal, id="raw_matches_wrapped"), + pytest.param(_op_aux_coord_resolved, assert_linequal, id="aux_coord_resolved"), + pytest.param(_op_clean_constraint, assert_conequal, id="clean_constraint"), +] + + +@pytest.mark.parametrize("builder, comparator", _EQUIVALENT_OPS) +def test_semantics_invariant_ops_agree(builder: Any, comparator: Any) -> None: + """Mode-invariant operations must be byte-identical under both semantics.""" + legacy, v1 = _build_under_both(builder) + comparator(legacy, v1) diff --git a/test/test_linear_expression.py b/test/test_linear_expression.py index 5ffd7de1..1f135e86 100644 --- a/test/test_linear_expression.py +++ b/test/test_linear_expression.py @@ -19,8 +19,8 @@ from xarray.testing import assert_equal from linopy import ( - EvolvingAPIWarning, LinearExpression, + LinopySemanticsWarning, Model, QuadraticExpression, Variable, @@ -296,6 +296,7 @@ def test_linear_expression_multi_indexed(u: Variable) -> None: assert isinstance(expr, LinearExpression) +@pytest.mark.legacy def test_multiply_expression_by_multiindex_level_constant(u: Variable) -> None: """ Expression over a MultiIndex dim times a single-level constant. @@ -307,7 +308,7 @@ def test_multiply_expression_by_multiindex_level_constant(u: Variable) -> None: """ by_level1 = xr.DataArray([10.0, 20.0], coords={"level1": [1, 2]}, dims=["level1"]) - with pytest.warns(EvolvingAPIWarning, match=r"broadcasting level subset"): + with pytest.warns(LinopySemanticsWarning, match=r"broadcasting level subset"): expr = (1 * u) * by_level1 coeffs = expr.coeffs.squeeze("_term") @@ -317,6 +318,15 @@ def test_multiply_expression_by_multiindex_level_constant(u: Variable) -> None: assert coeffs.sel(dim_3=(2, "b")).item() == 20.0 +@pytest.mark.v1 +def test_multiply_expression_by_mi_level_constant_raises_v1(u: Variable) -> None: + """v1: the implicit level projection in arithmetic must be explicit.""" + by_level1 = xr.DataArray([10.0, 20.0], coords={"level1": [1, 2]}, dims=["level1"]) + + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + (1 * u) * by_level1 + + def test_linear_expression_with_errors(m: Model, x: Variable) -> None: with pytest.raises(TypeError): x / x @@ -417,9 +427,12 @@ def test_linear_expression_substraction( assert res.data.notnull().all().to_array().all() +@pytest.mark.legacy def test_linear_expression_sum( x: Variable, y: Variable, z: Variable, v: Variable ) -> None: + # Legacy-only: ``v.loc[:9] + v.loc[10:]`` merges disjoint coords + # (forbidden by v1 Β§8). expr = 10 * x + y + z res = expr.sum("dim_0") @@ -439,9 +452,12 @@ def test_linear_expression_sum( assert len(expr.coords["dim_2"]) == 10 +@pytest.mark.legacy def test_linear_expression_sum_with_const( x: Variable, y: Variable, z: Variable, v: Variable ) -> None: + # Legacy-only: ``v.loc[:9] + v.loc[10:]`` merges disjoint coords + # (forbidden by v1 Β§8). expr = 10 * x + y + z + 10 res = expr.sum("dim_0") @@ -608,7 +624,17 @@ def test_linear_expression_multiplication_invalid( expr / x +@pytest.mark.legacy class TestCoordinateAlignment: + """ + Documents legacy positional/left-join alignment and silent NaN fill. + + The whole block is legacy-only: v1 raises on these cases (see + `test_legacy_violations.py` and convention.md Β§5/Β§8). Once later slices + flesh out v1's equivalent coverage, individual tests here can be + reclassified or removed. + """ + @pytest.fixture(params=["da", "series"]) def subset(self, request: Any) -> xr.DataArray | pd.Series: if request.param == "da": @@ -2033,15 +2059,18 @@ def test_merge(x: Variable, y: Variable, z: Variable) -> None: assert isinstance(res, LinearExpression) # now concat with same length of terms - expr1 = z.sel(dim_0=0).sum("dim_1") - expr2 = z.sel(dim_0=1).sum("dim_1") + # ``drop=True`` so the scalar ``dim_0`` coord doesn't survive each .sel + # and trip Β§11's aux-coord-conflict check (the two picks pin dim_0=0 + # vs dim_0=1). + expr1 = z.sel(dim_0=0, drop=True).sum("dim_1") + expr2 = z.sel(dim_0=1, drop=True).sum("dim_1") res = merge([expr1, expr2], dim="dim_1", cls=LinearExpression) assert res.nterm == 3 # now with different length of terms - expr1 = z.sel(dim_0=0, dim_1=slice(0, 1)).sum("dim_1") - expr2 = z.sel(dim_0=1).sum("dim_1") + expr1 = z.sel(dim_0=0, dim_1=slice(0, 1), drop=True).sum("dim_1") + expr2 = z.sel(dim_0=1, drop=True).sum("dim_1") res = merge([expr1, expr2], dim="dim_1", cls=LinearExpression) assert res.nterm == 3 @@ -2269,15 +2298,6 @@ def test_variable_names() -> None: assert expr.nterm == 2 assert expr.variable_names == {"a", "b"} - mask = xr.DataArray(False, coords=[time]) - expr = a + (b * 1).where(mask) - assert expr.nterm == 2 - assert expr.variable_names == {"a"} - - expr = (b * 1).where(mask) - assert expr.nterm == 1 - assert expr.variable_names == set() - expr = LinearExpression.from_constant(model=m, constant=5) assert expr.nterm == 0 assert expr.variable_names == set() @@ -2291,6 +2311,39 @@ def test_variable_names() -> None: assert expr.variable_names == {"a"} +@pytest.mark.legacy +def test_variable_names_masked_addend_legacy() -> None: + # Legacy: a fully masked addend's terms turn dead; the live variable remains. + m = Model() + time = pd.Index(range(3), name="time") + a = m.add_variables(name="a", coords=[time]) + b = m.add_variables(name="b", coords=[time]) + mask = xr.DataArray(False, coords=[time]) + + expr = a + (b * 1).where(mask) + assert expr.nterm == 2 + assert expr.variable_names == {"a"} + + expr = (b * 1).where(mask) + assert expr.nterm == 1 + assert expr.variable_names == set() + + +@pytest.mark.v1 +def test_variable_names_masked_addend_v1() -> None: + # v1 (section 6): absence propagates β€” a sum with a fully masked addend + # has no live terms anywhere. + m = Model() + time = pd.Index(range(3), name="time") + a = m.add_variables(name="a", coords=[time]) + b = m.add_variables(name="b", coords=[time]) + mask = xr.DataArray(False, coords=[time]) + + expr = a + (b * 1).where(mask) + assert expr.variable_names == set() + assert expr.isnull().all() + + def test_nterm() -> None: m = Model() time = pd.Index(range(3), name="time") @@ -2309,10 +2362,35 @@ def test_nterm() -> None: expr = a + b.where(all_false) assert expr.nterm == 2 - expr = expr.simplify() + +@pytest.mark.legacy +def test_nterm_simplify_collapses_fully_masked_addend() -> None: + # Legacy: a fully masked addend contributes dead terms that simplify() drops. + m = Model() + time = pd.Index(range(3), name="time") + all_false = xr.DataArray(False, coords=[time]) + a = m.add_variables(name="a", coords=[time]) + b = m.add_variables(name="b", coords=[time]) + + expr = (a + b.where(all_false)).simplify() assert expr.nterm == 1 +@pytest.mark.v1 +def test_nterm_simplify_absent_expression_has_no_terms() -> None: + # v1 (section 6): absence propagates β€” a fully masked addend absorbs the + # whole sum, so nothing is left to simplify. + m = Model() + time = pd.Index(range(3), name="time") + all_false = xr.DataArray(False, coords=[time]) + a = m.add_variables(name="a", coords=[time]) + b = m.add_variables(name="b", coords=[time]) + + expr = (a + b.where(all_false)).simplify() + assert expr.nterm == 0 + assert expr.isnull().all() + + class TestJoinParameter: @pytest.fixture def m2(self) -> Model: @@ -2335,9 +2413,12 @@ def c(self, m2: Model) -> Variable: return m2.variables["c"] class TestAddition: + @pytest.mark.legacy def test_add_join_none_preserves_default( self, a: Variable, b: Variable ) -> None: + # Legacy-only: a and b have different coords on dim ``i``; + # under v1 both arithmetic forms raise (see convention.md Β§8). result_default = a.to_linexpr() + b.to_linexpr() result_none = a.to_linexpr().add(b.to_linexpr(), join=None) assert_linequal(result_default, result_none) @@ -2596,9 +2677,11 @@ def test_div_constant_outer_fill_values(self, a: Variable) -> None: assert result.coeffs.squeeze().sel(i=0).item() == pytest.approx(1.0) class TestQuadratic: + @pytest.mark.legacy def test_quadratic_add_constant_join_inner( self, a: Variable, b: Variable ) -> None: + # Legacy-only: a*b has misaligned coords on ``i`` (Β§8 raises in v1). quad = a.to_linexpr() * b.to_linexpr() const = xr.DataArray([10, 20, 30], dims=["i"], coords={"i": [1, 2, 3]}) result = quad.add(const, join="inner") @@ -2610,9 +2693,11 @@ def test_quadratic_add_expr_join_inner(self, a: Variable) -> None: result = quad.add(const, join="inner") assert list(result.indexes["i"]) == [0, 1] + @pytest.mark.legacy def test_quadratic_mul_constant_join_inner( self, a: Variable, b: Variable ) -> None: + # Legacy-only: a*b has misaligned coords on ``i`` (Β§8 raises in v1). quad = a.to_linexpr() * b.to_linexpr() const = xr.DataArray([2, 3, 4], dims=["i"], coords={"i": [1, 2, 3]}) result = quad.mul(const, join="inner") diff --git a/test/test_optimization.py b/test/test_optimization.py index 1e771b22..c97ef28b 100644 --- a/test/test_optimization.py +++ b/test/test_optimization.py @@ -597,6 +597,7 @@ def test_duplicated_variables( assert all(np.isclose(model_with_duplicated_variables.solution["x"], 5, rtol=tol)) +@pytest.mark.legacy @pytest.mark.parametrize("solver,io_api,explicit_coordinate_names", params) def test_non_aligned_variables( model_with_non_aligned_variables: Model, @@ -604,6 +605,10 @@ def test_non_aligned_variables( io_api: str, explicit_coordinate_names: bool, ) -> None: + """ + Legacy-only: var+var on the same dim with different coords (see + convention.md Β§8). Under v1, the model construction itself raises. + """ status, condition = model_with_non_aligned_variables.solve( solver, io_api=io_api, explicit_coordinate_names=explicit_coordinate_names ) @@ -985,6 +990,7 @@ def test_modified_model( assert (modified_model.solution.y == 10).all() +@pytest.mark.legacy @pytest.mark.parametrize("solver,io_api,explicit_coordinate_names", params) def test_masked_variable_model( masked_variable_model: Model, @@ -992,6 +998,14 @@ def test_masked_variable_model( io_api: str, explicit_coordinate_names: bool, ) -> None: + """ + Legacy-only: asserts that ``x + y >= 10`` with ``y`` masked still + binds ``x >= 10`` at the masked slots β€” which only works because + legacy collapses the absent ``y`` to 0. Under v1 Β§6 the absence in + ``y`` propagates into the constraint and the constraint is dropped + at the masked slots, so ``x`` is free to be 0 there. The v1 way to + express the legacy intent is ``x + y.fillna(0) >= 10``. + """ masked_variable_model.solve( solver, io_api=io_api, explicit_coordinate_names=explicit_coordinate_names ) @@ -1006,6 +1020,70 @@ def test_masked_variable_model( assert_equal(x.add(y).solution, x.solution + y.solution.fillna(0)) +@pytest.mark.v1 +@pytest.mark.parametrize("solver,io_api,explicit_coordinate_names", params) +def test_masked_variable_model_v1_drops_constraint( + masked_variable_model: Model, + solver: str, + io_api: str, + explicit_coordinate_names: bool, +) -> None: + """ + v1 counterpart of ``test_masked_variable_model``. Under Β§6 the + absence of ``y`` at the last two slots propagates into ``x + y`` + and from there into the constraint, so the constraint drops at + those slots β€” ``x`` is no longer pinned to 10 there and the + objective ``2x + y`` drives it to 0 where it's still bound. + + Pin two things together: + 1. Model structure: con0 is masked at the absent slots (its label + is -1, no row emitted to the solver). This is the v1 invariant + that distinguishes us from legacy and is solver-independent. + 2. Solver outcome on the bound slots: ``x[:8]`` solves to 0 (the + constraint binds via ``y[:8] = 10``). ``x[-2:]`` is solver- + dependent β€” some solvers presolve away free variables and the + solution comes back as NaN β€” so we don't pin it here. + """ + con = masked_variable_model.constraints["con0"] + assert (con.labels.values[-2:] == -1).all() + assert (con.labels.values[:-2] != -1).all() + + masked_variable_model.solve( + solver, io_api=io_api, explicit_coordinate_names=explicit_coordinate_names + ) + x = masked_variable_model.variables.x + y = masked_variable_model.variables.y + tol = GPU_SOL_TOL if solver in gpu_solvers else CPU_SOL_TOL + assert y.solution[-2:].isnull().all() + assert (np.isclose(x.solution[:-2], 0, atol=tol)).all() + + +@pytest.mark.v1 +@pytest.mark.parametrize("solver,io_api,explicit_coordinate_names", params) +def test_masked_variable_model_v1_fillna_binds( + solver: str, + io_api: str, + explicit_coordinate_names: bool, +) -> None: + """ + Β§7 escape hatch under v1: ``x + y.fillna(0) >= 10`` revives the + masked slots as a present zero, so the constraint binds and the + legacy outcome (``x[-2:] == 10``) is recovered. The placement of + ``fillna`` is the caller's explicit statement of intent. + """ + m = Model() + lower = pd.Series(0, range(10)) + x = m.add_variables(lower, name="x") + mask = pd.Series([True] * 8 + [False, False]) + y = m.add_variables(lower, name="y", mask=mask) + m.add_constraints(x + y.fillna(0), GREATER_EQUAL, 10) + m.add_constraints(y, GREATER_EQUAL, 0) + m.add_objective(2 * x + y) + m.solve(solver, io_api=io_api, explicit_coordinate_names=explicit_coordinate_names) + tol = GPU_SOL_TOL if solver in gpu_solvers else CPU_SOL_TOL + assert (np.isclose(m.variables.x.solution[-2:], 10, rtol=tol)).all() + + @pytest.mark.parametrize("solver,io_api,explicit_coordinate_names", params) def test_masked_constraint_model( masked_constraint_model: Model, @@ -1273,6 +1351,7 @@ def test_auto_mask_variable_model( assert y.solution[:-2].notnull().all() +@pytest.mark.legacy @pytest.mark.parametrize("solver,io_api,explicit_coordinate_names", params) def test_auto_mask_constraint_model( auto_mask_constraint_model: Model, @@ -1280,7 +1359,11 @@ def test_auto_mask_constraint_model( io_api: str, explicit_coordinate_names: bool, ) -> None: - """Test that auto_mask=True correctly masks constraints with NaN RHS.""" + """ + Test that auto_mask=True correctly masks constraints with NaN RHS. + + Legacy-only: v1 forbids NaN constraint RHS (see convention.md Β§5/Β§12). + """ auto_mask_constraint_model.solve( solver, io_api=io_api, explicit_coordinate_names=explicit_coordinate_names ) diff --git a/test/test_piecewise_constraints.py b/test/test_piecewise_constraints.py index 72b57265..d02e87a7 100644 --- a/test/test_piecewise_constraints.py +++ b/test/test_piecewise_constraints.py @@ -1406,7 +1406,13 @@ def test_broadcast_points_dim_order_follows_exprs(self) -> None: # =========================================================================== +@pytest.mark.legacy class TestNaNMasking: + """ + NaN-as-masking patterns β€” legacy-only; v1 requires explicit ``mask=`` + or ``.where()`` (see convention.md Β§4, Β§5). + """ + def test_nan_masks_lambda_labels(self) -> None: """NaN in y_points produces masked labels in SOS2 formulation.""" m = Model() @@ -2056,8 +2062,13 @@ def test_expression_name_fallback(self) -> None: ) assert f"pwl0{PWL_LAMBDA_SUFFIX}" in m.variables + @pytest.mark.legacy def test_incremental_with_nan_mask(self) -> None: - """Incremental method with trailing NaN creates masked delta vars.""" + """ + Incremental method with trailing NaN creates masked delta vars. + + Legacy-only: NaN-as-mask in user input (see convention.md Β§5). + """ m = Model() gens = pd.Index(["a", "b"], name="gen") x = m.add_variables(coords=[gens], name="x") @@ -2334,6 +2345,7 @@ def test_convexity_invariant_to_x_direction(self) -> None: assert f_asc.method != "lp" assert f_desc.method != "lp" + @pytest.mark.legacy def test_lp_per_entity_nan_padding( self, nan_padded_pwl_model: Callable[[Method], Model] ) -> None: @@ -2341,12 +2353,15 @@ def test_lp_per_entity_nan_padding( Per-entity NaN-padded breakpoints with method='lp': padded segments must be masked out so they don't create spurious ``y ≀ 0`` constraints (bug-2 regression). + + Legacy-only: NaN-as-mask in user input (see convention.md Β§5). """ m = nan_padded_pwl_model("lp") m.solve() # f_b(10) on chord (5,10)β†’(15,15) is 12.5 assert abs(float(m.solution.sel({"entity": "b"})["y"]) - 12.5) < 1e-3 + @pytest.mark.legacy @pytest.mark.skipif(not _SOS_PATHS, reason="No SOS-capable solver installed") @pytest.mark.parametrize(("solver", "io_api"), _SOS_PATHS) def test_sos2_per_entity_nan_padding( @@ -2528,9 +2543,14 @@ def test_lp_domain_bound_infeasible_when_x_out_of_range(self) -> None: status, _ = m.solve() assert status != "ok" + @pytest.mark.legacy @pytest.mark.skipif(not _any_solvers, reason="no solver available") def test_lp_domain_uses_paired_valid_breakpoints(self) -> None: - """A trailing NaN in y must also shrink the LP x-domain.""" + """ + A trailing NaN in y must also shrink the LP x-domain. + + Legacy-only: NaN-as-mask in user input (see convention.md Β§5). + """ m = Model() x = m.add_variables(lower=0, upper=2, name="x") y = m.add_variables(lower=0, upper=10, name="y") diff --git a/test/test_variable.py b/test/test_variable.py index 8a900089..80883c3c 100644 --- a/test/test_variable.py +++ b/test/test_variable.py @@ -925,19 +925,28 @@ def test_mismatched_multiindex_raises( with pytest.raises(ValueError, match="MultiIndex.*does not match"): model.add_variables(upper=bound, coords=[midx], name="x") + @pytest.mark.legacy def test_single_level_bound_broadcasts( self, model: "Model", midx: pd.MultiIndex ) -> None: bound = DataArray([5, 6], dims=["l1"], coords={"l1": [0, 1]}) - # Implicit level projection is deprecated (scenario B) β€” warns until - # the v1 convention makes it an error. + # Implicit level projection is legacy-only (scenario B): warns under + # legacy semantics, raises under v1. with pytest.warns( - linopy.EvolvingAPIWarning, match=r"broadcasting level subset" + linopy.LinopySemanticsWarning, match=r"broadcasting level subset" ): var = model.add_variables(upper=bound, coords=[midx], name="x") assert var.dims == ("multi",) assert (var.data.upper == [5, 5, 6, 6]).all() + @pytest.mark.v1 + def test_single_level_bound_raises_v1( + self, model: "Model", midx: pd.MultiIndex + ) -> None: + bound = DataArray([5, 6], dims=["l1"], coords={"l1": [0, 1]}) + with pytest.raises(ValueError, match=r"not supported under the v1 convention"): + model.add_variables(upper=bound, coords=[midx], name="x") + def test_incomplete_level_bound_raises( self, model: "Model", midx: pd.MultiIndex ) -> None: