diff --git a/flixopt/components.py b/flixopt/components.py index 06313d7f6..1aab01a08 100644 --- a/flixopt/components.py +++ b/flixopt/components.py @@ -17,7 +17,13 @@ from .elements import Component, ComponentModel, Flow from .features import InvestmentModel, PiecewiseModel from .interface import InvestParameters, PiecewiseConversion, StatusParameters -from .modeling import BoundingPatterns, _scalar_safe_isel, _scalar_safe_isel_drop, _scalar_safe_reduce +from .modeling import ( + BoundingPatterns, + _scalar_safe_isel, + _scalar_safe_isel_drop, + _scalar_safe_reduce, + _set_constraint_lhs, +) from .structure import FlowSystemModel, VariableCategory, register_class_for_io if TYPE_CHECKING: @@ -849,7 +855,10 @@ def create_transmission_equation(self, name: str, in_flow: Flow, out_flow: Flow) ) if (self.element.absolute_losses is not None) and np.any(self.element.absolute_losses != 0): - con_transmission.lhs += in_flow.submodel.status.status * self.element.absolute_losses + _set_constraint_lhs( + con_transmission, + con_transmission.lhs + in_flow.submodel.status.status * self.element.absolute_losses, + ) return con_transmission diff --git a/flixopt/elements.py b/flixopt/elements.py index 446ef4bd7..760883fa1 100644 --- a/flixopt/elements.py +++ b/flixopt/elements.py @@ -16,7 +16,7 @@ from .core import PlausibilityError from .features import InvestmentModel, StatusModel from .interface import InvestParameters, StatusParameters -from .modeling import BoundingPatterns, ModelingPrimitives, ModelingUtilitiesAbstract +from .modeling import BoundingPatterns, ModelingPrimitives, ModelingUtilitiesAbstract, _set_constraint_lhs from .structure import ( Element, ElementModel, @@ -1033,7 +1033,7 @@ def _do_modeling(self): ) # Σ(inflows) + virtual_supply = Σ(outflows) + virtual_demand - eq_bus_balance.lhs += self.virtual_supply - self.virtual_demand + _set_constraint_lhs(eq_bus_balance, eq_bus_balance.lhs + self.virtual_supply - self.virtual_demand) # Add penalty shares as temporal effects (time-dependent) from .effects import PENALTY_EFFECT_LABEL diff --git a/flixopt/features.py b/flixopt/features.py index e85636435..cf214ea15 100644 --- a/flixopt/features.py +++ b/flixopt/features.py @@ -10,7 +10,7 @@ import linopy import numpy as np -from .modeling import BoundingPatterns, ModelingPrimitives, ModelingUtilities +from .modeling import BoundingPatterns, ModelingPrimitives, ModelingUtilities, _set_constraint_lhs from .structure import FlowSystemModel, Submodel, VariableCategory if TYPE_CHECKING: @@ -664,7 +664,9 @@ def _do_modeling(self): # Add it to the total (cluster_weight handles cluster representation, defaults to 1.0) # Sum over all temporal dimensions (time, and cluster if present) weighted_per_timestep = self.total_per_timestep * self._model.weights.get('cluster', 1.0) - self._eq_total.lhs -= weighted_per_timestep.sum(dim=self._model.temporal_dims) + _set_constraint_lhs( + self._eq_total, self._eq_total.lhs - weighted_per_timestep.sum(dim=self._model.temporal_dims) + ) def add_share( self, @@ -694,7 +696,7 @@ def add_share( raise ValueError('Cannot add share with scenario-dim to a model without scenario-dim') if name in self.shares: - self.share_constraints[name].lhs -= expression + _set_constraint_lhs(self.share_constraints[name], self.share_constraints[name].lhs - expression) else: # Temporal shares (with 'time' dim) are segment totals that need division category = VariableCategory.SHARE if 'time' in dims else None @@ -710,6 +712,6 @@ def add_share( ) if 'time' not in dims: - self._eq_total.lhs -= self.shares[name] + _set_constraint_lhs(self._eq_total, self._eq_total.lhs - self.shares[name]) else: - self._eq_total_per_timestep.lhs -= self.shares[name] + _set_constraint_lhs(self._eq_total_per_timestep, self._eq_total_per_timestep.lhs - self.shares[name]) diff --git a/flixopt/modeling.py b/flixopt/modeling.py index ff84c808f..e2bc59662 100644 --- a/flixopt/modeling.py +++ b/flixopt/modeling.py @@ -76,6 +76,22 @@ def _scalar_safe_reduce(data: xr.DataArray | Any, dim: str, method: str = 'mean' return data +def _set_constraint_lhs(constraint: linopy.Constraint, lhs) -> None: + """Replace a constraint's LHS, compatibly across linopy versions. + + Incremental accumulators (e.g. share totals) must mutate a constraint after + it is created. linopy >= 0.8 exposes ``Constraint.update(lhs=...)`` and + deprecates the ``.lhs`` setter (flixopt escalates that deprecation to an + error); linopy < 0.8 has only the setter, which is not deprecated there. + Dispatch on whichever API the installed linopy provides so the same call + works — and warns on neither — under both. + """ + if hasattr(constraint, 'update'): + constraint.update(lhs=lhs) + else: + constraint.lhs = lhs + + def _xr_allclose(a: xr.DataArray, b: xr.DataArray, rtol: float = 1e-5, atol: float = 1e-8) -> bool: """Check if two DataArrays are element-wise equal within tolerance. diff --git a/pyproject.toml b/pyproject.toml index d9f6bc18b..59f903b27 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -36,7 +36,7 @@ dependencies = [ "pandas >= 2.0.0, < 3", "xarray >=2024.2.0, <2026.5", # CalVer: allow through next calendar year # Optimization and data handling - "linopy >=0.5.1, <0.8", # Widened from patch pin to minor range + "linopy >=0.5.1, <0.10", # Widened from patch pin to minor range "netcdf4 >=1.6.1, <1.7.5", # 1.7.4 missing wheels, revert to < 2 later # Utilities "pyyaml >= 6.0.0, < 7", diff --git a/tests/test_clustering/test_expansion_regression.py b/tests/test_clustering/test_expansion_regression.py index 065c33ddd..eb9f96ef8 100644 --- a/tests/test_clustering/test_expansion_regression.py +++ b/tests/test_clustering/test_expansion_regression.py @@ -93,8 +93,13 @@ def test_expanded_storage(self, system_with_storage, solver_fixture): fs_e = fs_c.transform.expand() sol = fs_e.solution - # Storage dispatch varies by solver — check charge_state is non-trivial - assert float(np.nansum(sol['S|charge_state'].values)) > 0 + # Gas price and boiler efficiency are flat, so storage earns nothing and its + # dispatch is degenerate: every cycling depth (including none) is optimal, and + # which one the solver lands on shifts with solver and linopy version. Assert + # the expansion produced a usable charge_state, not one arbitrary optimum. + charge_state = sol['S|charge_state'].values + assert charge_state.shape == (N_HOURS + 1,) + assert np.isfinite(charge_state).all() # Net discharge should be ~0 (balanced storage) assert float(np.nansum(sol['S|netto_discharge'].values)) == pytest.approx(0, abs=1e-4) diff --git a/tests/test_math/test_clustering.py b/tests/test_math/test_clustering.py index aaa37923c..f5dfc0de4 100644 --- a/tests/test_math/test_clustering.py +++ b/tests/test_math/test_clustering.py @@ -336,11 +336,15 @@ def test_storage_cyclic_charge_discharge_pattern(self, optimize): discharge_fr = fs.solution['Battery(discharge)|flow_rate'].values[:, :4] assert_allclose(discharge_fr, [[0, 50, 0, 50], [0, 50, 0, 50]], atol=1e-5) + # Both clusters carry identical data, so the absolute SOC offset is degenerate: + # any level in [50, 100] is optimal. Assert what the model actually determines -- + # the charge/discharge pattern and the cyclic wrap -- not one arbitrary offset. charge_state = fs.solution['Battery|charge_state'] assert charge_state.dims == ('cluster', 'time') - cs_c0 = charge_state.isel(cluster=0).values[:5] - cs_c1 = charge_state.isel(cluster=1).values[:5] - assert_allclose(cs_c0, [50, 50, 0, 50, 0], atol=1e-5) - assert_allclose(cs_c1, [100, 100, 50, 100, 50], atol=1e-5) + for cluster in (0, 1): + cs = charge_state.isel(cluster=cluster).values[:5] + assert_allclose(np.diff(cs), [0, -50, 50, -50], atol=1e-5) + assert_allclose(cs[0], cs[3], atol=1e-5) # cyclic wrap within the cluster + assert 50 - 1e-5 <= cs[0] <= 100 + 1e-5 assert_allclose(fs.solution['objective'].item(), 100.0, rtol=1e-5)