From e6203f51cdf8f3034ee2077898be12d0ffb16b20 Mon Sep 17 00:00:00 2001 From: LukasFrankenQ Date: Mon, 4 Sep 2023 21:07:53 +0100 Subject: [PATCH] pre-testrun version added to prepare_sector_network --- rules/build_sector.smk | 2 + scripts/build_egs_potentials.py | 22 ++++-- scripts/prepare_sector_network.py | 126 ++++++++++++++++++++++++++++++ 3 files changed, 144 insertions(+), 6 deletions(-) diff --git a/rules/build_sector.smk b/rules/build_sector.smk index 4655ebba..25efd7db 100644 --- a/rules/build_sector.smk +++ b/rules/build_sector.smk @@ -783,6 +783,8 @@ 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", solar_thermal_total=RESOURCES + "solar_thermal_total_elec_s{simpl}_{clusters}.nc" if config["sector"]["solar_thermal"] diff --git a/scripts/build_egs_potentials.py b/scripts/build_egs_potentials.py index aecfd72a..61feee92 100644 --- a/scripts/build_egs_potentials.py +++ b/scripts/build_egs_potentials.py @@ -23,7 +23,7 @@ import json import pandas as pd import geopandas as gpd -from shapely.geometry import Polygon, Point +from shapely.geometry import Polygon def prepare_egs_data(egs_file): @@ -108,13 +108,23 @@ if __name__ == "__main__": egs_data.index = egs_data.geometry.astype(str) egs_shapes = egs_data.geometry - network_shapes = gpd.read_file(snakemake.input.shapes).set_index("name", drop=True).set_crs(epsg=4326) - egs_shapes = egs_data.geometry + network_shapes = ( + gpd.read_file(snakemake.input.shapes) + .set_index("name", drop=True) + .set_crs(epsg=4326) + ) + + overlap_matrix = ( + pd.DataFrame( + index=network_shapes.index, + columns=(egs_shapes := egs_data.geometry).astype(str).values) + ) - overlap_matrix = pd.DataFrame(index=network_shapes.index, columns=egs_shapes.index) - for name, polygon in network_shapes.geometry.items(): - overlap_matrix.loc[name] = egs_shapes.intersection(polygon).area / egs_shapes.area + overlap_matrix.loc[name] = ( + egs_shapes + .intersection(polygon).area + ) / egs_shapes.area overlap_matrix.to_csv(snakemake.output["egs_overlap"]) diff --git a/scripts/prepare_sector_network.py b/scripts/prepare_sector_network.py index 11406bff..de023f30 100644 --- a/scripts/prepare_sector_network.py +++ b/scripts/prepare_sector_network.py @@ -3280,6 +3280,121 @@ def set_temporal_aggregation(n, opts, solver_name): return n +def add_egs_potential( + n, + egs_potentials, + egs_overlap, + costs, +): + """ + Adds EGS potential to model. + + Built in scripts/build_egs_potentials.py + """ + + n.add( + "Carrier", + "geothermal heat", + nice_name="Geothermal Heat", + color=snakemake.config["plotting"]["tech_colors"]["geothermal"], + co2_emissions=costs.at["geothermal", "CO2 intensity"], + ) + + config = snakemake.config + + overlap = pd.read_csv(egs_overlap, index_col=0) + indicator_matrix = ( + pd.DataFrame( + np.ceil(overlap.values), + index=overlap.index, + columns=overlap.columns + ) + ) + + egs_potentials = pd.read_csv(egs_potentials, index_col=0) + + Nyears = n.snapshot_weightings.generators.sum() / 8760 + dr = config["costs"]["fill_values"]["discount rate"] + lt = costs.at["geothermal", "lifetime"] + + egs_annuity = calculate_annuity(lt, dr) + + egs_potentials["capital_cost"] = ( + (egs_annuity + 0.02 / 1.02) * egs_potentials["CAPEX"] * Nyears + ) + + efficiency = costs.at["geothermal", "efficiency electricity"] + + # p_nom conversion GW -> MW + # capital_cost conversion Euro/kW -> Euro/MW + egs_potentials["p_nom_max"] = egs_potentials["p_nom_max"] * 1000.0 + egs_potentials["capital_cost"] = egs_potentials["capital_cost"] * 1000.0 + + n.add( + "Bus", + "EU geothermal heat", + carrier="geothermal heat", + unit="MWh_th", + location="EU", + ) + + n.add( + "Generator", + "EU geothermal heat", + bus="EU geothermal heat", + carrier="geothermal heat", + p_nom_max=np.inf, + p_nom_extendable=True, + ) + + for bus, overlaps in overlap.iterrows(): + if not overlaps.sum(): + continue + + overlaps.index = np.arange(len(overlaps)) + + overlaps = overlaps.loc[overlaps > 0.] + indicators = indicator_matrix.loc[bus] + indicators.index = np.arange(len(indicators)) + + bus_egs = egs_potentials.loc[overlaps.loc[overlaps > 0.].index] + + if not len(bus_egs): + continue + + bus_egs["p_nom_max"] = bus_egs["p_nom_max"].multiply(overlaps) + bus_egs = bus_egs.loc[bus_egs.p_nom_max > 0.] + + if config["sector"]["enhanced_geothermal_best_only"]: + bus_egs = bus_egs.sort_values(by="capital_cost").iloc[:1] + appendix = "" + else: + appendix = " " + bus_egs.index + + bus_egs.index = np.arange(len(bus_egs)).astype(str) + + p_nom_max = bus_egs["p_nom_max"] + p_nom_max.index = f"{bus} enhanced geothermal" + appendix + + capital_cost = bus_egs["capital_cost"] + capital_cost.index = f"{bus} enhanced geothermal" + appendix + + n.madd( + "Link", + f"{bus} enhanced geothermal" + appendix, + bus0="EU geothermal heat", + bus1=bus, + bus2="co2 atmosphere", + carrier="geothermal heat", + p_nom_extendable=True, + p_nom_max=p_nom_max / efficiency, + capital_cost=capital_cost * efficiency, + efficiency=efficiency, + efficiency2=costs.at["geothermal", "CO2 intensity"] * efficiency, + ) + + + if __name__ == "__main__": if "snakemake" not in globals(): from _helpers import mock_snakemake @@ -3453,6 +3568,17 @@ if __name__ == "__main__": if options.get("cluster_heat_buses", False) and not first_year_myopic: cluster_heat_buses(n) + if options.get("enhanced_geothermal"): + logger.info("Adding Enhanced Geothermal Potential.") + + add_egs_potential( + n, + snakemake.input["egs_potentials"], + snakemake.input["egs_overlap"], + costs + ) + + n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards))) sanitize_carriers(n, snakemake.config)