diff --git a/doc/release_notes.rst b/doc/release_notes.rst index f9ce4a76..008eb6c4 100644 --- a/doc/release_notes.rst +++ b/doc/release_notes.rst @@ -7,6 +7,8 @@ Upcoming Version * Add documentation about `LinearExpression.where` with `drop=True`. Add `BaseExpression.variable_names` property. * Add ``BaseExpression.has_terms`` property: boolean array, true at slots with at least one live term (`#741 `_). * ``Variable.fix(value)`` now places ``value`` correctly on variables with named dimensions; previously array values could be misaligned. +* ``Variable.fix()`` now fixes a variable by collapsing its bounds (``lower = upper = value``) instead of adding a ``__fix__`` equality constraint; ``unfix()`` restores the original bounds (`#769 `_). A fix outside the current bounds now warns and overrides instead of raising, and its shadow price appears as the variable's reduced cost rather than a constraint dual. +* Binary variable bounds are now respected by the solver, so fixing a binary works (they were previously forced to ``[0, 1]``). **Features** diff --git a/examples/manipulating-models.ipynb b/examples/manipulating-models.ipynb index 6903386b..5762eda0 100644 --- a/examples/manipulating-models.ipynb +++ b/examples/manipulating-models.ipynb @@ -357,7 +357,7 @@ "id": "30", "metadata": {}, "source": [ - "Calling `unfix()` on all variables removes the fix constraints and `unrelax()` restores the integrality of `z`." + "Calling `unfix()` on all variables restores their original bounds and `unrelax()` restores the integrality of `z`." ] }, { diff --git a/linopy/constants.py b/linopy/constants.py index 11194407..bf68ef01 100644 --- a/linopy/constants.py +++ b/linopy/constants.py @@ -37,7 +37,9 @@ class PerformanceWarning(UserWarning): short_LESS_EQUAL: LESS_EQUAL, } -FIX_CONSTRAINT_PREFIX = "__fix__" +STASHED_LOWER = "_stashed_lower" +STASHED_UPPER = "_stashed_upper" +STASHED_ATTRS: list[str] = [STASHED_LOWER, STASHED_UPPER] TERM_DIM = "_term" STACKED_TERM_DIM = "_stacked_term" diff --git a/linopy/io.py b/linopy/io.py index 62f1ca13..77530dce 100644 --- a/linopy/io.py +++ b/linopy/io.py @@ -249,6 +249,9 @@ def bounds_to_file( list(m.variables.continuous) + list(m.variables.integers) + list(m.variables.semi_continuous) + + [ + n for n in m.variables.binaries if m.variables[n].fixed + ] # fixed binaries need bounds ) if not len(list(names)): return diff --git a/linopy/model.py b/linopy/model.py index 4fea3176..d2bf021a 100644 --- a/linopy/model.py +++ b/linopy/model.py @@ -1149,16 +1149,8 @@ def remove_variables(self, name: str) -> None: ------- None. """ - from linopy.constants import FIX_CONSTRAINT_PREFIX - variable = self.variables[name] - # Clean up fix constraint if present - fix_name = f"{FIX_CONSTRAINT_PREFIX}{name}" - if fix_name in self.constraints: - self.constraints.remove(fix_name) - - # Clean up relaxed registry if present self._relaxed_registry.pop(name, None) to_remove = [ diff --git a/linopy/solvers.py b/linopy/solvers.py index 5b36a7d5..23618735 100644 --- a/linopy/solvers.py +++ b/linopy/solvers.py @@ -1268,12 +1268,6 @@ def _build_solver_model( [integrality_map[v] for v in vtypes[int_mask]], dtype=np.int32 ) h.changeColsIntegrality(len(labels), labels, integrality) - if len(model.binaries): - labels = np.arange(len(vtypes))[vtypes == "B"] - n = len(labels) - h.changeColsBounds( - n, labels, np.zeros_like(labels), np.ones_like(labels) - ) c = M.c h.changeColsCost(len(c), np.arange(len(c), dtype=np.int32), c) @@ -2827,9 +2821,6 @@ def _build_solver_model( if len(model.binaries.labels) + len(model.integers.labels) > 0: idx = [i for (i, v) in enumerate(M.vtypes) if v in ["B", "I"]] task.putvartypelist(idx, [mosek.variabletype.type_int] * len(idx)) - if len(model.binaries.labels) > 0: - bidx = [i for (i, v) in enumerate(M.vtypes) if v == "B"] - task.putvarboundlistconst(bidx, mosek.boundkey.ra, 0.0, 1.0) if len(model.constraints) > 0: if set_names: diff --git a/linopy/variables.py b/linopy/variables.py index 3f6ffebf..0eacfdc4 100644 --- a/linopy/variables.py +++ b/linopy/variables.py @@ -56,10 +56,12 @@ ) from linopy.config import options from linopy.constants import ( - FIX_CONSTRAINT_PREFIX, HELPER_DIMS, SOS_DIM_ATTR, SOS_TYPE_ATTR, + STASHED_ATTRS, + STASHED_LOWER, + STASHED_UPPER, TERM_DIM, ) from linopy.types import ( @@ -1013,7 +1015,7 @@ def flat(self) -> DataFrame: ------- df : pandas.DataFrame """ - ds = self.data + ds = self.data.drop_vars(STASHED_ATTRS, errors="ignore") def mask_func(data: pd.DataFrame) -> pd.Series: return data["labels"] != -1 @@ -1033,7 +1035,8 @@ def to_polars(self) -> pl.DataFrame: ------- pl.DataFrame """ - df = to_polars(self.data) + ds = self.data.drop_vars(STASHED_ATTRS, errors="ignore") + df = to_polars(ds) df = filter_nulls_polars(df) check_has_nulls_polars(df, name=f"{self.type} {self.name}") return df @@ -1341,10 +1344,16 @@ def fix( overwrite: bool = True, ) -> None: """ - Fix the variable to a given value by adding an equality constraint. + Fix the variable to a given value by collapsing its bounds. + + Sets ``lower = upper = value``. If no value is given, the current solution value is used. + A fix value outside the variable's current bounds emits a warning, but + does not cause infeasibilities (the bounds are overridden). Fixing a + binary variable to anything other than 0 or 1 raises. + Parameters ---------- value : float/array_like, optional @@ -1354,8 +1363,9 @@ def fix( Integer and binary variables are always rounded to 0 decimal places. Default is 8. overwrite : bool, optional - If True (default), overwrite an existing fix constraint for this - variable. If False, raise an error if the variable is already fixed. + If True (default), re-fix a variable that is already fixed to the + new value (the originally stashed bounds are kept). If False, raise + an error if the variable is already fixed. """ if value is None: try: @@ -1368,50 +1378,72 @@ def fix( ) raise ValueError(msg) from None + is_fixed = self.fixed + is_binary = self.attrs["binary"] + is_integer = self.attrs["integer"] + + if is_fixed and not overwrite: + msg = ( + f"Variable '{self.name}' is already fixed. Use " + "overwrite=True to replace the existing fix value." + ) + raise ValueError(msg) + value = broadcast_to_coords( value, self.coords, label=f"fix() for variable '{self.name}'" ) - if self.attrs.get("integer") or self.attrs.get("binary"): + if is_binary and not (np.isclose(value, 0) | np.isclose(value, 1)).all(): + msg = ( + f"Cannot fix binary variable '{self.name}' to a value " + "other than 0 or 1." + ) + raise ValueError(msg) + + if is_integer or is_binary: value = value.round(0) else: value = value.round(decimals) - if (value < self.lower).any() or (value > self.upper).any(): - msg = ( - f"Fix values for variable '{self.name}' are outside the " - "variable bounds." - ) - raise ValueError(msg) + if is_fixed: + lower, upper = self.data[STASHED_LOWER], self.data[STASHED_UPPER] + else: + lower, upper = self.data.lower, self.data.upper - constraint_name = f"{FIX_CONSTRAINT_PREFIX}{self.name}" + if not is_binary and ((value < lower).any() or (value > upper).any()): + warn( + f"Fix values for variable '{self.name}' lie outside its current " + "bounds; the bounds are overridden by the fix value.", + UserWarning, + stacklevel=2, + ) - if constraint_name in self.model.constraints: - if not overwrite: - msg = ( - f"Variable '{self.name}' is already fixed. Use " - "overwrite=True to replace the existing fix constraint." - ) - raise ValueError(msg) - self.model.remove_constraints(constraint_name) + if not is_fixed: + self._data = assign_multiindex_safe( + self.data, + **{STASHED_LOWER: lower, STASHED_UPPER: upper}, + ) - self.model.add_constraints(self, "=", value, name=constraint_name) + self.lower = value + self.upper = value def unfix(self) -> None: """ - Remove the fix constraint for this variable. + Unfix the variable, restoring the bounds it had before :meth:`fix`. """ - constraint_name = f"{FIX_CONSTRAINT_PREFIX}{self.name}" - if constraint_name in self.model.constraints: - self.model.remove_constraints(constraint_name) + if not self.fixed: + return + + self.lower = self.data[STASHED_LOWER] + self.upper = self.data[STASHED_UPPER] + self._data = self.data.drop_vars(STASHED_ATTRS) @property def fixed(self) -> bool: """ Return whether the variable is currently fixed. """ - constraint_name = f"{FIX_CONSTRAINT_PREFIX}{self.name}" - return constraint_name in self.model.constraints + return all(attr in self.data for attr in STASHED_ATTRS) class AtIndexer: @@ -1727,7 +1759,7 @@ def fix( decimals : int, optional Number of decimal places to round continuous variables to. overwrite : bool, optional - If True, overwrite existing fix constraints. + If True, re-fix variables that are already fixed. """ for var in self.data.values(): var.fix(value=value, decimals=decimals, overwrite=overwrite) diff --git a/test/test_fix_relax.py b/test/test_fix_relax.py index 666d4cfd..ad5c8de2 100644 --- a/test/test_fix_relax.py +++ b/test/test_fix_relax.py @@ -1,5 +1,6 @@ """Tests for Variable.fix(), Variable.unfix(), and Variable.fixed.""" +import warnings from pathlib import Path import numpy as np @@ -55,58 +56,92 @@ def test_fix_scalar_dtypes(self, model_with_solution: Model, value: object) -> N m = model_with_solution m.variables["x"].fix(value=value) assert m.variables["x"].fixed - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}x"] - np.testing.assert_almost_equal(con.rhs.item(), 5.0) + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 5.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 5.0) @pytest.mark.parametrize("value", ARRAY_VALUES) def test_fix_array_dtypes(self, model_with_solution: Model, value: object) -> None: m = model_with_solution m.variables["y"].fix(value=value) assert m.variables["y"].fixed - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}y"] - np.testing.assert_array_almost_equal(con.rhs.values, [2.5, -1.5]) + np.testing.assert_array_almost_equal(m.variables["y"].lower.values, [2.5, -1.5]) + np.testing.assert_array_almost_equal(m.variables["y"].upper.values, [2.5, -1.5]) def test_fix_uses_solution(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["x"].fix() assert m.variables["x"].fixed - assert f"{FIX_CONSTRAINT_PREFIX}x" in m.constraints + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 3.14159265) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 3.14159265) + + def test_fix_does_not_add_constraint(self, model_with_solution: Model) -> None: + m = model_with_solution + m.variables["x"].fix(value=5.0) + assert len(m.constraints) == 0 def test_fix_rounds_binary(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["z"].fix() - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}z"] - np.testing.assert_equal(con.rhs.item(), 1.0) + np.testing.assert_equal(m.variables["z"].lower.item(), 1.0) + np.testing.assert_equal(m.variables["z"].upper.item(), 1.0) + + def test_fix_binary_to_zero(self, model_with_solution: Model) -> None: + m = model_with_solution + m.variables["z"].fix(value=0) + np.testing.assert_equal(m.variables["z"].lower.item(), 0.0) + np.testing.assert_equal(m.variables["z"].upper.item(), 0.0) + + @pytest.mark.parametrize("value", [5, 0.4, 0.6, -1]) + def test_fix_binary_outside_domain_raises( + self, model_with_solution: Model, value: float + ) -> None: + m = model_with_solution + with pytest.raises(ValueError, match="binary variable"): + m.variables["z"].fix(value=value) def test_fix_rounds_integer(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["w"].fix() - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}w"] - np.testing.assert_equal(con.rhs.item(), 42.0) + np.testing.assert_equal(m.variables["w"].lower.item(), 42.0) + np.testing.assert_equal(m.variables["w"].upper.item(), 42.0) def test_fix_rounds_continuous(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["x"].fix(decimals=4) - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}x"] - np.testing.assert_almost_equal(con.rhs.item(), 3.1416, decimal=4) + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 3.1416, decimal=4) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 3.1416, decimal=4) - def test_fix_raises_above_upper_bound(self, model_with_solution: Model) -> None: + def test_fix_above_upper_bound_warns_and_overrides( + self, model_with_solution: Model + ) -> None: m = model_with_solution - with pytest.raises(ValueError, match="outside the variable bounds"): + with pytest.warns(UserWarning, match="outside its current"): m.variables["x"].fix(value=11.0) + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 11.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 11.0) - def test_fix_raises_below_lower_bound(self, model_with_solution: Model) -> None: + def test_fix_below_lower_bound_warns_and_overrides( + self, model_with_solution: Model + ) -> None: m = model_with_solution - with pytest.raises(ValueError, match="outside the variable bounds"): + with pytest.warns(UserWarning, match="outside its current"): m.variables["x"].fix(value=-1.0) + np.testing.assert_almost_equal(m.variables["x"].lower.item(), -1.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), -1.0) + + def test_fix_within_bounds_does_not_warn(self, model_with_solution: Model) -> None: + m = model_with_solution + with warnings.catch_warnings(): + warnings.simplefilter("error", UserWarning) + m.variables["x"].fix(value=5.0) def test_fix_small_overshoot_rounded_within_bounds( self, model_with_solution: Model ) -> None: m = model_with_solution m.variables["x"].fix(value=10.0000000001) - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}x"] - np.testing.assert_almost_equal(con.rhs.item(), 10.0) + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 10.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 10.0) def test_fix_raises_if_already_fixed_no_overwrite( self, model_with_solution: Model @@ -116,28 +151,54 @@ def test_fix_raises_if_already_fixed_no_overwrite( with pytest.raises(ValueError, match="already fixed"): m.variables["x"].fix(value=5.0, overwrite=False) - def test_fix_overwrite_replaces_existing(self, model_with_solution: Model) -> None: + def test_fix_overwrite_keeps_original_stashed_bounds( + self, model_with_solution: Model + ) -> None: m = model_with_solution m.variables["x"].fix(value=3.0) m.variables["x"].fix(value=5.0, overwrite=True) - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}x"] - np.testing.assert_almost_equal(con.rhs.item(), 5.0) + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 5.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 5.0) + m.variables["x"].unfix() + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 0.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 10.0) def test_fix_multidimensional(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["y"].fix() assert m.variables["y"].fixed - con = m.constraints[f"{FIX_CONSTRAINT_PREFIX}y"] - np.testing.assert_array_almost_equal(con.rhs.values, [2.71828, -1.41421]) + np.testing.assert_array_almost_equal( + m.variables["y"].lower.values, [2.71828, -1.41421] + ) + np.testing.assert_array_almost_equal( + m.variables["y"].upper.values, [2.71828, -1.41421] + ) class TestVariableUnfix: - def test_unfix_removes_constraint(self, model_with_solution: Model) -> None: + def test_unfix_restores_bounds(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["x"].fix(value=5.0) m.variables["x"].unfix() assert not m.variables["x"].fixed - assert f"{FIX_CONSTRAINT_PREFIX}x" not in m.constraints + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 0.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 10.0) + + def test_unfix_restores_multidimensional_bounds( + self, model_with_solution: Model + ) -> None: + m = model_with_solution + m.variables["y"].fix() + m.variables["y"].unfix() + np.testing.assert_array_almost_equal(m.variables["y"].lower.values, [-5, -5]) + np.testing.assert_array_almost_equal(m.variables["y"].upper.values, [5, 5]) + + def test_unfix_restores_binary_bounds(self, model_with_solution: Model) -> None: + m = model_with_solution + m.variables["z"].fix() + m.variables["z"].unfix() + np.testing.assert_equal(m.variables["z"].lower.item(), 0.0) + np.testing.assert_equal(m.variables["z"].upper.item(), 1.0) def test_unfix_noop_if_not_fixed(self, model_with_solution: Model) -> None: m = model_with_solution @@ -158,9 +219,7 @@ class TestUnfixDoesNotUnrelaxIndependently: def test_unfix_on_relaxed_only_variable(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["z"].relax() - # unfix should be a no-op — no fix constraint exists m.variables["z"].unfix() - # relaxation should still be in effect assert m.variables["z"].relaxed assert not m.variables["z"].attrs["binary"] @@ -182,6 +241,9 @@ def test_unfix_does_not_unrelax(self, model_with_solution: Model) -> None: m.variables["z"].relax() m.variables["z"].unfix() assert not m.variables["z"].fixed + # unfix restores the original binary bounds regardless of relaxation + np.testing.assert_equal(m.variables["z"].lower.item(), 0.0) + np.testing.assert_equal(m.variables["z"].upper.item(), 1.0) # relaxation is independent — still in effect assert m.variables["z"].relaxed assert not m.variables["z"].attrs["binary"] @@ -238,6 +300,8 @@ def test_unfix_all(self, model_with_solution: Model) -> None: m.variables.unfix() for name in m.variables: assert not m.variables[name].fixed + np.testing.assert_almost_equal(m.variables["x"].lower.item(), 0.0) + np.testing.assert_almost_equal(m.variables["x"].upper.item(), 10.0) def test_fix_integers_only(self, model_with_solution: Model) -> None: m = model_with_solution @@ -424,7 +488,7 @@ def test_remove_fixed_variable(self, model_with_solution: Model) -> None: m = model_with_solution m.variables["x"].fix(value=5.0) m.remove_variables("x") - assert f"{FIX_CONSTRAINT_PREFIX}x" not in m.constraints + assert "x" not in m.variables def test_remove_relaxed_variable(self, model_with_solution: Model) -> None: m = model_with_solution @@ -432,10 +496,45 @@ def test_remove_relaxed_variable(self, model_with_solution: Model) -> None: m.variables["z"].relax() m.remove_variables("z") assert "z" not in m._relaxed_registry - assert f"{FIX_CONSTRAINT_PREFIX}z" not in m.constraints + assert "z" not in m.variables class TestFixIO: + def test_fixed_bounds_survive_netcdf( + self, model_with_solution: Model, tmp_path: Path + ) -> None: + m = model_with_solution + m.variables["x"].fix(value=5.0) + m.variables["y"].fix() + + path = tmp_path / "model.nc" + m.to_netcdf(path) + + from linopy.io import read_netcdf + + m2 = read_netcdf(path) + assert m2.variables["x"].fixed + assert m2.variables["y"].fixed + np.testing.assert_almost_equal(m2.variables["x"].lower.item(), 5.0) + np.testing.assert_almost_equal(m2.variables["x"].upper.item(), 5.0) + + def test_unfix_after_roundtrip_restores_bounds( + self, model_with_solution: Model, tmp_path: Path + ) -> None: + m = model_with_solution + m.variables["x"].fix(value=5.0) + + path = tmp_path / "model.nc" + m.to_netcdf(path) + + from linopy.io import read_netcdf + + m2 = read_netcdf(path) + m2.variables["x"].unfix() + assert not m2.variables["x"].fixed + np.testing.assert_almost_equal(m2.variables["x"].lower.item(), 0.0) + np.testing.assert_almost_equal(m2.variables["x"].upper.item(), 10.0) + def test_relaxed_registry_survives_netcdf( self, model_with_solution: Model, tmp_path: Path ) -> None: @@ -452,9 +551,8 @@ def test_relaxed_registry_survives_netcdf( m2 = read_netcdf(path) assert m2._relaxed_registry == {"z": "binary", "w": "integer"} - # Fix constraints should also survive - assert f"{FIX_CONSTRAINT_PREFIX}z" in m2.constraints - assert f"{FIX_CONSTRAINT_PREFIX}w" in m2.constraints + assert m2.variables["z"].fixed + assert m2.variables["w"].fixed def test_empty_registry_netcdf( self, model_with_solution: Model, tmp_path: Path diff --git a/test/test_optimization.py b/test/test_optimization.py index a2912c6f..1e771b22 100644 --- a/test/test_optimization.py +++ b/test/test_optimization.py @@ -737,6 +737,50 @@ def test_milp_binary_model( ).all() +FIXED_VAR_CASES = [ + pytest.param("continuous", {}, 7.0, 100, id="continuous-lower-raised"), + pytest.param("continuous", {}, 3.0, -100, id="continuous-upper-lowered"), + pytest.param("integer", {"integer": True}, 7.0, 100, id="integer-lower-raised"), + pytest.param("integer", {"integer": True}, 3.0, -100, id="integer-upper-lowered"), + pytest.param("binary", {"binary": True}, 1.0, 100, id="binary-1"), + pytest.param("binary", {"binary": True}, 0.0, -100, id="binary-0"), +] + + +@pytest.mark.parametrize("kind,var_kwargs,fixval,coef", FIXED_VAR_CASES) +@pytest.mark.parametrize( + "solver,io_api,explicit_coordinate_names", + [p for p in params if p[0] not in ["mindopt"]], +) +def test_fixed_variable_is_held( + solver: str, + io_api: str, + explicit_coordinate_names: bool, + kind: str, + var_kwargs: dict, + fixval: float, + coef: float, +) -> None: + if kind in ("integer", "binary") and solver not in feasible_mip_solvers: + pytest.skip(f"{solver} does not support MIP") + + m = Model() + if "binary" in var_kwargs: + v = m.add_variables(name="v", **var_kwargs) + else: + v = m.add_variables(lower=0, upper=10, name="v", **var_kwargs) + x = m.add_variables(lower=0, upper=10, name="x") + m.add_constraints(x >= 2, name="c") + m.add_objective(x + coef * v) + v.fix(fixval) + + status, condition = m.solve( + solver, io_api=io_api, explicit_coordinate_names=explicit_coordinate_names + ) + assert condition == "optimal" + assert float(m.solution.v) == pytest.approx(fixval) + + @pytest.mark.parametrize( "solver,io_api,explicit_coordinate_names", [p for p in params if p[0] not in ["mindopt"]],