diff --git a/Snakefile b/Snakefile index 90fcfb56..cf829322 100644 --- a/Snakefile +++ b/Snakefile @@ -280,6 +280,23 @@ else: build_biomass_transport_costs_output = {} +if config["sector"].get("sequestration_potential", False): + rule build_sequestration_potentials: + input: + sequestration_potential="data/complete_map_2020_unit_Mt.geojson", + regions_onshore=pypsaeur("resources/regions_onshore_elec_s{simpl}_{clusters}.geojson"), + regions_offshore=pypsaeur("resources/regions_offshore_elec_s{simpl}_{clusters}.geojson"), + output: + sequestration_potential="resources/co2_sequestration_potential_elec_s{simpl}_{clusters}.csv" + threads: 1 + resources: mem_mb=4000 + benchmark: "benchmarks/build_sequestration_potentials_s{simpl}_{clusters}" + script: "scripts/build_sequestration_potentials.py" + build_sequestration_potentials_output = rules.build_sequestration_potentials.output +else: + build_sequestration_potentials_output = {} + + rule build_salt_cavern_potentials: input: salt_caverns="data/h2_salt_caverns_GWh_per_sqkm.geojson", @@ -512,7 +529,8 @@ rule prepare_sector_network: solar_thermal_rural="resources/solar_thermal_rural_elec_s{simpl}_{clusters}.nc" if config["sector"]["solar_thermal"] else [], **build_retro_cost_output, **build_biomass_transport_costs_output, - **gas_infrastructure + **gas_infrastructure, + **build_sequestration_potentials_output output: RDIR + '/prenetworks/elec_s{simpl}_{clusters}_lv{lv}_{opts}_{sector_opts}_{planning_horizons}.nc' threads: 1 resources: mem_mb=2000 diff --git a/config.default.yaml b/config.default.yaml index ee1c5059..5fa072d4 100644 --- a/config.default.yaml +++ b/config.default.yaml @@ -250,9 +250,11 @@ sector: dac: true co2_vent: false SMR: true + sequestration_potential: true # geological co2 storage potential co2_sequestration_potential: 200 #MtCO2/a sequestration potential for Europe co2_sequestration_cost: 10 #EUR/tCO2 for sequestration of CO2 - co2_network: false + co2_spatial: false + co2network: false cc_fraction: 0.9 # default fraction of CO2 captured with post-combustion capture hydrogen_underground_storage: true hydrogen_underground_storage_locations: @@ -275,6 +277,7 @@ sector: gas_network_connectivity_upgrade: 1 # https://networkx.org/documentation/stable/reference/algorithms/generated/networkx.algorithms.connectivity.edge_augmentation.k_edge_augmentation.html#networkx.algorithms.connectivity.edge_augmentation.k_edge_augmentation gas_distribution_grid: true gas_distribution_grid_cost_factor: 1.0 #multiplies cost in data/costs.csv + biomass_spatial: false # biomass transport between nodes biomass_transport: false # biomass transport between nodes conventional_generation: # generator : carrier OCGT: gas diff --git a/scripts/build_sequestration_potentials.py b/scripts/build_sequestration_potentials.py new file mode 100644 index 00000000..cf8ab6d6 --- /dev/null +++ b/scripts/build_sequestration_potentials.py @@ -0,0 +1,48 @@ +import pandas as pd +import geopandas as gpd + +def area(gdf): + """Returns area of GeoDataFrame geometries in square kilometers.""" + return gdf.to_crs(epsg=3035).area.div(1e6) + + +def allocate_sequestration_potential(gdf, regions, attr='conservative estimate Mt', threshold=3): + gdf = gdf.loc[gdf[attr] > threshold, [attr, "geometry"]] + gdf["area_sqkm"] = area(gdf) + overlay = gpd.overlay(regions, gdf, keep_geom_type=True) + overlay["share"] = area(overlay) / overlay["area_sqkm"] + adjust_cols = overlay.columns.difference({"name", "area_sqkm", "geometry", "share"}) + overlay[adjust_cols] = overlay[adjust_cols].multiply(overlay["share"], axis=0) + gdf_regions = overlay.groupby("name").sum() + gdf_regions.drop(["area_sqkm", "share"], axis=1, inplace=True) + return gdf_regions.squeeze() + + +if __name__ == "__main__": + if 'snakemake' not in globals(): + from helper import mock_snakemake + snakemake = mock_snakemake( + 'build_sequestration_potentials', + simpl='', + clusters="181" + ) + + # TODO move to config.yaml + threshold = 3 + include_onshore = False + + gdf = gpd.read_file(snakemake.input.sequestration_potential) + + regions = gpd.read_file(snakemake.input.regions_offshore) + if include_onshore: + onregions = gpd.read_file(snakemake.input.regions_onshore) + regions = pd.concat([regions, onregions]).dissolve(by='name').reset_index() + + attr = snakemake.config['sector']["sequestration_potential"] + kwargs = dict(attr=attr, threshold=threshold) if isinstance(attr, str) else {} + + s = allocate_sequestration_potential(gdf, regions, **kwargs) + + s = s.where(s>threshold).dropna() + + s.to_csv(snakemake.output.sequestration_potential) diff --git a/scripts/helper.py b/scripts/helper.py index e6ddfd4a..62ae33c0 100644 --- a/scripts/helper.py +++ b/scripts/helper.py @@ -138,6 +138,6 @@ def parse(l): def update_config_with_sector_opts(config, sector_opts): for o in sector_opts.split("-"): - if o.startswith("CF:"): + if o.startswith("CF+"): l = o.split("+")[1:] update_config(config, parse(l)) \ No newline at end of file diff --git a/scripts/prepare_sector_network.py b/scripts/prepare_sector_network.py index b6c052be..76b71fc8 100644 --- a/scripts/prepare_sector_network.py +++ b/scripts/prepare_sector_network.py @@ -44,7 +44,7 @@ def define_spatial(nodes, options): spatial.biomass = SimpleNamespace() - if options["biomass_transport"]: + if options.get("biomass_spatial", options["biomass_transport"]): spatial.biomass.nodes = nodes + " solid biomass" spatial.biomass.locations = nodes spatial.biomass.industry = nodes + " solid biomass for industry" @@ -61,7 +61,7 @@ def define_spatial(nodes, options): spatial.co2 = SimpleNamespace() - if options["co2_network"]: + if options.get("co2_spatial", options["co2network"]): spatial.co2.nodes = nodes + " co2 stored" spatial.co2.locations = nodes spatial.co2.vents = nodes + " co2 vent" @@ -88,8 +88,11 @@ def define_spatial(nodes, options): spatial.gas.locations = ["EU"] spatial.gas.biogas = ["EU biogas"] spatial.gas.industry = ["gas for industry"] - spatial.gas.industry_cc = ["gas for industry CC"] spatial.gas.biogas_to_gas = ["EU biogas to gas"] + if options.get("co2_spatial", options["co2network"]): + spatial.gas.industry_cc = nodes + " gas for industry CC" + else: + spatial.gas.industry_cc = ["gas for industry CC"] spatial.gas.df = pd.DataFrame(vars(spatial.gas), index=nodes) @@ -507,10 +510,18 @@ def add_co2_tracking(n, options): unit="t_co2" ) + if options["sequestration_potential"]: + # TODO make configurable options + upper_limit = 25000 # Mt + annualiser = 25 # TODO research suitable value + e_nom_max = pd.read_csv(snakemake.input.sequestration_potential, index_col=0).squeeze() + e_nom_max = e_nom_max.reindex(spatial.co2.locations).fillna(0.).clip(upper=upper_limit).mul(1e6) / annualiser # t + e_nom_max = e_nom_max.rename(index=lambda x: x + " co2 stored") + n.madd("Store", spatial.co2.nodes, e_nom_extendable=True, - e_nom_max=np.inf, + e_nom_max=e_nom_max, capital_cost=options['co2_sequestration_cost'], carrier="co2 stored", bus=spatial.co2.nodes @@ -1218,7 +1229,7 @@ def add_storage_and_grids(n, costs): bus0=spatial.coal.nodes, bus1=spatial.nodes, bus2="co2 atmosphere", - bus3="co2 stored", + bus3=spatial.co2.nodes, marginal_cost=costs.at['coal', 'efficiency'] * costs.at['coal', 'VOM'], #NB: VOM is per MWel capital_cost=costs.at['coal', 'efficiency'] * costs.at['coal', 'fixed'] + costs.at['biomass CHP capture', 'fixed'] * costs.at['coal', 'CO2 intensity'], #NB: fixed cost is per MWel p_nom_extendable=True, @@ -1828,7 +1839,7 @@ def add_biomass(n, costs): else: biogas_potentials_spatial = biomass_potentials["biogas"].sum() - if options["biomass_transport"]: + if options.get("biomass_spatial", options["biomass_transport"]): solid_biomass_potentials_spatial = biomass_potentials["solid biomass"].rename(index=lambda x: x + " solid biomass") else: solid_biomass_potentials_spatial = biomass_potentials["solid biomass"].sum() @@ -2052,7 +2063,7 @@ def add_industry(n, costs): unit="MWh_LHV" ) - if options["biomass_transport"]: + if options.get("biomass_spatial", options["biomass_transport"]): p_set = industrial_demand.loc[spatial.biomass.locations, "solid biomass"].rename(index=lambda x: x + " solid biomass for industry") / 8760 else: p_set = industrial_demand["solid biomass"].sum() / 8760 @@ -2775,7 +2786,7 @@ if __name__ == "__main__": if "noH2network" in opts: remove_h2_network(n) - if options["co2_network"]: + if options["co2network"]: add_co2_network(n, costs) solver_name = snakemake.config["solving"]["solver"]["name"]