diff --git a/config/config.default.yaml b/config/config.default.yaml index 62959b6324..7c3183a234 100644 --- a/config/config.default.yaml +++ b/config/config.default.yaml @@ -419,6 +419,13 @@ renewable: eia_norm_year: false eia_correct_by_capacity: false eia_approximate_missing: false + # land-eligibility settings for alternative carbon dioxide removal technologies, + # feeding determine_availability_matrix.py the same way onwind/solar do above + biochar: + corine: [12, 13, 14] # CORINE Land Cover codes for non-irrigated arable land (potential biochar feedstock land) + natura: false + cutout: default + excluder_resolution: 250 # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#conventional conventional: @@ -837,6 +844,9 @@ sector: methanation: true coal_cc: false dac: true + biochar: + enable: false + heat_output: false co2_vent: false heat_vent: urban central: true @@ -1070,6 +1080,11 @@ industry: reference_year: 2023 oil_refining_emissions: 0.013 +biochar: + application_per_sqkm: 1500 # tonnes of biochar per km² (≈15 t/ha) + max_land_usage: 0.2 # maximum fraction of suitable CORINE area used for biochar + number_years: 25 # amortisation period (years) for biochar CO2 storage + # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#costs costs: year: 2050 diff --git a/config/plotting.default.yaml b/config/plotting.default.yaml index e5156d485a..f796cd2532 100644 --- a/config/plotting.default.yaml +++ b/config/plotting.default.yaml @@ -738,3 +738,5 @@ plotting: import NH3: '#e2ed74' import oil: '#93eda2' import methanol: '#87d0e6' + biochar heat: '#228B22' + co2 biochar: '#005000' diff --git a/rules/build_sector.smk b/rules/build_sector.smk index 59b6540fea..4b5502eef8 100755 --- a/rules/build_sector.smk +++ b/rules/build_sector.smk @@ -920,6 +920,58 @@ rule build_biomass_potentials: scripts("build_biomass_potentials.py") +rule determine_carbon_dioxide_removal_availability_matrix: + input: + corine=ancient(rules.retrieve_corine.output["tif_file"]), + regions=resources("regions_onshore_base_s_{clusters}.geojson"), + cutout=lambda w: input_cutout( + w, config_provider("renewable", w.technology, "cutout")(w) + ), + output: + resources( + "availability_matrix_carbon_dioxide_removal_{clusters}_{technology}.nc" + ), + log: + logs( + "determine_carbon_dioxide_removal_availability_matrix_{clusters}_{technology}.log" + ), + wildcard_constraints: + technology="biochar", + threads: config["atlite"].get("nprocesses", 4) + resources: + mem_mb=config["atlite"].get("nprocesses", 4) * 5000, + params: + renewable=config_provider("renewable"), + plot_availability_matrix=config_provider("atlite", "plot_availability_matrix"), + message: + "Determining availability matrix for {wildcards.clusters} clusters and {wildcards.technology} carbon dioxide removal technology" + script: + scripts("determine_availability_matrix.py") + + +rule build_available_land: + input: + availability_matrix=resources( + "availability_matrix_carbon_dioxide_removal_{clusters}_{technology}.nc" + ), + regions=resources("regions_onshore_base_s_{clusters}.geojson"), + cutout=lambda w: input_cutout( + w, config_provider("renewable", w.technology, "cutout")(w) + ), + output: + csv_file=resources("{technology}_available_land_s_{clusters}.csv"), + log: + logs("build_available_land_{clusters}_{technology}.log"), + wildcard_constraints: + technology="biochar", + resources: + mem_mb=8000, + params: + renewable=config_provider("renewable"), + script: + scripts("build_available_land.py") + + rule build_biomass_transport_costs: input: sc1="data/biomass_transport_costs_supplychain1.csv", @@ -1689,6 +1741,11 @@ rule prepare_sector_network: if config_provider("sector", "district_heating", "ates", "enable")(w) else [] ), + biochar_potentials=lambda w: ( + resources("biochar_available_land_s_{clusters}.csv") + if config_provider("sector", "biochar", "enable")(w) + else [] + ), output: resources( "networks/base_s_{clusters}_{opts}_{sector_opts}_{planning_horizons}.nc" diff --git a/scripts/build_available_land.py b/scripts/build_available_land.py new file mode 100644 index 0000000000..a08f3d903b --- /dev/null +++ b/scripts/build_available_land.py @@ -0,0 +1,71 @@ +# SPDX-FileCopyrightText: Contributors to PyPSA-Eur +# +# SPDX-License-Identifier: MIT +""" +Aggregate a per-bus land-eligibility availability matrix (as produced by +``determine_availability_matrix.py`` / ``determine_rock_weathering_availability_matrix.py``) +into eligible land area per network node, in km². + +Shared by biochar, afforestation and rock weathering: the per-bus fractional +eligibility (0-1 per cutout grid cell) is multiplied by the cutout grid cell +area and summed per bus, recovering the same eligible-area-in-km² number the +previous per-technology scripts (``build_corine_potentials.py``, +``build_EW_potentials.py``) each computed independently via +``atlite.gis.shape_availability``. + +Output: CSV with columns ``area [sqkm]`` (total node area) and +``potential [sqkm]`` (eligible area after exclusions), indexed by node name. +The downstream consumers (``add_biochar``, ``build_afforestation_potentials``, +``add_rock_weathering``) apply their own technology-specific rate (t/km², +growth rate, etc.) to the ``potential [sqkm]`` column. +""" + +import logging + +import geopandas as gpd +import pandas as pd +import xarray as xr + +from scripts._helpers import configure_logging, load_cutout, set_scenario_config + +logger = logging.getLogger(__name__) + + +if __name__ == "__main__": + if "snakemake" not in globals(): + from scripts._helpers import mock_snakemake + + snakemake = mock_snakemake( + "build_available_land", clusters="39", technology="biochar" + ) + configure_logging(snakemake) + set_scenario_config(snakemake) + + technology = snakemake.wildcards.technology + + availability = xr.open_dataarray(snakemake.input.availability_matrix) + cutout = load_cutout(snakemake.input.cutout) + regions = ( + gpd.read_file(snakemake.input.regions).set_index("name").rename_axis("bus") + ) + + cell_area = cutout.grid.to_crs(3035).area / 1e6 + cell_area = xr.DataArray( + cell_area.values.reshape(cutout.shape), + [cutout.coords["y"], cutout.coords["x"]], + ) + + potential = (availability * cell_area).sum(dim=["x", "y"]).to_pandas() + area = regions.geometry.to_crs(3035).area / 1e6 + + df = pd.DataFrame({"area [sqkm]": area, "potential [sqkm]": potential}) + df.index.name = "node" + + logger.info( + "Total %s potential: %.0f sqkm (of %.0f sqkm total node area)", + technology, + df["potential [sqkm]"].sum(), + df["area [sqkm]"].sum(), + ) + + df.to_csv(snakemake.output.csv_file) diff --git a/scripts/prepare_sector_network.py b/scripts/prepare_sector_network.py index 91efa0bc0d..95dbedd3e9 100755 --- a/scripts/prepare_sector_network.py +++ b/scripts/prepare_sector_network.py @@ -1247,6 +1247,150 @@ def add_dac(n, costs): ) +def add_biochar(n, costs): + """ + Add biochar (CDR via biomass pyrolysis and stable-carbon soil storage) + to the network as Bus, Store, and Link components, all sharing a single + "co2 biochar" carrier; optionally exports pyrolysis waste heat to urban + central heating. + + A single pyrolysis Link per node converts solid biomass (bus2) and + electricity (bus3) into stored biochar-carbon (bus1, drawn from the + atmosphere via bus0) and, if enabled, district heat (bus4, with an + overflow Store for heat that cannot be absorbed by the local heat + network). The store's capacity is dimensioned from the eligible CORINE + land area (see ``determine_availability_matrix.py`` and + ``build_available_land.py``) multiplied by a per-km2 application rate + and the stable-carbon fraction (``co2_per_tonne``), amortised over an + assumed storage horizon. + + Parameters + ---------- + n : pypsa.Network + The PyPSA network container object + costs : pd.DataFrame + Costs and parameters for different technologies. Must contain a + 'biochar pyrolysis' entry with 'biomass-input', 'yield-biochar', + 'electricity-input', 'heat-output', 'capital_cost', and 'VOM' + parameters + + Returns + ------- + None + Modifies the network object in-place by adding the biochar Bus, + Store, and Link (and, if enabled, heat-output Bus/Link/Store) + + Notes + ----- + Reads ``snakemake.input.biochar_potentials`` (eligible area in km2 per + node) and ``snakemake.config["biochar"]["application_per_sqkm"]`` / + ``["max_land_usage"]`` / ``["number_years"]``. Heat output is gated by + ``sector.heating`` and ``sector.biochar.heat_output``. + """ + logger.info("Adding biochar.") + + biochar_potentials = pd.read_csv(snakemake.input.biochar_potentials).set_index( + "node" + ) + + n.add("Carrier", "co2 biochar") + + n.add( + "Bus", + spatial.nodes + " co2 biochar", + location=spatial.nodes, + carrier="co2 biochar", + unit="t_co2", + ) + + co2_per_tonne = ( + 1 + / costs.at["biochar pyrolysis", "biomass-input"] + * 1 + / costs.at["biochar pyrolysis", "yield-biochar"] + ) # tCO2 / t_biochar + + n.add( + "Store", + spatial.nodes + " co2 biochar", + bus=spatial.nodes + " co2 biochar", + carrier="co2 biochar", + e_nom_extendable=True, + e_nom_max=( + biochar_potentials["potential [sqkm]"].values + * co2_per_tonne + * snakemake.config["biochar"]["application_per_sqkm"] + * snakemake.config["biochar"]["max_land_usage"] + / snakemake.config["biochar"]["number_years"] + ), + ) + + if len(spatial.biomass.nodes) == 1: + biomass_buses = spatial.biomass.nodes + else: + biomass_buses = spatial.nodes + " solid biomass" + + # heat output into urban central heat system + if ( + snakemake.config["sector"]["heating"] + and snakemake.config["sector"]["biochar"]["heat_output"] + ): + logger.info("Adding biochar heat output to urban central heat system.") + biochar_heat_buses = [] + for node in spatial.nodes: + if (node + " urban central heat") in n.buses.index: + biochar_heat_bus = node + " biochar heat" + biochar_heat_buses.append(biochar_heat_bus) + n.add("Bus", biochar_heat_bus, carrier="biochar heat") + n.add( + "Link", + biochar_heat_bus, + bus0=biochar_heat_bus, + bus1=node + " urban central heat", + p_nom_extendable=True, + carrier="biochar heat", + ) + biochar_heat_bus_waste = node + " biochar heat waste" + n.add("Bus", biochar_heat_bus_waste, carrier="biochar heat") + n.add( + "Store", + biochar_heat_bus_waste, + bus=biochar_heat_bus_waste, + e_nom_extendable=True, + carrier="biochar heat", + ) + n.add( + "Link", + biochar_heat_bus_waste, + bus0=biochar_heat_bus, + bus1=biochar_heat_bus_waste, + p_nom_extendable=True, + carrier="biochar heat", + ) + else: + biochar_heat_buses.append(None) + else: + biochar_heat_buses = [None] + + n.add( + "Link", + spatial.nodes + " biochar", + bus0="co2 atmosphere", + bus1=spatial.nodes + " co2 biochar", + bus2=biomass_buses, + bus3=spatial.nodes, + bus4=biochar_heat_buses, + carrier="co2 biochar", + capital_cost=costs.at["biochar pyrolysis", "capital_cost"], + marginal_cost=costs.at["biochar pyrolysis", "VOM"], + efficiency=1.0, + efficiency2=-costs.at["biochar pyrolysis", "biomass-input"], + efficiency3=-costs.at["biochar pyrolysis", "electricity-input"], + efficiency4=costs.at["biochar pyrolysis", "heat-output"], + p_nom_extendable=True, + ) + + def add_co2limit(n, options, co2_totals_file, countries, nyears, limit): """ Add a global CO2 emissions constraint to the network. @@ -6518,6 +6662,9 @@ def add_import_options( if options["dac"]: add_dac(n, costs) + if options.get("biochar", {}).get("enable"): + add_biochar(n, costs) + if not options["electricity_transmission_grid"]: decentral(n)