Fix OPF time_elapsed calculation

doc/whatsnew/v0-3-1.rst

@@ -165,3 +165,14 @@ Changes
   that has no voltage issues at all before the reduction. Correcting 35 of 2421 lines
   brings the voltage band back from 0.818–1.136 to 0.895–1.040 and lets the expansion
   run through (openego/eDisGo#767)
+* Fixed the time step length handed to the OPF. ``time_elapsed`` drives every
+  inter-time-step coupling of the PowerModels problem — battery, thermal storage, DSM
+  and electromobility energy balances — and was computed as
+  ``int((snapshots[1] - snapshots[0]).seconds / 3600)``. Both parts were wrong:
+  ``.seconds`` is the seconds *component* of the difference rather than its total, so a
+  step of exactly one day came out as 0 and one of 25 h as 1 h; and the ``int()``
+  truncated any sub-hourly step to 0, which degenerates the couplings to
+  ``energy_2 == energy_1`` and decouples power from energy entirely. The value is now
+  derived from ``total_seconds()`` and returned as a float. The snapshots are also
+  checked: fewer than two, non-equidistant or non-increasing steps raise instead of
+  silently producing a meaningless global value (openego/eDisGo#693)

edisgo/io/powermodels_io.py

@@ -191,8 +191,8 @@ def to_powermodels(
     pm = _init_pm()
     timesteps = len(psa_net.snapshots)  # number of considered timesteps
     pm["name"] = f"ding0_{edisgo_object.topology.id}_t_{timesteps}"
-    pm["time_elapsed"] = int(
-        (psa_net.snapshots[1] - psa_net.snapshots[0]).seconds / 3600
+    pm["time_elapsed"] = _get_time_elapsed_in_hours(
+        psa_net.snapshots
     )  # length of timesteps in hours
     pm["baseMVA"] = s_base
     pm["source_version"] = 2
@@ -2084,6 +2084,42 @@ def _build_component_timeseries(
     pm["time_series"][kind] = pm_comp
 
 
+def _get_time_elapsed_in_hours(snapshots):
+    """
+    Calculate time elapsed in hours between two consecutive snapshots.
+
+    Parameters
+    ----------
+    snapshots : :pandas:`pandas.DatetimeIndex<DatetimeIndex>`
+        DatetimeIndex of snapshots.
+
+    Returns
+    -------
+    float
+        Time elapsed in hours between two consecutive snapshots.
+    """
+    if len(snapshots) < 2:
+        raise ValueError(
+            "At least two snapshots are required to determine time_elapsed "
+            "for the PowerModels OPF."
+        )
+
+    snapshot_deltas = snapshots.to_series().diff().dropna()
+
+    if not (snapshot_deltas == snapshot_deltas.iloc[0]).all():
+        raise ValueError(
+            "PowerModels OPF requires equidistant snapshots because "
+            "inter-timestep couplings use one global time_elapsed value."
+        )
+
+    time_elapsed = snapshot_deltas.iloc[0].total_seconds() / 3600
+
+    if time_elapsed <= 0:
+        raise ValueError("Snapshot time step must be positive.")
+
+    return time_elapsed
+
+
 def _mapping(
     psa_net,
     edisgo_obj,

tests/io/test_powermodels_io.py

@@ -4,6 +4,7 @@
 
 from edisgo import EDisGo
 from edisgo.io import powermodels_io
+from edisgo.io.powermodels_io import _get_time_elapsed_in_hours
 from edisgo.tools.tools import aggregate_district_heating_components
 
 
@@ -427,3 +428,61 @@ def test_no_check_for_opf_version_below_three(self, edisgo_obj):
         self._set_overlying_grid(edisgo_obj, full[:1], [1.0])
         pm, hv_flex_dict = powermodels_io.to_powermodels(edisgo_obj, opf_version=2)
         assert pm["opf_version"] == 2
+
+
+# test _get_time_elapsed_in_hours for inter-timestep couplings in
+# Julia to simulate
+# 1) timesteps <1h, 1h, 2h, >24h and
+# 2) whether too few, non-equidistant or
+#    negative timestamps raise an error
+def _snapshots(*timestamps):
+    return pd.to_datetime(timestamps)
+
+
+# 1) test various snapshot intervals and expected time elapsed in hours
+@pytest.mark.parametrize(
+    ("snapshots", "expected"),
+    [
+        (pd.date_range("2035-01-01", periods=3, freq="15min"), 0.25),
+        (pd.date_range("2035-01-01", periods=3, freq="h"), 1.0),
+        (pd.date_range("2035-01-01", periods=3, freq="2h"), 2.0),
+        (_snapshots("2035-01-01 00:00", "2035-01-02 01:00"), 25.0),
+    ],
+)
+# test expected time elapsed in hours for above defined snapshot intervals
+def test_get_time_elapsed_in_hours(snapshots, expected):
+    assert _get_time_elapsed_in_hours(snapshots) == pytest.approx(expected)
+
+
+# 2) test three other cases:
+# 1. not enough snapshots (1 snapshot)
+# 2. non-equidistant snapshots (15min, 1h, 1h 15min)
+# 3. negative time elapsed (snapshots in reverse order)
+@pytest.mark.parametrize(
+    ("snapshots", "error_message"),
+    [
+        (
+            _snapshots("2035-01-01 00:00"),
+            "At least two snapshots",
+        ),
+        (
+            _snapshots(
+                "2035-01-01 00:00",
+                "2035-01-01 00:15",
+                "2035-01-01 01:15",
+            ),
+            "equidistant",
+        ),
+        (
+            _snapshots(
+                "2035-01-01 01:00",
+                "2035-01-01 00:00",
+            ),
+            "positive",
+        ),
+    ],
+)
+# test that ValueError is raised for the above three cases
+def test_get_time_elapsed_in_hours_raises(snapshots, error_message):
+    with pytest.raises(ValueError, match=error_message):
+        _get_time_elapsed_in_hours(snapshots)