fix: ensure correct dtype of "reverse" column (#1525)

* fix: ensure correct dtype of "reverse" column

* apply boolean dtype to reverse column in perfect-foresight

* reinsert objective check

* fix: improve selection of bidirectionals

* fix: relax condition on sanitizing "reversed" column

* revert some unneeded changes

* debug: print dataframes in summary (to be reverted)

* filter empty dfs in make summary

* use highs-simplex instead of highs-default for ci

---------

Co-authored-by: Lukas Trippe <[email protected]>

config/config.default.yaml

@@ -963,6 +963,12 @@ solving:
       ipm_optimality_tolerance: 1e-4
       parallel: "on"
       random_seed: 123
+    highs-simplex:
+      solver: "simplex"
+      parallel: "on"
+      primal_feasibility_tolerance: 1e-5
+      dual_feasibility_tolerance: 1e-5
+      random_seed: 123
     gurobi-default:
       threads: 8
       method: 2 # barrier

config/test/config.electricity.yaml

@@ -80,10 +80,10 @@ lines:
 solving:
   solver:
     name: highs
-    options: highs-default
+    options: highs-simplex
 
   check_objective:
-    enable: false
+    enable: true
     expected_value: 3.8120188094e+07
 
 plotting:

config/test/config.myopic.yaml

@@ -85,7 +85,7 @@ industry:
 solving:
   solver:
     name: highs
-    options: highs-default
+    options: highs-simplex
   mem: 4000
 
 plotting:

config/test/config.overnight.yaml

@@ -79,12 +79,12 @@ industry:
 solving:
   solver:
     name: highs
-    options: highs-default
+    options: highs-simplex
   mem: 4000
 
   check_objective:
-    enable: false
-    expected_value: 7.0847670388e+08
+    enable: true
+    expected_value: 6.96e+08
 
 plotting:
   map:

config/test/config.perfect.yaml

@@ -84,12 +84,12 @@ industry:
 solving:
   solver:
     name: highs
-    options: highs-default
+    options: highs-simplex
   mem: 4000
 
   check_objective:
-    enable: false
-    expected_value: 1.4427662256e+10
+    enable: true
+    expected_value: 1.3778753459e+10
 
 plotting:
   map:

config/test/config.scenarios.yaml

@@ -58,4 +58,4 @@ renewable:
 solving:
   solver:
     name: highs
-    options: highs-default
+    options: highs-simplex

scripts/_helpers.py

@@ -17,6 +17,7 @@
 
 import fiona
 import pandas as pd
+import pypsa
 import pytz
 import requests
 import yaml
@@ -895,6 +896,24 @@ def get_snapshots(snapshots, drop_leap_day=False, freq="h", **kwargs):
     return time
 
 
+def sanitize_custom_columns(n: pypsa.Network):
+    """
+    Sanitize non-standard columns used throughout the workflow.
+
+    Parameters
+    ----------
+        n (pypsa.Network): The network object.
+
+    Returns
+    -------
+        None
+    """
+    if "reversed" in n.links.columns:
+        # Replace NA values with default value False
+        n.links.loc[n.links.reversed.isna(), "reversed"] = False
+        n.links.reversed = n.links.reversed.astype(bool)
+
+
 def rename_techs(label: str) -> str:
     """
     Rename technology labels for better readability.

scripts/add_existing_baseyear.py

@@ -17,6 +17,7 @@
 import xarray as xr
 from _helpers import (
     configure_logging,
+    sanitize_custom_columns,
     set_scenario_config,
     update_config_from_wildcards,
 )
@@ -620,25 +621,6 @@ def add_heating_capacities_installed_before_baseyear(
             )
 
 
-def set_defaults(n):
-    """
-    Set default values for missing values in the network.
-
-    Parameters
-    ----------
-        n (pypsa.Network): The network object.
-
-    Returns
-    -------
-        None
-    """
-    if "Link" in n.components:
-        if "reversed" in n.links.columns:
-            # Replace NA values with default value False
-            n.links.loc[n.links.reversed.isna(), "reversed"] = False
-            n.links.reversed = n.links.reversed.astype(bool)
-
-
 # %%
 if __name__ == "__main__":
     if "snakemake" not in globals():
@@ -706,13 +688,12 @@ def set_defaults(n):
         )
 
     # Set defaults for missing missing values
-    set_defaults(n)
 
     if options.get("cluster_heat_buses", False):
         cluster_heat_buses(n)
 
     n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
 
+    sanitize_custom_columns(n)
     sanitize_carriers(n, snakemake.config)
-
     n.export_to_netcdf(snakemake.output[0])

scripts/make_summary.py

@@ -338,7 +338,7 @@ def calculate_supply(n, label, supply):
             for end in [col[3:] for col in c.df.columns if col[:3] == "bus"]:
                 items = c.df.index[c.df["bus" + end].map(bus_map).fillna(False)]
 
-                if len(items) == 0:
+                if len(items) == 0 or c.pnl["p" + end].empty:
                     continue
 
                 # lots of sign compensation for direction and to do maximums
@@ -390,7 +390,7 @@ def calculate_supply_energy(n, label, supply_energy):
             for end in [col[3:] for col in c.df.columns if col[:3] == "bus"]:
                 items = c.df.index[c.df[f"bus{str(end)}"].map(bus_map).fillna(False)]
 
-                if len(items) == 0:
+                if len(items) == 0 or c.pnl["p" + end].empty:
                     continue
 
                 s = (-1) * c.pnl["p" + end][items].multiply(

scripts/prepare_perfect_foresight.py

@@ -12,9 +12,11 @@
 import pypsa
 from _helpers import (
     configure_logging,
+    sanitize_custom_columns,
     set_scenario_config,
     update_config_from_wildcards,
 )
+from add_electricity import sanitize_carriers
 from add_existing_baseyear import add_build_year_to_new_assets
 from pypsa.descriptors import expand_series
 from six import iterkeys
@@ -574,4 +576,6 @@ def update_heat_pump_efficiency(n: pypsa.Network, years: list[int]):
     update_heat_pump_efficiency(n=n, years=years)
 
     # export network
+    sanitize_custom_columns(n)
+    sanitize_carriers(n, snakemake.config)
     n.export_to_netcdf(snakemake.output[0])

scripts/prepare_sector_network.py

@@ -4378,8 +4378,8 @@ def lossy_bidirectional_links(n, carrier, efficiencies={}):
     rev_links["reversed"] = True
     rev_links.index = rev_links.index.map(lambda x: x + "-reversed")
 
+    n.links["reversed"] = n.links.get("reversed", False)
     n.links = pd.concat([n.links, rev_links], sort=False)
-    n.links["reversed"] = n.links["reversed"].fillna(False).infer_objects(copy=False)
     n.links["length_original"] = n.links["length_original"].fillna(n.links.length)
 
     # do compression losses after concatenation to take electricity consumption at bus0 in either direction

scripts/solve_network.py

@@ -784,33 +784,16 @@ def add_battery_constraints(n):
 
 
 def add_lossy_bidirectional_link_constraints(n):
-    if not n.links.p_nom_extendable.any() or "reversed" not in n.links.columns:
+    if not n.links.p_nom_extendable.any() or not any(n.links.get("reversed", [])):
         return
 
-    n.links["reversed"] = n.links.reversed.fillna(0).astype(bool)
     carriers = n.links.loc[n.links.reversed, "carrier"].unique()  # noqa: F841
-
-    forward_i = n.links.query(
-        "carrier in @carriers and ~reversed and p_nom_extendable"
+    backwards = n.links.query(
+        "carrier in @carriers and p_nom_extendable and reversed"
     ).index
-
-    def get_backward_i(forward_i):
-        return pd.Index(
-            [
-                (
-                    re.sub(r"-(\d{4})$", r"-reversed-\1", s)
-                    if re.search(r"-\d{4}$", s)
-                    else s + "-reversed"
-                )
-                for s in forward_i
-            ]
-        )
-
-    backward_i = get_backward_i(forward_i)
-
-    lhs = n.model["Link-p_nom"].loc[backward_i]
-    rhs = n.model["Link-p_nom"].loc[forward_i]
-
+    forwards = backwards.str.replace("-reversed", "")
+    lhs = n.model["Link-p_nom"].loc[backwards]
+    rhs = n.model["Link-p_nom"].loc[forwards]
     n.model.add_constraints(lhs == rhs, name="Link-bidirectional_sync")