Merge branch 'main' into 928-fix-discrepancies-between-db-and-models

Signed-off-by: GustaafL <41048720+GustaafL@users.noreply.github.com>

.vscode/settings.json

@@ -18,5 +18,6 @@
         "flexmeasures"
     ],
     "python.testing.unittestEnabled": false,
-    "python.testing.pytestEnabled": true
+    "python.testing.pytestEnabled": true,
+    "python.analysis.autoImportCompletions": true
 }

documentation/changelog.rst

@@ -12,15 +12,17 @@ New features
 -------------
 
 * Better navigation experience through listings (sensors / assets / users / accounts) in the :abbr:`UI (user interface)`, by heading to the selected entity upon a click (or CTRL + click) anywhere within a row [see `PR #923 <https://github.com/FlexMeasures/flexmeasures/pull/923>`_]
-* New flexmeasures configuration setting `FLEXMEASURES_ENFORCE_SECURE_CONTENT_POLICY` for upgrading insecure `http` requests to secured requests `https` [see `PR #920 <https://github.com/FlexMeasures/flexmeasures/pull/920>`_]
+* Define device-level power constraints as sensors to create schedules with changing power limits. [see `PR #897 <https://github.com/FlexMeasures/flexmeasures/pull/897>`_]
+* Allow to provide external storage usage or gain components using the ``soc-usage`` and ``soc-gain`` fields of the `flex-model` [see `PR #906 <https://github.com/FlexMeasures/flexmeasures/pull/906>`_]
 
 Infrastructure / Support
 ----------------------
 
 * Align database and models of `annotations`, `data_sources`, and `timed_belief` [see `PR #929 <https://github.com/FlexMeasures/flexmeasures/pull/929>`_]
-* Remove obsolete database tables `price`, `power`, `market`, `market_type`, `weather`, `asset`, and `weather_sensor` [see `PR #921 <https://github.com/FlexMeasures/flexmeasures/pull/921>`_]
 * New documentation section on constructing a flex model for :abbr:`V2G (vehicle-to-grid)` [see `PR #885 <https://github.com/FlexMeasures/flexmeasures/pull/885>`_]
 * Allow charts in plugins to show currency codes (such as EUR) as currency symbols (€) [see `PR #922 <https://github.com/FlexMeasures/flexmeasures/pull/922>`_]
+* Remove obsolete database tables `price`, `power`, `market`, `market_type`, `weather`, `asset`, and `weather_sensor` [see `PR #921 <https://github.com/FlexMeasures/flexmeasures/pull/921>`_]
+* New flexmeasures configuration setting `FLEXMEASURES_ENFORCE_SECURE_CONTENT_POLICY` for upgrading insecure `http` requests to secured requests `https` [see `PR #920 <https://github.com/FlexMeasures/flexmeasures/pull/920>`_]
 
 Bugfixes
 -----------

documentation/concepts/device_scheduler.rst

@@ -48,7 +48,8 @@ Symbol                              Variable in the Code
 :math:`P^{ems}_{min}(j)`              ems_derivative_min                                 Minimum flow of the EMS during time period :math:`j`.
 :math:`P^{ems}_{max}(j)`              ems_derivative_max                                 Maximum flow of the EMS during time period :math:`j`.
 :math:`Commitment(c,j)`               commitment_quantity                                Commitment c (at EMS level) over time step :math:`j`.
-:math:`M`                             M                                                  Large constant number, upper bound of :math:`Power_{up}(d,j)` and :math:`|Power_{down}(d,j)|`
+:math:`M`                             M                                                  Large constant number, upper bound of :math:`Power_{up}(d,j)` and :math:`|Power_{down}(d,j)|`.
+:math:`D(d,j)`                        stock_delta                                        Explicit energy gain or loss of device :math:`d` during time period :math:`j`.
 ================================ ================================================ ==============================================================================================================  
 
 
@@ -86,8 +87,7 @@ change of :math:`Stock(d,j)`, taking into account conversion efficiencies but no
 .. math::
   :name: stock
 
-    \Delta Stock(d,j) = \frac{P_{down}(d,j)}{\eta_{down}(d,j) } + P_{up}(d,j)  \cdot \eta_{up}(d,j)
-
+    \Delta Stock(d,j) = \frac{P_{down}(d,j)}{\eta_{down}(d,j) } + P_{up}(d,j)  \cdot \eta_{up}(d,j) + D(d,j)
 
 
 .. math:: 

flexmeasures/api/v3_0/sensors.py

