included Lisas PR suggestions

config/config.default.yaml

@@ -577,13 +577,15 @@ sector:
       solar: 3
       offwind-ac: 3
       offwind-dc: 3
-  enhanced_geothermal: false
+  enhanced_geothermal:
-  enhanced_geothermal_optimism: false # if true, egs costs are reducing towards 2050 according to Aghahosseini et al., (2020)
+    enable: false
-  enhanced_geothermal_performant: true # if true, adds only the cheapest patch of EGS potential to each region
+    optimism: false # if true, egs costs are reducing towards 2050 according to Aghahosseini et al., (2020)
-  enhanced_geothermal_flexible: true # if true, adds a storage unit simulating flexible operation of EGS, see Ricks et al. 2023
+    performant: true # if true, adds only the cheapest patch of EGS potential to each region
-  enhanced_geothermal_var_cf: true # if true, adds time-dependent capacity factor to EGS, see Ricks et al. 2023
+    flexible: true # if true, adds a storage unit simulating flexible operation of EGS, see Ricks et al. 2023
-  enhanced_geothermal_reservoir_max_hours: 240 # relavant for flexible EGS, see Ricks et al. 2023
+    var_cf: true # if true, adds time-dependent capacity factor to EGS, see Ricks et al. 2023
-  enhanced_geothermal_reservoir_max_boost: 0.25 # share of generation that can be added by flexible EGS, see Ricks et al. 2023
+    max_hours: 240 # relavant for flexible EGS, see Ricks et al. 2023
+    max_boost: 0.25 # share of generation that can be added by flexible EGS, see Ricks et al. 2023
+    sustainability_factor: 0.0025 # share of capacity that is replenished from earth's core
 
 # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#industry
 industry:

envs/environment.yaml

@@ -20,7 +20,7 @@ dependencies:
 - openpyxl!=3.1.1
 - pycountry
 - seaborn
-- snakemake-minimal>=8.5
+- snakemake-minimal>=8.5,<8.6
 - memory_profiler
 - yaml
 - pytables

rules/build_sector.smk

@@ -816,36 +816,6 @@ def input_profile_offwind(w):
     }
 
 
-rule build_egs_potentials:
-    params:
-        snapshots=config_provider("snapshots"),
-    input:
-        egs_cost="data/egs_costs.json",
-        regions=resources("regions_onshore_elec_s{simpl}_{clusters}.geojson"),
-        air_temperature=(
-            resources("temp_air_total_elec_s{simpl}_{clusters}.nc")
-            if config_provider("sector", "enhanced_geothermal_var_cf")
-            else []
-        ),
-    output:
-        egs_potentials=resources("egs_potentials_s{simpl}_{clusters}.csv"),
-        egs_overlap=resources("egs_overlap_s{simpl}_{clusters}.csv"),
-        egs_capacity_factors=(
-            resources("egs_capacity_factors_s{simpl}_{clusters}.csv")
-            if config_provider("sector", "enhanced_geothermal_var_cf")
-            else []
-        ),
-    threads: 1
-    resources:
-        mem_mb=2000,
-    log:
-        logs("build_egs_potentials_s{simpl}_{clusters}.log"),
-    conda:
-        "../envs/environment.yaml"
-    script:
-        "../scripts/build_egs_potentials.py"
-
-
 rule prepare_sector_network:
     params:
         time_resolution=config_provider("clustering", "temporal", "resolution_sector"),
@@ -946,13 +916,6 @@ rule prepare_sector_network:
         cop_air_total=resources("cop_air_total_elec_s{simpl}_{clusters}.nc"),
         cop_air_rural=resources("cop_air_rural_elec_s{simpl}_{clusters}.nc"),
         cop_air_urban=resources("cop_air_urban_elec_s{simpl}_{clusters}.nc"),
