diff --git a/documentation/changelog.rst b/documentation/changelog.rst index f6f68c5005..9bc9c68c5f 100644 --- a/documentation/changelog.rst +++ b/documentation/changelog.rst @@ -101,6 +101,7 @@ Infrastructure / Support Bugfixes ----------- +* The big-M bounding commitment deviations now accounts for the committed quantity, so a committed quantity far beyond the devices' flow limits no longer renders the problem infeasible when a non-convex cost curve activates the commitment-sign constraints; the bound is also tighter (covering one time step rather than the whole horizon, for flow commitments), which improves MIP numerics [see `PR #2411 `_] * In a multi-device flex-model, a device without a stock (e.g. a converter port or curtailable generator) silently disabled constraint validation for all devices after it; validation now covers every device, and also newly checks that each device's power bounds do not contradict each other, so a contradictory hard bound fails with a clear per-time-step message instead of a bare solver infeasibility [see `PR #2252 `_] * The scheduler now rejects a commitment that no constraint would bind — a stock commitment naming no device or known stock group, or a commodity commitment for a commodity that no commitment maps devices to — instead of silently dropping it from the problem, or letting a favourably priced deviation make the problem unbounded; the error names the commitment [see `PR #2410 `_ and `PR #2413 `_] * Show icons for more asset types in the UI's asset structure view, which previously fell back to a question mark: the ``wind``, ``process`` and ``heat-storage`` types that FlexMeasures seeds by default, and EV infrastructure under its various names (such as ``one-way_evse``, ``two-way_evse``, ``evse``, ``charging_station`` and ``charging_hub``) and building services equipment (``hvac``, ``ahu``, ``dhw``, ``heatpump``, ``chiller``, ``lighting`` and ``other-loads``). Asset type names are now matched ignoring case and separators, so an asset type named ``charge-point`` gets the same icon as ``chargepoint`` [see `PR #2391 `_] diff --git a/flexmeasures/data/models/planning/scheduling_problem.py b/flexmeasures/data/models/planning/scheduling_problem.py index 5bb27c847f..6ac0e38986 100644 --- a/flexmeasures/data/models/planning/scheduling_problem.py +++ b/flexmeasures/data/models/planning/scheduling_problem.py @@ -24,7 +24,11 @@ from flask import current_app from pandas.tseries.frequencies import to_offset -from flexmeasures.data.models.planning import Commitment, FlowCommitment +from flexmeasures.data.models.planning import ( + Commitment, + FlowCommitment, + StockCommitment, +) from flexmeasures.data.models.planning.utils import initialize_df, initialize_series infinity = float("inf") @@ -503,9 +507,61 @@ def prepare_scheduling_problem( # noqa C901 bigM_columns = ["derivative max", "derivative min", "derivative equals"] # Compute a good value for our Big-Ms # Md is used to constrain the search space for device power - # Mc is used to constrain the search space for commitment deviations Md = np.nanmax([np.nanmax(d[bigM_columns].abs()) for d in device_constraints]) - Mc = np.nansum([np.nansum(d[bigM_columns].abs()) for d in device_constraints]) + + # Mc is used to constrain the search space for commitment deviations. + # A deviation makes up the gap between a committed quantity and an aggregated device flow (or stock change), + # so it is bounded by the largest absolute committed quantity plus the devices' summed flow limits: + # summed at any one time step for flow commitments, + # and summed over the whole horizon for stock commitments (a stock change accumulates flows since the start). + # A device without flow limits at some time step contributes nothing there (as it did before this bound was tightened), + # its power being unbounded anyway. + per_device_step_limits = [ + d[bigM_columns].astype(float).abs().max(axis=1).fillna(0).to_numpy() + for d in device_constraints + ] + per_step_total = sum(per_device_step_limits) + Mc = float(np.max(per_step_total)) + has_stock_commitment = any( + "class" in c.columns and (c["class"] == StockCommitment).any() + for c in commitments + ) + if has_stock_commitment: + # A stock change is not a raw flow: + # it passes through the derivative efficiencies (P_up * eta_up + P_down / eta_down) and includes the explicit stock delta, + # so each device's flow limits are scaled by its worst-case conversion gain, and its stock deltas are added, + # before summing over the horizon. + horizon_stock_change_limit = 0.0 + for d, flow_limits in zip(device_constraints, per_device_step_limits): + gain = np.ones(len(flow_limits)) + if "derivative up efficiency" in d.columns: + gain = np.maximum( + gain, + d["derivative up efficiency"].astype(float).fillna(1).to_numpy(), + ) + if "derivative down efficiency" in d.columns: + gain = np.maximum( + gain, + 1 + / d["derivative down efficiency"] + .astype(float) + .fillna(1) + .to_numpy(), + ) + deltas = ( + d["stock delta"].astype(float).fillna(0).abs().to_numpy() + if "stock delta" in d.columns + else 0.0 + ) + horizon_stock_change_limit += float(np.sum(flow_limits * gain + deltas)) + Mc = max(Mc, horizon_stock_change_limit) + if commitments: + quantities = np.abs( + np.concatenate([c["quantity"].to_numpy(dtype=float) for c in commitments]) + ) + finite_quantities = quantities[np.isfinite(quantities)] + if len(finite_quantities): + Mc += float(np.max(finite_quantities)) # Both Md and Mc have to be 1 MW, at least Md = max(Md, 1) diff --git a/flexmeasures/data/models/planning/tests/test_big_m.py b/flexmeasures/data/models/planning/tests/test_big_m.py new file mode 100644 index 0000000000..b0b18a90df --- /dev/null +++ b/flexmeasures/data/models/planning/tests/test_big_m.py @@ -0,0 +1,144 @@ +"""Tests for the big-M values bounding the scheduler's search space.""" + +from __future__ import annotations + +from datetime import timedelta + +import numpy as np +import pandas as pd +import pytest + +from flexmeasures.data.models.planning import FlowCommitment, StockCommitment +from flexmeasures.data.models.planning.linear_optimization import device_scheduler +from flexmeasures.data.models.planning.scheduling_problem import ( + prepare_scheduling_problem, +) +from flexmeasures.data.models.planning.utils import initialize_df + +#: Run every test in this module under both scheduler backends (see conftest). +RUN_UNDER_EACH_SOLVER = True + +COLUMNS = [ + "equals", + "max", + "min", + "efficiency", + "derivative equals", + "derivative max", + "derivative min", + "derivative down efficiency", + "derivative up efficiency", + "stock delta", +] + +START = pd.Timestamp("2020-01-01T00:00:00") +END = pd.Timestamp("2020-01-01T04:00:00") +RESOLUTION = timedelta(hours=1) + + +def make_device_constraints(power_capacity: float) -> pd.DataFrame: + device_constraints = initialize_df(COLUMNS, START, END, RESOLUTION) + device_constraints["derivative max"] = power_capacity + device_constraints["derivative min"] = -power_capacity + return device_constraints + + +def make_index() -> pd.DatetimeIndex: + return initialize_df(COLUMNS, START, END, RESOLUTION).index + + +def test_Mc_covers_flow_limits_per_step_plus_the_committed_quantity(): + """For flow commitments, a deviation spans at most the committed quantity plus one time step's summed flow limits.""" + commitment = FlowCommitment( + name="energy", + quantity=-100, + upwards_deviation_price=1, + downwards_deviation_price=-1, + index=make_index(), + ) + problem = prepare_scheduling_problem( + device_constraints=[make_device_constraints(0.5), make_device_constraints(2)], + ems_constraints=initialize_df(COLUMNS, START, END, RESOLUTION), + commitments=[commitment], + ) + assert problem.Md == 2 + assert problem.Mc == 0.5 + 2 + 100 + + +def test_Mc_covers_the_horizon_for_stock_commitments(): + """A stock commitment's deviation accumulates flows since the start, so Mc must cover the whole horizon.""" + index = make_index() + commitment = StockCommitment( + name="soc", + quantity=0.5, + upwards_deviation_price=1, + downwards_deviation_price=-1, + device=pd.Series(0, index=index), + index=index, + ) + problem = prepare_scheduling_problem( + device_constraints=[make_device_constraints(0.5), make_device_constraints(2)], + ems_constraints=initialize_df(COLUMNS, START, END, RESOLUTION), + commitments=[commitment], + ) + # 4 time steps of 0.5 + 2 flow limits each, plus the committed quantity. + assert problem.Mc == 4 * 2.5 + 0.5 + + +def test_Mc_scales_stock_changes_by_conversion_efficiencies_and_stock_deltas(): + """A stock change passes through the derivative efficiencies and includes the stock delta, unlike a raw flow.""" + index = make_index() + commitment = StockCommitment( + name="soc", + quantity=0.5, + upwards_deviation_price=1, + downwards_deviation_price=-1, + device=pd.Series(0, index=index), + index=index, + ) + device_constraints = make_device_constraints(0.5) + device_constraints["derivative up efficiency"] = 2 + device_constraints["derivative down efficiency"] = 0.5 + device_constraints["stock delta"] = 0.25 + problem = prepare_scheduling_problem( + device_constraints=[device_constraints], + ems_constraints=initialize_df(COLUMNS, START, END, RESOLUTION), + commitments=[commitment], + ) + # 4 time steps of 0.5 flow scaled by the worst-case conversion gain max(2, 1/0.5) plus a 0.25 stock delta, + # plus the committed quantity. + assert problem.Mc == 4 * (0.5 * 2 + 0.25) + 0.5 + + +def test_Mc_is_at_least_one(): + problem = prepare_scheduling_problem( + device_constraints=[make_device_constraints(0.001)], + ems_constraints=initialize_df(COLUMNS, START, END, RESOLUTION), + ) + assert problem.Mc == 1 + + +def test_large_committed_quantity_remains_feasible_under_a_non_convex_cost_curve(): + """A committed quantity far beyond the devices' flow limits must not be cut off by Mc. + + The non-convex prices (summed upwards price below summed downwards price) activate the commitment-sign constraints, + in which Mc caps the deviations. + Before Mc accounted for the committed quantity, the required deviation exceeded Mc and the problem was infeasible. + """ + commitment = FlowCommitment( + name="energy", + quantity=-100, + upwards_deviation_price=-1, + downwards_deviation_price=1, + index=make_index(), + ) + schedule, costs, results, model = device_scheduler( + device_constraints=[make_device_constraints(0.5)], + ems_constraints=initialize_df(COLUMNS, START, END, RESOLUTION), + commitments=[commitment], + ) + assert results.solver.termination_condition == "optimal" + # The upwards deviation earns 1 per unit, so the device consumes at full power. + np.testing.assert_allclose(schedule[0].values, 0.5, atol=1e-6) + # Each of the 4 steps deviates upwards by 100.5 at price -1. + assert costs == pytest.approx(-4 * 100.5) diff --git a/flexmeasures/data/tests/test_scheduling_simultaneous.py b/flexmeasures/data/tests/test_scheduling_simultaneous.py index ebf4c4cd8e..a83e0f0b5c 100644 --- a/flexmeasures/data/tests/test_scheduling_simultaneous.py +++ b/flexmeasures/data/tests/test_scheduling_simultaneous.py @@ -134,13 +134,6 @@ def test_create_simultaneous_jobs( # Define expected costs based on resolution expected_total_cost = -3.2775 - # Aggregate (unscoped) commitment semantics (issue #2379): the sample commitment - # rewarding supply binds the site's *aggregate* flow, so it stays inactive while - # the site is net-consuming and does not bias the per-device dispatch. Under the - # earlier per-device binding it wrongly rewarded the battery's supply, shifting the - # EV/battery split (EV costs were €2.3125); the total cost is unchanged either way. - expected_ev_costs = 2.2375 - expected_battery_costs = expected_total_cost - expected_ev_costs # Check costs np.testing.assert_approx_equal( @@ -149,6 +142,15 @@ def test_create_simultaneous_jobs( 4, f"Total costs should be €{expected_total_cost}, got €{total_cost}", ) + # Fairness benchmark: how the total cost is split between the devices. + # The optimum is degenerate (both devices face the same prices), + # so the split records which vertex the solver lands on, and model changes may move it. + # When one does, update the expectation deliberately, as a record of the change's fairness impact. + # PR #2411 (tightened Mc, changing the sign-constraint coefficients) moved the EV costs from €2.2375 to €2.3125, + # which happens to match the split seen under the old per-device commitment binding of issue #2379; + # that semantics is now pinned directly by the commitment-cost assertion below, rather than by the split. + expected_ev_costs = 2.3125 + expected_battery_costs = expected_total_cost - expected_ev_costs np.testing.assert_approx_equal( ev_costs, expected_ev_costs, @@ -161,6 +163,19 @@ def test_create_simultaneous_jobs( 4, f"Battery costs should be €{expected_battery_costs}, got €{battery_costs}", ) + # Aggregate (unscoped) commitment semantics (issue #2379): + # the sample commitment rewarding supply binds the site's *aggregate* flow, + # so it collects a reward only where the site as a whole net-produces. + # Under the earlier per-device binding it wrongly rewarded the battery's own supply, + # also while the site was net-consuming. + np.testing.assert_approx_equal( + job.meta["scheduler_info"]["commitment_costs"][ + "a sample commitment rewarding supply" + ], + -0.038, + 4, + "the aggregate flow commitment should be rewarded for the site's aggregate net supply only", + ) np.testing.assert_approx_equal( job.meta["scheduler_info"]["commitment_costs"]["electricity net energy"], expected_total_cost,