@@ -295,6 +295,9 @@ def trigger_schedule(
         and aggregate production should be priced by sensor 10,
         where the aggregate power flow in the EMS is described by the sum over sensors 13, 14 and 15
         (plus the flexible sensor being optimized, of course).
+
+        The battery consumption power capacity is limited by sensor 42 and the production capacity is constant (30 kW).
+
         Note that, if forecasts for sensors 13, 14 and 15 are not available, a schedule cannot be computed.
 
         .. code-block:: json
@@ -325,7 +328,9 @@ def trigger_schedule(
                     "soc-max": 25,
                     "roundtrip-efficiency": 0.98,
                     "storage-efficiency": 0.9999,
-                    "power-capacity": "25kW"
+                    "power-capacity": "25kW",
+                    "consumption-capacity" : {"sensor" : 42},
+                    "production-capacity" : "30 kW"
                 },
                 "flex-context": {
                     "consumption-price-sensor": 9,

flexmeasures/conftest.py

@@ -778,6 +778,33 @@ def create_test_battery_assets(
     )
     db.session.add(test_battery_sensor_no_prices)
 
+    test_battery_dynamic_power_capacity = GenericAsset(
+        name="Test battery with dynamic power capacity",
+        owner=setup_accounts["Prosumer"],
+        generic_asset_type=battery_type,
+        latitude=10,
+        longitude=100,
+        attributes=dict(
+            capacity_in_mw=10,
+            max_soc_in_mwh=20,
+            min_soc_in_mwh=0,
+            soc_in_mwh=2.0,
+            market_id=setup_markets["epex_da"].id,
+        ),
+    )
+    test_battery_dynamic_capacity_power_sensor = Sensor(
+        name="power",
+        generic_asset=test_battery_dynamic_power_capacity,
+        event_resolution=timedelta(minutes=15),
+        unit="MW",
+        attributes=dict(
+            capacity_in_mw=10,
+            production_capacity="8 MW",
+            consumption_capacity="0.5 MW",
+        ),
+    )
+    db.session.add(test_battery_dynamic_capacity_power_sensor)
+
     test_small_battery = GenericAsset(
         name="Test small battery",
         owner=setup_accounts["Prosumer"],
@@ -817,6 +844,7 @@ def create_test_battery_assets(
         "Test battery": test_battery,
         "Test battery with no known prices": test_battery_no_prices,
         "Test small battery": test_small_battery,
+        "Test battery with dynamic power capacity": test_battery_dynamic_power_capacity,
     }
 
 
@@ -1075,3 +1103,67 @@ def error_generator():
     @roles_accepted(ADMIN_ROLE)
     def vips_only():
         return jsonify({"message": "Nothing bad happened."}), 200
+
+
+@pytest.fixture(scope="module")
+def capacity_sensors(db, add_battery_assets, setup_sources):
+    battery = add_battery_assets["Test battery with dynamic power capacity"]
+    production_capacity_sensor = Sensor(
+        name="production capacity",
+        generic_asset=battery,
+        unit="kW",
+        event_resolution="PT15M",
+        attributes={"consumption_is_positive": True},
+    )
+    consumption_capacity_sensor = Sensor(
+        name="consumption capacity",
+        generic_asset=battery,
+        unit="kW",
+        event_resolution="PT15M",
+        attributes={"consumption_is_positive": True},
+    )
+
+    db.session.add_all([production_capacity_sensor, consumption_capacity_sensor])
+    db.session.flush()
+
+    time_slots = pd.date_range(
+        datetime(2015, 1, 2), datetime(2015, 1, 2, 7, 45), freq="15T"
+    ).tz_localize("Europe/Amsterdam")
+
+    values = [200] * 4 * 4 + [300] * 4 * 4
+
+    beliefs = [
+        TimedBelief(
+            event_start=dt,
+            belief_horizon=parse_duration("PT0M"),
+            event_value=val,
+            sensor=production_capacity_sensor,
+            source=setup_sources["Seita"],
+        )
+        for dt, val in zip(time_slots, values)
+    ]
+    db.session.add_all(beliefs)
+    db.session.commit()
+
+    time_slots = pd.date_range(
+        datetime(2015, 1, 2), datetime(2015, 1, 2, 7, 45), freq="15T"
+    ).tz_localize("Europe/Amsterdam")
+
+    values = [250] * 4 * 4 + [150] * 4 * 4
+
+    beliefs = [
+        TimedBelief(
+            event_start=dt,
+            belief_horizon=parse_duration("PT0M"),
+            event_value=val,
+            sensor=consumption_capacity_sensor,
+            source=setup_sources["Seita"],
+        )
+        for dt, val in zip(time_slots, values)
+    ]
+    db.session.add_all(beliefs)
+    db.session.commit()
+
+    yield dict(
+        production=production_capacity_sensor, consumption=consumption_capacity_sensor
+    )

flexmeasures/data/models/planning/linear_optimization.py

@@ -52,6 +52,7 @@ def device_scheduler(  # noqa C901
         derivative equals: exact amount of flow (we do this by clamping derivative min and derivative max)
         derivative down efficiency: conversion efficiency of flow out of a device (flow out : stock decrease)
         derivative up efficiency: conversion efficiency of flow into a device (stock increase : flow in)
+        stock delta: predefined stock delta to apply to the storage device. Positive values cause an increase and negative values a decrease
     EMS constraints are on an EMS level. Handled constraints (listed by column name):
         derivative max: maximum flow
         derivative min: minimum flow
@@ -102,6 +103,14 @@ def device_scheduler(  # noqa C901
     # M has to be 1 MW, at least
     M = max(M, 1)
 
+    for d in range(len(device_constraints)):
+        if "stock delta" not in device_constraints[d].columns:
+            device_constraints[d]["stock delta"] = 0
+        else:
+            device_constraints[d]["stock delta"] = device_constraints[d][
+                "stock delta"
+            ].fillna(0)
+
     # Turn prices per commitment into prices per commitment flow
     if len(commitment_downwards_deviation_price) != 0:
         if all(
@@ -223,6 +232,9 @@ def device_derivative_up_efficiency(m, d, j):
             return 1
         return eff
 
+    def device_stock_delta(m, d, j):
+        return device_constraints[d]["stock delta"].iloc[j]
+
     model.up_price = Param(model.c, model.j, initialize=price_up_select)
     model.down_price = Param(model.c, model.j, initialize=price_down_select)
     model.commitment_quantity = Param(
@@ -245,6 +257,7 @@ def device_derivative_up_efficiency(m, d, j):
     model.device_derivative_up_efficiency = Param(
         model.d, model.j, initialize=device_derivative_up_efficiency
     )
+    model.stock_delta = Param(model.d, model.j, initialize=device_stock_delta)
 
     # Add variables
     model.ems_power = Var(model.d, model.j, domain=Reals, initialize=0)
@@ -268,6 +281,7 @@ def device_bounds(m, d, j):
             (
                 m.device_power_down[d, k] / m.device_derivative_down_efficiency[d, k]
                 + m.device_power_up[d, k] * m.device_derivative_up_efficiency[d, k]
+                + m.stock_delta[d, k]
             )
             for k in range(0, j + 1)
         ]

flexmeasures/data/models/planning/storage.py

@@ -9,7 +9,6 @@
 import numpy as np
 from flask import current_app
 
-
 from flexmeasures.data.models.planning import Scheduler, SchedulerOutputType
 from flexmeasures.data.models.planning.linear_optimization import device_scheduler
 from flexmeasures.data.models.planning.utils import (
@@ -19,13 +18,14 @@
     initialize_df,
     get_power_values,
     fallback_charging_policy,
+    get_continuous_series_sensor_or_quantity,
 )
 from flexmeasures.data.models.planning.exceptions import InfeasibleProblemException
 from flexmeasures.data.schemas.scheduling.storage import StorageFlexModelSchema
 from flexmeasures.data.schemas.scheduling import FlexContextSchema
 from flexmeasures.utils.time_utils import get_max_planning_horizon
 from flexmeasures.utils.coding_utils import deprecated
-from flexmeasures.utils.unit_utils import ur
+from flexmeasures.utils.unit_utils import ur, convert_units
 
 
 def check_and_convert_power_capacity(
@@ -63,6 +63,7 @@ class MetaStorageScheduler(Scheduler):
         "derivative min",
         "derivative down efficiency",
         "derivative up efficiency",
+        "stock delta",
     ]
 
     def compute_schedule(self) -> pd.Series | None:
@@ -198,14 +199,71 @@ def _prepare(self, skip_validation: bool = False) -> tuple:  # noqa: C901
             soc_min,
         )
 
+        consumption_capacity = self.flex_model.get("consumption_capacity")
+        production_capacity = self.flex_model.get("production_capacity")
+
         if sensor.get_attribute("is_strictly_non_positive"):
             device_constraints[0]["derivative min"] = 0
         else:
-            device_constraints[0]["derivative min"] = power_capacity_in_mw * -1
+            device_constraints[0]["derivative min"] = (
+                -1
+            ) * get_continuous_series_sensor_or_quantity(
+                quantity_or_sensor=production_capacity,
+                actuator=sensor,
+                unit=sensor.unit,
+                query_window=(start, end),
+                resolution=resolution,
+                beliefs_before=belief_time,
+                fallback_attribute="production_capacity",
+                max_value=convert_units(power_capacity_in_mw, "MW", sensor.unit),
+            )
         if sensor.get_attribute("is_strictly_non_negative"):
             device_constraints[0]["derivative max"] = 0
         else:
-            device_constraints[0]["derivative max"] = power_capacity_in_mw
+            device_constraints[0][
+                "derivative max"
+            ] = get_continuous_series_sensor_or_quantity(
+                quantity_or_sensor=consumption_capacity,
+                actuator=sensor,
+                unit=sensor.unit,
+                query_window=(start, end),
+                resolution=resolution,
+                beliefs_before=belief_time,
+                fallback_attribute="consumption_capacity",
+                max_value=convert_units(power_capacity_in_mw, "MW", sensor.unit),
+            )
+
+        soc_gain = self.flex_model.get("soc_gain", [])
+        soc_usage = self.flex_model.get("soc_usage", [])
+
+        all_stock_delta = []
+
+        for is_usage, soc_delta in zip([False, True], [soc_gain, soc_usage]):
+            for component in soc_delta:
+                stock_delta_series = get_continuous_series_sensor_or_quantity(
+                    quantity_or_sensor=component,
+                    actuator=sensor,
+                    unit="MW",
+                    query_window=(start, end),
+                    resolution=resolution,
+                    beliefs_before=belief_time,
+                )
+
+                # example: 4 MW sustained over 15 minutes gives 1 MWh
+                stock_delta_series *= resolution / timedelta(
+                    hours=1
+                )  # MW -> MWh / resolution
+
+                if is_usage:
+                    stock_delta_series *= -1
+
+                all_stock_delta.append(stock_delta_series)
+
+        if len(all_stock_delta) > 0:
+            all_stock_delta = pd.concat(all_stock_delta, axis=1)
+
+            device_constraints[0]["stock delta"] = all_stock_delta.sum(1)
+            device_constraints[0]["stock delta"] *= timedelta(hours=1) / resolution
 
         # Apply round-trip efficiency evenly to charging and discharging
         device_constraints[0]["derivative down efficiency"] = (

flexmeasures/data/models/planning/tests/conftest.py

@@ -219,6 +219,59 @@ def process(db, building, setup_sources) -> dict[str, Sensor]:
     return _process
 
 
+@pytest.fixture(scope="module")
+def add_stock_delta(db, add_battery_assets, setup_sources) -> dict[str, Sensor]:
+    """
+    Different usage forecast sensors are defined:
+        - "delta fails": the usage forecast exceeds the maximum power.
+        - "delta": the usage forecast can be fulfilled just right. This coincides with the schedule resolution.
+        - "delta hourly": the event resolution is changed to test that the schedule is still feasible.
+                          This has a greater resolution.
+        - "delta 5min": the event resolution is reduced even more. This sensor has a resolution smaller than that used
+                        for the scheduler.
+    """
+
+    battery = add_battery_assets["Test battery"]
+    capacity = battery.get_attribute("capacity_in_mw")
+    sensors = {}
+    sensor_specs = [
+        ("delta fails", timedelta(minutes=15), capacity * 1.2),
+        ("delta", timedelta(minutes=15), capacity),
+        ("delta hourly", timedelta(hours=1), capacity),
+        ("delta 5min", timedelta(minutes=5), capacity),
+    ]
+
+    for name, resolution, value in sensor_specs:
+        # 1 days of test data
+        time_slots = initialize_index(
+            start=pd.Timestamp("2015-01-01").tz_localize("Europe/Amsterdam"),
+            end=pd.Timestamp("2015-01-02").tz_localize("Europe/Amsterdam"),
+            resolution=resolution,
+        )
+
+        stock_delta_sensor = Sensor(
+            name=name,
+            unit="MW",
+            event_resolution=resolution,
+            generic_asset=battery,
+        )
+        db.session.add(stock_delta_sensor)
+        db.session.flush()
+
+        stock_gain = [value] * len(time_slots)
+
+        add_as_beliefs(
+            db,
+            stock_delta_sensor,
+            stock_gain,
+            time_slots,
+            setup_sources["Seita"],
+        )
+        sensors[name] = stock_delta_sensor
+
+    return sensors
+
+
 def add_as_beliefs(db, sensor, values, time_slots, source):
     beliefs = [
         TimedBelief(