-        egs_potentials=RESOURCES + "egs_potentials_s{simpl}_{clusters}.csv",
-        egs_overlap=RESOURCES + "egs_overlap_s{simpl}_{clusters}.csv",
-        egs_capacity_factors=(
-            RESOURCES + "egs_capacity_factors_s{simpl}_{clusters}.csv"
-            if config["sector"]["enhanced_geothermal_var_cf"]
-            else []
-        ),
         solar_thermal_total=lambda w: (
             resources("solar_thermal_total_elec_s{simpl}_{clusters}.nc")
             if config_provider("sector", "solar_thermal")(w)

scripts/build_egs_potentials.py

@@ -93,6 +93,10 @@ def get_capacity_factors(network_regions_file, air_temperatures_file):
     in Decarbonized Elec Systems.
     """
 
+    # these values are taken from the paper's
+    # Supplementary Figure 20 from https://zenodo.org/records/7093330
+    # and relate deviations of the ambient temperature from the year-average
+    # ambient temperature to EGS capacity factors.
     delta_T = [-15, -10, -5, 0, 5, 10, 15, 20]
     cf = [1.17, 1.13, 1.07, 1, 0.925, 0.84, 0.75, 0.65]
 
@@ -118,6 +122,10 @@ def get_capacity_factors(network_regions_file, air_temperatures_file):
     snapshots = pd.date_range(freq="h", **snakemake.params.snapshots)
     capacity_factors = pd.DataFrame(index=snapshots)
 
+    # bespoke computation of capacity factors for each bus.
+    # Considering the respective temperatures, we compute
+    # the deviation from the average temperature and relate it
+    # to capacity factors based on the data from above.
     for bus in index:
         temp = air_temp.sel(name=bus).to_dataframe()["temperature"]
         capacity_factors[bus] = np.interp((temp - temp.mean()).values, x, y)
@@ -135,18 +143,20 @@ if __name__ == "__main__":
             clusters=37,
         )
 
-    sustainability_factor = 0.0025
+    egs_config = snakemake.params["sector"]["enhanced_geothermal"]
+    costs_config = snakemake.params["costs"]
+
+    sustainability_factor = egs_config["sustainability_factor"]
     # the share of heat that is replenished from the earth's core.
     # we are not constraining ourselves to the sustainable share, but
     # inversely apply it to our underlying data, which refers to the
-    # sustainable heat.
+    # sustainable heat. Source: Relative magnitude of sustainable heat vs
-
+    # nonsustainable heat in the paper "From hot rock to useful energy..."
-    config = snakemake.config
 
     egs_data = prepare_egs_data(snakemake.input.egs_cost)
 
-    if config["sector"]["enhanced_geothermal_optimism"]:
+    if egs_config["optimism"]:
-        egs_data = egs_data[(year := config["costs"]["year"])]
+        egs_data = egs_data[(year := costs_config["year"])]
         logger.info(
             f"EGS optimism! Building EGS potentials with costs estimated for {year}."
         )

scripts/prepare_sector_network.py

@@ -3582,7 +3582,8 @@ def add_enhanced_geothermal(
             "'add_enhanced_geothermal' not implemented for multiple geothermal nodes."
         )
 
-    config = snakemake.config
+    egs_config = snakemake.params["sector"]["enhanced_geothermal"]
+    costs_config = snakemake.config["costs"]
 
     # matrix defining the overlap between gridded geothermal potential estimation, and bus regions
     overlap = pd.read_csv(egs_overlap, index_col=0)
@@ -3590,7 +3591,7 @@ def add_enhanced_geothermal(
     egs_potentials = pd.read_csv(egs_potentials, index_col=0)
 
     Nyears = n.snapshot_weightings.generators.sum() / 8760
-    dr = config["costs"]["fill_values"]["discount rate"]
+    dr = costs_config["fill_values"]["discount rate"]
     lt = costs.at["geothermal", "lifetime"]
     FOM = costs.at["geothermal", "FOM"]
 
@@ -3673,7 +3674,7 @@ def add_enhanced_geothermal(
         bus_egs["p_nom_max"] = bus_egs["p_nom_max"].multiply(bus_overlap)
         bus_egs = bus_egs.loc[bus_egs.p_nom_max > 0.0]
 
-        if config["sector"]["enhanced_geothermal_performant"]:
+        if egs_config["performant"]:
             bus_egs = bus_egs.sort_values(by="capital_cost").iloc[:1]
             appendix = pd.Index([""])
         else:
@@ -3718,7 +3719,7 @@ def add_enhanced_geothermal(
             p_nom_extendable=True,
             carrier="geothermal organic rankine cycle",
             capital_cost=plant_capital_cost * efficiency_orc,
-            efficiency=efficiency_orc if not as_chp else efficiency_orc * 2.0,
+            efficiency=efficiency_orc,
         )
 
         if as_chp and bus + " urban central heat" in n.buses.index:
@@ -3731,7 +3732,7 @@ def add_enhanced_geothermal(
                 capital_cost=plant_capital_cost
                 * efficiency_orc
                 * costs.at["geothermal", "district heating cost"],
-                efficiency=efficiency_dh * 2.0,
+                efficiency=efficiency_dh,
                 p_nom_extendable=True,
             )
         elif as_chp and not bus + " urban central heat" in n.buses.index:
@@ -3739,15 +3740,13 @@ def add_enhanced_geothermal(
                 efficiency_orc
             )
 
-        if snakemake.params.sector["enhanced_geothermal_flexible"]:
+        if egs_config["flexible"]:
             # this StorageUnit represents flexible operation using the geothermal reservoir.
             # Hence, it is counter-intuitive to install it at the surface bus,
             # this is however the more lean and computationally efficient solution.
 
-            max_hours = snakemake.params.sector[
+            max_hours = egs_config["reservoir_max_hours"]
-                "enhanced_geothermal_reservoir_max_hours"
+            boost = egs_config["reservoir_max_boost"]
-            ]
-            boost = snakemake.params.sector["enhanced_geothermal_reservoir_max_boost"]
 
             n.add(
                 "StorageUnit",
@@ -3912,7 +3911,7 @@ if __name__ == "__main__":
     if options.get("cluster_heat_buses", False) and not first_year_myopic:
         cluster_heat_buses(n)
 
-    if options.get("enhanced_geothermal", False):
+    if options["enhanced_geothermal"].get("enable", False):
         logger.info("Adding Enhanced Geothermal Potential.")
         add_enhanced_geothermal(
             n, snakemake.input["egs_potentials"], snakemake.input["egs_overlap"], costs

scripts/solve_network.py

@@ -816,7 +816,7 @@ def add_geothermal_chp_constraint(n):
 def add_flexible_egs_constraint(n):
     well_index = n.links.loc[n.links.carrier == "geothermal heat"].index
     storage_index = n.storage_units.loc[
-        n.storage_units.carrier == "geothermal heat  "
+        n.storage_units.carrier == "geothermal heat"
     ].index
 
     p_nom_rhs = n.model["Link-p_nom"].loc[well_index]