diff --git a/borg-it b/borg-it old mode 100755 new mode 100644 diff --git a/config/config.default.yaml b/config/config.default.yaml index 6ba8039e..83371c37 100644 --- a/config/config.default.yaml +++ b/config/config.default.yaml @@ -67,7 +67,6 @@ snapshots: # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#enable enable: retrieve: auto - prepare_links_p_nom: false retrieve_databundle: true retrieve_cost_data: true build_cutout: false @@ -355,7 +354,6 @@ biomass: - Secondary Forestry residues - woodchips - Sawdust - Residues from landscape care - - Municipal waste not included: - Sugar from sugar beet - Rape seed @@ -369,6 +367,25 @@ biomass: biogas: - Manure solid, liquid - Sludge + municipal solid waste: + - Municipal waste + share_unsustainable_use_retained: + 2020: 1 + 2025: 0.66 + 2030: 0.33 + 2035: 0 + 2040: 0 + 2045: 0 + 2050: 0 + share_sustainable_potential_available: + 2020: 0 + 2025: 0.33 + 2030: 0.66 + 2035: 1 + 2040: 1 + 2045: 1 + 2050: 1 + # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#solar-thermal solar_thermal: @@ -409,6 +426,22 @@ sector: 2045: 0.8 2050: 1.0 district_heating_loss: 0.15 + forward_temperature: 90 #C + return_temperature: 50 #C + heat_source_cooling: 6 #K + heat_pump_cop_approximation: + refrigerant: ammonia + heat_exchanger_pinch_point_temperature_difference: 5 #K + isentropic_compressor_efficiency: 0.8 + heat_loss: 0.0 + heat_pump_sources: + urban central: + - air + urban decentral: + - air + rural: + - air + - ground cluster_heat_buses: true heat_demand_cutout: default bev_dsm_restriction_value: 0.75 @@ -491,7 +524,7 @@ sector: aviation_demand_factor: 1. HVC_demand_factor: 1. time_dep_hp_cop: true - heat_pump_sink_T: 55. + heat_pump_sink_T_individual_heating: 55. reduce_space_heat_exogenously: true reduce_space_heat_exogenously_factor: 2020: 0.10 # this results in a space heat demand reduction of 10% @@ -613,6 +646,7 @@ sector: biomass_to_liquid: false electrobiofuels: false biosng: false + municipal_solid_waste: false limit_max_growth: enable: false # allowing 30% larger than max historic growth @@ -640,6 +674,7 @@ sector: max_amount: 1390 # TWh upstream_emissions_factor: .1 #share of solid biomass CO2 emissions at full combustion + # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#industry industry: St_primary_fraction: @@ -733,7 +768,7 @@ industry: # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#costs costs: year: 2030 - version: v0.9.0 + version: v0.9.1 social_discountrate: 0.02 fill_values: FOM: 0 @@ -1036,6 +1071,7 @@ plotting: biogas: '#e3d37d' biomass: '#baa741' solid biomass: '#baa741' + municipal solid waste: '#91ba41' solid biomass import: '#d5ca8d' solid biomass transport: '#baa741' solid biomass for industry: '#7a6d26' @@ -1050,6 +1086,7 @@ plotting: services rural biomass boiler: '#c6cf98' services urban decentral biomass boiler: '#dde5b5' biomass to liquid: '#32CD32' + unsustainable bioliquids: '#32CD32' electrobiofuels: 'red' BioSNG: '#123456' # power transmission diff --git a/data/links_p_nom.csv b/data/links_p_nom.csv index bd7a4c95..56a99e52 100644 --- a/data/links_p_nom.csv +++ b/data/links_p_nom.csv @@ -5,7 +5,7 @@ Cross-Channel,France - Echingen 50°41′48″N 1°38′21″E / 50.69667 Volgograd-Donbass,Russia - Volzhskaya 48°49′34″N 44°40′20″E / 48.82611°N 44.67222°E,Ukraine - Mikhailovskaya 48°39′13″N 38°33′56″E / 48.65361°N 38.56556°E,475(0/475),400,750.0,1965,Merc/Thyr,Shut down in 2014,[1],44.672222222222224,48.82611111111111,38.565555555555555,48.65361111111111 Konti-Skan 1,Denmark - Vester Hassing 57°3′46″N 10°5′24″E / 57.06278°N 10.09000°E,Sweden - Stenkullen 57°48′15″N 12°19′13″E / 57.80417°N 12.32028°E,176(87/89),250,250.0,1965,Merc,Replaced in August 2006 by modern converters using thyristors,[1],10.09,57.062777777777775,12.320277777777777,57.80416666666667 SACOI 1a,Italy - Suvereto 43°3′10″N 10°41′42″E / 43.05278°N 10.69500°E ( before 1992: Italy - San Dalmazio 43°15′43″N 10°55′05″E / 43.26194°N 10.91806°E),"France- Lucciana 42°31′40″N 9°26′59″E / 42.52778°N 9.44972°E",483(365/118),200,200.0,1965,Merc,"Replaced in 1986 by Thyr- multiterminal scheme",[1],10.695,43.05277777777778,9.449722222222222,42.52777777777778 -SACOI 1b,"France- Lucciana 42°31′40″N 9°26′59″E / 42.52778°N 9.44972°E", "Codrongianos- Italy 40°39′7″N 8°42′48″E / 40.65194°N 8.71333°E",483(365/118),200,200.0,1965,Merc,"Replaced in 1986 by Thyr- multiterminal scheme",[1],9.449722222222222,42.52777777777778,8.679351,40.65765 +SACOI 1b,"France- Lucciana 42°31′40″N 9°26′59″E / 42.52778°N 9.44972°E","Codrongianos- Italy 40°39′7″N 8°42′48″E / 40.65194°N 8.71333°E",483(365/118),200,200.0,1965,Merc,"Replaced in 1986 by Thyr- multiterminal scheme",[1],9.449722222222222,42.52777777777778,8.679351,40.65765 Kingsnorth,UK - Kingsnorth 51°25′11″N 0°35′46″E / 51.41972°N 0.59611°E,UK - London-Beddington 51°22′23″N 0°7′38″W / 51.37306°N 0.12722°W,85(85/0),266,320.0,1975,Merc,Bipolar scheme Supplier: English Electric Shut down in 1987,[33],0.5961111111111111,51.41972222222222,-0.1272222222222222,51.37305555555555 Skagerrak 1 + 2,Denmark - Tjele 56°28′44″N 9°34′1″E / 56.47889°N 9.56694°E,Norway - Kristiansand 58°15′36″N 7°53′55″E / 58.26000°N 7.89861°E,230(130/100),250,500.0,1977,Thyr,Supplier: STK(Nexans) Control system upgrade by ABB in 2007,[34][35][36],9.566944444444445,56.47888888888889,7.898611111111111,58.26 Gotland 2,Sweden - Västervik 57°43′41″N 16°38′51″E / 57.72806°N 16.64750°E,Sweden - Yigne 57°35′13″N 18°11′44″E / 57.58694°N 18.19556°E,99.5(92.9/6.6),150,130.0,1983,Thyr,Supplier: ABB,,16.6475,57.72805555555556,18.195555555555554,57.58694444444444 @@ -23,7 +23,7 @@ Visby-Nas,Sweden - Nas 57°05′58″N 18°14′27″E / 57.09944°N 18.24 SwePol,Poland - Wierzbięcin 54°30′8″N 16°53′28″E / 54.50222°N 16.89111°E,Sweden - Stärnö 56°9′11″N 14°50′29″E / 56.15306°N 14.84139°E,245(245/0),450,600.0,2000,Thyr,Supplier: ABB,[38],16.891111111111112,54.50222222222222,14.841388888888888,56.153055555555554 Tjæreborg,Denmark - Tjæreborg/Enge 55°26′52″N 8°35′34″E / 55.44778°N 8.59278°E,Denmark - Tjæreborg/Substation 55°28′07″N 8°33′36″E / 55.46861°N 8.56000°E,4.3(4.3/0),9,7.0,2000,IGBT,Interconnection to wind power generating stations,,8.592777777777778,55.44777777777778,8.56,55.46861111111111 Italy-Greece,Greece - Arachthos 39°11′00″N 20°57′48″E / 39.18333°N 20.96333°E,Italy - Galatina 40°9′53″N 18°7′49″E / 40.16472°N 18.13028°E,310(200/110),400,500.0,2001,Thyr,,,20.963333333333335,39.18333333333333,18.130277777777778,40.164722222222224 -Moyle,UK - Auchencrosh 55°04′10″N 4°58′50″W / 55.06944°N 4.98056°W,UK - N. Ireland- Ballycronan More 54°50′34″N 5°46′11″W / 54.84278°N 5.76972°W,63.5(63.5/0),250,2501.0,2001,Thyr,"Supplier: Siemens- Nexans",[39],-4.980555555555556,55.06944444444444,-5.769722222222223,54.842777777777776 +Moyle,UK - Auchencrosh 55°04′10″N 4°58′50″W / 55.06944°N 4.98056°W,UK - N. Ireland- Ballycronan More 54°50′34″N 5°46′11″W / 54.84278°N 5.76972°W,63.5(63.5/0),250,500.0,2001,Thyr,"Supplier: Siemens- Nexans",[39],-4.980555555555556,55.06944444444444,-5.769722222222223,54.842777777777776 HVDC Troll,Norway - Kollsnes 60°33′01″N 4°50′26″E / 60.55028°N 4.84056°E,Norway - Offshore platform Troll A 60°40′00″N 3°40′00″E / 60.66667°N 3.66667°E,70(70/0),60,80.0,2004,IGBT,Power supply for offshore gas compressor Supplier: ABB,[40],4.8405555555555555,60.55027777777778,3.6666666666666665,60.666666666666664 Estlink,Finland - Espoo 60°12′14″N 24°33′06″E / 60.20389°N 24.55167°E,Estonia - Harku 59°23′5″N 24°33′37″E / 59.38472°N 24.56028°E,105(105/0),150,350.0,2006,IGBT,Supplier: ABB,[40],24.551666666666666,60.20388888888889,24.560277777777777,59.38472222222222 NorNed,Netherlands - Eemshaven 53°26′4″N 6°51′57″E / 53.43444°N 6.86583°E,Norway - Feda 58°16′58″N 6°51′55″E / 58.28278°N 6.86528°E,580(580/0),450,700.0,2008,Thyr,"Supplier: ABB- Nexans",[40],6.865833333333334,53.434444444444445,6.865277777777778,58.28277777777778 diff --git a/doc/conf.py b/doc/conf.py index 5c4b3b89..f0d1ca37 100644 --- a/doc/conf.py +++ b/doc/conf.py @@ -53,7 +53,6 @@ extensions = [ autodoc_mock_imports = [ "atlite", "snakemake", - "pycountry", "rioxarray", "country_converter", "tabula", @@ -341,4 +340,6 @@ texinfo_documents = [ # Example configuration for intersphinx: refer to the Python standard library. -intersphinx_mapping = {"https://docs.python.org/": None} +intersphinx_mapping = { + "https://docs.python.org/": ("https://docs.python.org/3", None), +} diff --git a/doc/configtables/biomass.csv b/doc/configtables/biomass.csv index f5b4841f..865d247e 100644 --- a/doc/configtables/biomass.csv +++ b/doc/configtables/biomass.csv @@ -5,3 +5,5 @@ classes ,,, -- solid biomass,--,Array of biomass comodity,The comodity that are included as solid biomass -- not included,--,Array of biomass comodity,The comodity that are not included as a biomass potential -- biogas,--,Array of biomass comodity,The comodity that are included as biogas +share_unsustainable_use_retained,--,Dictionary with planning horizons as keys., Share of unsustainable biomass use retained using primary production of Eurostat data as reference +share_sustainable_potential_available,--,Dictionary with planning horizons as keys., Share determines phase-in of ENSPRESO biomass potentials diff --git a/doc/configtables/enable.csv b/doc/configtables/enable.csv index c74d0eff..0268319e 100644 --- a/doc/configtables/enable.csv +++ b/doc/configtables/enable.csv @@ -1,6 +1,5 @@ ,Unit,Values,Description enable,str or bool,"{auto, true, false}","Switch to include (true) or exclude (false) the retrieve_* rules of snakemake into the workflow; 'auto' sets true|false based on availability of an internet connection to prevent issues with snakemake failing due to lack of internet connection." -prepare_links_p_nom,bool,"{true, false}","Switch to retrieve current HVDC projects from `Wikipedia `_" retrieve_databundle,bool,"{true, false}","Switch to retrieve databundle from zenodo via the rule :mod:`retrieve_databundle` or whether to keep a custom databundle located in the corresponding folder." retrieve_cost_data,bool,"{true, false}","Switch to retrieve technology cost data from `technology-data repository `_." build_cutout,bool,"{true, false}","Switch to enable the building of cutouts via the rule :mod:`build_cutout`." diff --git a/doc/configtables/sector.csv b/doc/configtables/sector.csv index c96aa629..e0deb8ca 100644 --- a/doc/configtables/sector.csv +++ b/doc/configtables/sector.csv @@ -9,6 +9,18 @@ district_heating,--,,`prepare_sector_network.py `_ to one to save memory. ,,, bev_dsm_restriction _value,--,float,Adds a lower state of charge (SOC) limit for battery electric vehicles (BEV) to manage its own energy demand (DSM). Located in `build_transport_demand.py `_. Set to 0 for no restriction on BEV DSM @@ -72,7 +84,7 @@ boilers,--,"{true, false}",Add option for transforming gas into heat using gas b resistive_heaters,--,"{true, false}",Add option for transforming electricity into heat using resistive heaters (independently from gas boilers) oil_boilers,--,"{true, false}",Add option for transforming oil into heat using boilers biomass_boiler,--,"{true, false}",Add option for transforming biomass into heat using boilers -overdimension_individual_heating,--,"float",Add option for overdimensioning individual heating systems by a certain factor. This allows them to cover heat demand peaks e.g. 10% higher than those in the data with a setting of 1.1. +overdimension_individual_heating,--,float,Add option for overdimensioning individual heating systems by a certain factor. This allows them to cover heat demand peaks e.g. 10% higher than those in the data with a setting of 1.1. chp,--,"{true, false}",Add option for using Combined Heat and Power (CHP) micro_chp,--,"{true, false}",Add option for using Combined Heat and Power (CHP) for decentral areas. solar_thermal,--,"{true, false}",Add option for using solar thermal to generate heat. @@ -139,6 +151,7 @@ biogas_upgrading_cc,--,"{true, false}",Add option to capture CO2 from biomass up conventional_generation,,,Add a more detailed description of conventional carriers. Any power generation requires the consumption of fuel from nodes representing that fuel. biomass_to_liquid,--,"{true, false}",Add option for transforming solid biomass into liquid fuel with the same properties as oil biosng,--,"{true, false}",Add option for transforming solid biomass into synthesis gas with the same properties as natural gas +municipal_solid_waste,--,"{true, false}",Add option for municipal solid waste limit_max_growth,,, -- enable,--,"{true, false}",Add option to limit the maximum growth of a carrier -- factor,p.u.,float,The maximum growth factor of a carrier (e.g. 1.3 allows 30% larger than max historic growth) diff --git a/doc/preparation.rst b/doc/preparation.rst index 4585f4db..83f9781c 100644 --- a/doc/preparation.rst +++ b/doc/preparation.rst @@ -41,11 +41,6 @@ Rule ``build_cutout`` .. automodule:: build_cutout -Rule ``prepare_links_p_nom`` -=============================== - -.. automodule:: prepare_links_p_nom - .. _base: Rule ``base_network`` diff --git a/doc/release_notes.rst b/doc/release_notes.rst index d2453bb1..7404e2ef 100644 --- a/doc/release_notes.rst +++ b/doc/release_notes.rst @@ -10,9 +10,19 @@ Release Notes Upcoming Release ================ +* Added unsustainable biomass potentials for solid, gaseous, and liquid biomass. The potentials can be phased-out and/or + substituted by the phase-in of sustainable biomass types using the config parameters + ``biomass: share_unsustainable_use_retained`` and ``biomass: share_sustainable_potential_available``. + +* The rule ``prepare_links_p_nom`` was removed since it was outdated and not used. + +* Changed heat pump COP approximation for central heating to be based on `Jensen et al. (2018) `__ and a default forward temperature of 90C. This is more realistic for district heating than the previously used approximation method. + +* split solid biomass potentials into solid biomass and municipal solid waste. Add option to use municipal solid waste. This option is only activated in combination with the flag ``waste_to_energy`` + * Add option to import solid biomass -* Add option to produce electrobiofuels (flag ``electrobiofuels`) from solid biomass and hydrogen, as a combination of BtL and Fischer-Tropsch to make more use of the biogenic carbon +* Add option to produce electrobiofuels (flag ``electrobiofuels``) from solid biomass and hydrogen, as a combination of BtL and Fischer-Tropsch to make more use of the biogenic carbon * Add flag ``sector: fossil_fuels`` in config to remove the option of importing fossil fuels diff --git a/doc/requirements.txt b/doc/requirements.txt index a1cd0a5c..dca414fc 100644 --- a/doc/requirements.txt +++ b/doc/requirements.txt @@ -17,7 +17,6 @@ tabula-py # cartopy scikit-learn -pycountry pyyaml seaborn memory_profiler diff --git a/envs/environment.yaml b/envs/environment.yaml index febd6ea2..c8d8a633 100644 --- a/envs/environment.yaml +++ b/envs/environment.yaml @@ -18,7 +18,6 @@ dependencies: # Dependencies of the workflow itself - xlrd - openpyxl!=3.1.1 -- pycountry - seaborn - snakemake-minimal>=8.14 - memory_profiler diff --git a/rules/build_electricity.smk b/rules/build_electricity.smk index 18ff8230..34472f27 100644 --- a/rules/build_electricity.smk +++ b/rules/build_electricity.smk @@ -2,21 +2,6 @@ # # SPDX-License-Identifier: MIT -if config["enable"].get("prepare_links_p_nom", False): - - rule prepare_links_p_nom: - output: - "data/links_p_nom.csv", - log: - logs("prepare_links_p_nom.log"), - threads: 1 - resources: - mem_mb=1500, - conda: - "../envs/environment.yaml" - script: - "../scripts/prepare_links_p_nom.py" - rule build_electricity_demand: params: @@ -106,8 +91,8 @@ rule build_shapes: params: countries=config_provider("countries"), input: - naturalearth=ancient("data/bundle/naturalearth/ne_10m_admin_0_countries.shp"), - eez=ancient("data/bundle/eez/World_EEZ_v8_2014.shp"), + naturalearth=ancient("data/naturalearth/ne_10m_admin_0_countries_deu.shp"), + eez=ancient("data/eez/World_EEZ_v12_20231025_gpkg/eez_v12.gpkg"), nuts3=ancient("data/bundle/NUTS_2013_60M_SH/data/NUTS_RG_60M_2013.shp"), nuts3pop=ancient("data/bundle/nama_10r_3popgdp.tsv.gz"), nuts3gdp=ancient("data/bundle/nama_10r_3gdp.tsv.gz"), diff --git a/rules/build_sector.smk b/rules/build_sector.smk index d8309899..9136cf7c 100644 --- a/rules/build_sector.smk +++ b/rules/build_sector.smk @@ -217,13 +217,27 @@ rule build_temperature_profiles: rule build_cop_profiles: params: - heat_pump_sink_T=config_provider("sector", "heat_pump_sink_T"), + heat_pump_sink_T_decentral_heating=config_provider( + "sector", "heat_pump_sink_T_individual_heating" + ), + forward_temperature_central_heating=config_provider( + "sector", "district_heating", "forward_temperature" + ), + return_temperature_central_heating=config_provider( + "sector", "district_heating", "return_temperature" + ), + heat_source_cooling_central_heating=config_provider( + "sector", "district_heating", "heat_source_cooling" + ), + heat_pump_cop_approximation_central_heating=config_provider( + "sector", "district_heating", "heat_pump_cop_approximation" + ), + heat_pump_sources=config_provider("sector", "heat_pump_sources"), input: temp_soil_total=resources("temp_soil_total_elec_s{simpl}_{clusters}.nc"), temp_air_total=resources("temp_air_total_elec_s{simpl}_{clusters}.nc"), output: - cop_soil_total=resources("cop_soil_total_elec_s{simpl}_{clusters}.nc"), - cop_air_total=resources("cop_air_total_elec_s{simpl}_{clusters}.nc"), + cop_profiles=resources("cop_profiles_elec_s{simpl}_{clusters}.nc"), resources: mem_mb=20000, log: @@ -233,7 +247,7 @@ rule build_cop_profiles: conda: "../envs/environment.yaml" script: - "../scripts/build_cop_profiles.py" + "../scripts/build_cop_profiles/run.py" def solar_thermal_cutout(wildcards): @@ -331,7 +345,8 @@ rule build_biomass_potentials: "https://zenodo.org/records/10356004/files/ENSPRESO_BIOMASS.xlsx", keep_local=True, ), - nuts2="data/bundle/nuts/NUTS_RG_10M_2013_4326_LEVL_2.geojson", # https://gisco-services.ec.europa.eu/distribution/v2/nuts/download/#nuts21 + eurostat="data/eurostat/Balances-April2023", + nuts2="data/bundle/nuts/NUTS_RG_10M_2013_4326_LEVL_2.geojson", regions_onshore=resources("regions_onshore_elec_s{simpl}_{clusters}.geojson"), nuts3_population=ancient("data/bundle/nama_10r_3popgdp.tsv.gz"), swiss_cantons=ancient("data/ch_cantons.csv"), @@ -344,7 +359,7 @@ rule build_biomass_potentials: biomass_potentials=resources( "biomass_potentials_s{simpl}_{clusters}_{planning_horizons}.csv" ), - threads: 1 + threads: 8 resources: mem_mb=1000, log: @@ -940,7 +955,10 @@ rule prepare_sector_network: countries=config_provider("countries"), adjustments=config_provider("adjustments", "sector"), emissions_scope=config_provider("energy", "emissions"), + biomass=config_provider("biomass"), RDIR=RDIR, + heat_pump_sources=config_provider("sector", "heat_pump_sources"), + heat_systems=config_provider("sector", "heat_systems"), input: unpack(input_profile_offwind), **rules.cluster_gas_network.output, @@ -1020,8 +1038,7 @@ rule prepare_sector_network: ), temp_soil_total=resources("temp_soil_total_elec_s{simpl}_{clusters}.nc"), temp_air_total=resources("temp_air_total_elec_s{simpl}_{clusters}.nc"), - cop_soil_total=resources("cop_soil_total_elec_s{simpl}_{clusters}.nc"), - cop_air_total=resources("cop_air_total_elec_s{simpl}_{clusters}.nc"), + cop_profiles=resources("cop_profiles_elec_s{simpl}_{clusters}.nc"), solar_thermal_total=lambda w: ( resources("solar_thermal_total_elec_s{simpl}_{clusters}.nc") if config_provider("sector", "solar_thermal")(w) diff --git a/rules/retrieve.smk b/rules/retrieve.smk index 18b0ddd2..86c6b998 100644 --- a/rules/retrieve.smk +++ b/rules/retrieve.smk @@ -4,6 +4,7 @@ import requests from datetime import datetime, timedelta +from shutil import move, unpack_archive if config["enable"].get("retrieve", "auto") == "auto": config["enable"]["retrieve"] = has_internet_access() @@ -15,8 +16,6 @@ if config["enable"]["retrieve"] is False: if config["enable"]["retrieve"] and config["enable"].get("retrieve_databundle", True): datafiles = [ "je-e-21.03.02.xls", - "eez/World_EEZ_v8_2014.shp", - "naturalearth/ne_10m_admin_0_countries.shp", "NUTS_2013_60M_SH/data/NUTS_RG_60M_2013.shp", "nama_10r_3popgdp.tsv.gz", "nama_10r_3gdp.tsv.gz", @@ -53,6 +52,8 @@ if config["enable"]["retrieve"] and config["enable"].get("retrieve_databundle", log: "logs/retrieve_eurostat_data.log", retries: 2 + conda: + "../envs/retrieve.yaml" script: "../scripts/retrieve_eurostat_data.py" @@ -62,6 +63,8 @@ if config["enable"]["retrieve"] and config["enable"].get("retrieve_databundle", log: "logs/retrieve_eurostat_household_data.log", retries: 2 + conda: + "../envs/retrieve.yaml" script: "../scripts/retrieve_eurostat_household_data.py" @@ -211,6 +214,64 @@ if config["enable"]["retrieve"]: move(input[0], output[0]) +if config["enable"]["retrieve"]: + + rule retrieve_eez: + params: + zip="data/eez/World_EEZ_v12_20231025_gpkg.zip", + output: + gpkg="data/eez/World_EEZ_v12_20231025_gpkg/eez_v12.gpkg", + run: + import os + import requests + from uuid import uuid4 + + name = str(uuid4())[:8] + org = str(uuid4())[:8] + + response = requests.post( + "https://www.marineregions.org/download_file.php", + params={"name": "World_EEZ_v12_20231025_gpkg.zip"}, + data={ + "name": name, + "organisation": org, + "email": f"{name}@{org}.org", + "country": "Germany", + "user_category": "academia", + "purpose_category": "Research", + "agree": "1", + }, + ) + + with open(params["zip"], "wb") as f: + f.write(response.content) + output_folder = Path(params["zip"]).parent + unpack_archive(params["zip"], output_folder) + os.remove(params["zip"]) + + + +if config["enable"]["retrieve"]: + + # Download directly from naciscdn.org which is a redirect from naturalearth.com + # (https://www.naturalearthdata.com/downloads/10m-cultural-vectors/10m-admin-0-countries/) + # Use point-of-view (POV) variant of Germany so that Crimea is included. + rule retrieve_naturalearth_countries: + input: + storage( + "https://naciscdn.org/naturalearth/10m/cultural/ne_10m_admin_0_countries_deu.zip" + ), + params: + zip="data/naturalearth/ne_10m_admin_0_countries_deu.zip", + output: + countries="data/naturalearth/ne_10m_admin_0_countries_deu.shp", + run: + move(input[0], params["zip"]) + output_folder = Path(output["countries"]).parent + unpack_archive(params["zip"], output_folder) + os.remove(params["zip"]) + + if config["enable"]["retrieve"]: # Some logic to find the correct file URL # Sometimes files are released delayed or ahead of schedule, check which file is currently available diff --git a/rules/solve_myopic.smk b/rules/solve_myopic.smk index 21fb7169..8d7fa284 100644 --- a/rules/solve_myopic.smk +++ b/rules/solve_myopic.smk @@ -9,6 +9,7 @@ rule add_existing_baseyear: sector=config_provider("sector"), existing_capacities=config_provider("existing_capacities"), costs=config_provider("costs"), + heat_pump_sources=config_provider("sector", "heat_pump_sources"), input: network=RESULTS + "prenetworks/elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}.nc", @@ -21,8 +22,7 @@ rule add_existing_baseyear: config_provider("scenario", "planning_horizons", 0)(w) ) ), - cop_soil_total=resources("cop_soil_total_elec_s{simpl}_{clusters}.nc"), - cop_air_total=resources("cop_air_total_elec_s{simpl}_{clusters}.nc"), + cop_profiles=resources("cop_profiles_elec_s{simpl}_{clusters}.nc"), existing_heating_distribution=resources( "existing_heating_distribution_elec_s{simpl}_{clusters}_{planning_horizons}.csv" ), @@ -69,6 +69,7 @@ rule add_brownfield: snapshots=config_provider("snapshots"), drop_leap_day=config_provider("enable", "drop_leap_day"), carriers=config_provider("electricity", "renewable_carriers"), + heat_pump_sources=config_provider("sector", "heat_pump_sources"), input: unpack(input_profile_tech_brownfield), simplify_busmap=resources("busmap_elec_s{simpl}.csv"), @@ -77,8 +78,7 @@ rule add_brownfield: + "prenetworks/elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}.nc", network_p=solved_previous_horizon, #solved network at previous time step costs=resources("costs_{planning_horizons}.csv"), - cop_soil_total=resources("cop_soil_total_elec_s{simpl}_{clusters}.nc"), - cop_air_total=resources("cop_air_total_elec_s{simpl}_{clusters}.nc"), + cop_profiles=resources("cop_profiles_elec_s{simpl}_{clusters}.nc"), output: RESULTS + "prenetworks-brownfield/elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}.nc", diff --git a/rules/solve_perfect.smk b/rules/solve_perfect.smk index 51cb3920..a06c6dfa 100644 --- a/rules/solve_perfect.smk +++ b/rules/solve_perfect.smk @@ -7,6 +7,7 @@ rule add_existing_baseyear: sector=config_provider("sector"), existing_capacities=config_provider("existing_capacities"), costs=config_provider("costs"), + heat_pump_sources=config_provider("sector", "heat_pump_sources"), input: network=RESULTS + "prenetworks/elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}.nc", @@ -19,8 +20,7 @@ rule add_existing_baseyear: config_provider("scenario", "planning_horizons", 0)(w) ) ), - cop_soil_total=resources("cop_soil_total_elec_s{simpl}_{clusters}.nc"), - cop_air_total=resources("cop_air_total_elec_s{simpl}_{clusters}.nc"), + cop_profiles=resources("cop_profiles_elec_s{simpl}_{clusters}.nc"), existing_heating_distribution=resources( "existing_heating_distribution_elec_s{simpl}_{clusters}_{planning_horizons}.csv" ), diff --git a/scripts/add_existing_baseyear.py b/scripts/add_existing_baseyear.py index d5333675..f67c38d9 100644 --- a/scripts/add_existing_baseyear.py +++ b/scripts/add_existing_baseyear.py @@ -24,6 +24,10 @@ from _helpers import ( from add_electricity import sanitize_carriers from prepare_sector_network import cluster_heat_buses, define_spatial, prepare_costs +from scripts.definitions.heat_sector import HeatSector +from scripts.definitions.heat_system import HeatSystem +from scripts.definitions.heat_system_type import HeatSystemType + logger = logging.getLogger(__name__) cc = coco.CountryConverter() idx = pd.IndexSlice @@ -416,14 +420,14 @@ def add_power_capacities_installed_before_baseyear(n, grouping_years, costs, bas def add_heating_capacities_installed_before_baseyear( - n, - baseyear, - grouping_years, - ashp_cop, - gshp_cop, - time_dep_hp_cop, - costs, - default_lifetime, + n: pypsa.Network, + baseyear: int, + grouping_years: list, + cop: dict, + time_dep_hp_cop: bool, + costs: pd.DataFrame, + default_lifetime: int, + existing_heating: pd.DataFrame, ): """ Parameters @@ -435,141 +439,158 @@ def add_heating_capacities_installed_before_baseyear( currently assumed heating capacities split between residential and services proportional to heating load in both 50% capacities in rural buses 50% in urban buses + cop: xr.DataArray + DataArray with time-dependent coefficients of performance (COPs) heat pumps. Coordinates are heat sources (see config), heat system types (see :file:`scripts/enums/HeatSystemType.py`), nodes and snapshots. + time_dep_hp_cop: bool + If True, time-dependent (dynamic) COPs are used for heat pumps """ logger.debug(f"Adding heating capacities installed before {baseyear}") - existing_heating = pd.read_csv( - snakemake.input.existing_heating_distribution, header=[0, 1], index_col=0 - ) + for heat_system in existing_heating.columns.get_level_values(0).unique(): + heat_system = HeatSystem(heat_system) - for name in existing_heating.columns.get_level_values(0).unique(): - name_type = "central" if name == "urban central" else "decentral" + nodes = pd.Index( + n.buses.location[n.buses.index.str.contains(f"{heat_system} heat")] + ) - nodes = pd.Index(n.buses.location[n.buses.index.str.contains(f"{name} heat")]) - - if (name_type != "central") and options["electricity_distribution_grid"]: + if (not heat_system == HeatSystem.URBAN_CENTRAL) and options[ + "electricity_distribution_grid" + ]: nodes_elec = nodes + " low voltage" else: nodes_elec = nodes - heat_pump_type = "air" if "urban" in name else "ground" - - # Add heat pumps - costs_name = f"decentral {heat_pump_type}-sourced heat pump" - - cop = {"air": ashp_cop, "ground": gshp_cop} - - if time_dep_hp_cop: - efficiency = cop[heat_pump_type][nodes] - else: - efficiency = costs.at[costs_name, "efficiency"] - - too_large_grouping_years = [gy for gy in grouping_years if gy >= int(baseyear)] - if too_large_grouping_years: - logger.warning( - f"Grouping years >= baseyear are ignored. Dropping {too_large_grouping_years}." - ) - valid_grouping_years = pd.Series( - [ - int(grouping_year) - for grouping_year in grouping_years - if int(grouping_year) + default_lifetime > int(baseyear) - and int(grouping_year) < int(baseyear) + too_large_grouping_years = [ + gy for gy in grouping_years if gy >= int(baseyear) ] - ) + if too_large_grouping_years: + logger.warning( + f"Grouping years >= baseyear are ignored. Dropping {too_large_grouping_years}." + ) + valid_grouping_years = pd.Series( + [ + int(grouping_year) + for grouping_year in grouping_years + if int(grouping_year) + default_lifetime > int(baseyear) + and int(grouping_year) < int(baseyear) + ] + ) - assert valid_grouping_years.is_monotonic_increasing + assert valid_grouping_years.is_monotonic_increasing - # get number of years of each interval - _years = valid_grouping_years.diff() - # Fill NA from .diff() with value for the first interval - _years[0] = valid_grouping_years[0] - baseyear + default_lifetime - # Installation is assumed to be linear for the past - ratios = _years / _years.sum() + # get number of years of each interval + _years = valid_grouping_years.diff() + # Fill NA from .diff() with value for the first interval + _years[0] = valid_grouping_years[0] - baseyear + default_lifetime + # Installation is assumed to be linear for the past + ratios = _years / _years.sum() for ratio, grouping_year in zip(ratios, valid_grouping_years): + # Add heat pumps + for heat_source in snakemake.params.heat_pump_sources[ + heat_system.system_type.value + ]: + costs_name = heat_system.heat_pump_costs_name(heat_source) - n.madd( - "Link", - nodes, - suffix=f" {name} {heat_pump_type} heat pump-{grouping_year}", - bus0=nodes_elec, - bus1=nodes + " " + name + " heat", - carrier=f"{name} {heat_pump_type} heat pump", - efficiency=efficiency, - capital_cost=costs.at[costs_name, "efficiency"] - * costs.at[costs_name, "fixed"], - p_nom=existing_heating.loc[nodes, (name, f"{heat_pump_type} heat pump")] - * ratio - / costs.at[costs_name, "efficiency"], - build_year=int(grouping_year), - lifetime=costs.at[costs_name, "lifetime"], - ) + efficiency = ( + cop.sel( + heat_system=heat_system.system_type.value, + heat_source=heat_source, + name=nodes, + ) + .to_pandas() + .reindex(index=n.snapshots) + if time_dep_hp_cop + else costs.at[costs_name, "efficiency"] + ) + + n.madd( + "Link", + nodes, + suffix=f" {heat_system} {heat_source} heat pump-{grouping_year}", + bus0=nodes_elec, + bus1=nodes + " " + heat_system.value + " heat", + carrier=f"{heat_system} {heat_source} heat pump", + efficiency=efficiency, + capital_cost=costs.at[costs_name, "efficiency"] + * costs.at[costs_name, "fixed"], + p_nom=existing_heating.loc[ + nodes, (heat_system.value, f"{heat_source} heat pump") + ] + * ratio + / costs.at[costs_name, "efficiency"], + build_year=int(grouping_year), + lifetime=costs.at[costs_name, "lifetime"], + ) # add resistive heater, gas boilers and oil boilers n.madd( "Link", nodes, - suffix=f" {name} resistive heater-{grouping_year}", + suffix=f" {heat_system} resistive heater-{grouping_year}", bus0=nodes_elec, - bus1=nodes + " " + name + " heat", - carrier=name + " resistive heater", - efficiency=costs.at[f"{name_type} resistive heater", "efficiency"], + bus1=nodes + " " + heat_system.value + " heat", + carrier=heat_system.value + " resistive heater", + efficiency=costs.at[ + heat_system.resistive_heater_costs_name, "efficiency" + ], capital_cost=( - costs.at[f"{name_type} resistive heater", "efficiency"] - * costs.at[f"{name_type} resistive heater", "fixed"] + costs.at[heat_system.resistive_heater_costs_name, "efficiency"] + * costs.at[heat_system.resistive_heater_costs_name, "fixed"] ), p_nom=( - existing_heating.loc[nodes, (name, "resistive heater")] + existing_heating.loc[nodes, (heat_system.value, "resistive heater")] * ratio - / costs.at[f"{name_type} resistive heater", "efficiency"] + / costs.at[heat_system.resistive_heater_costs_name, "efficiency"] ), build_year=int(grouping_year), - lifetime=costs.at[f"{name_type} resistive heater", "lifetime"], + lifetime=costs.at[heat_system.resistive_heater_costs_name, "lifetime"], ) n.madd( "Link", nodes, - suffix=f" {name} gas boiler-{grouping_year}", + suffix=f"{heat_system} gas boiler-{grouping_year}", bus0="EU gas" if "EU gas" in spatial.gas.nodes else nodes + " gas", - bus1=nodes + " " + name + " heat", + bus1=nodes + " " + heat_system.value + " heat", bus2="co2 atmosphere", - carrier=name + " gas boiler", - efficiency=costs.at[f"{name_type} gas boiler", "efficiency"], + carrier=heat_system.value + " gas boiler", + efficiency=costs.at[heat_system.gas_boiler_costs_name, "efficiency"], efficiency2=costs.at["gas", "CO2 intensity"], capital_cost=( - costs.at[f"{name_type} gas boiler", "efficiency"] - * costs.at[f"{name_type} gas boiler", "fixed"] + costs.at[heat_system.gas_boiler_costs_name, "efficiency"] + * costs.at[heat_system.gas_boiler_costs_name, "fixed"] ), p_nom=( - existing_heating.loc[nodes, (name, "gas boiler")] + existing_heating.loc[nodes, (heat_system.value, "gas boiler")] * ratio - / costs.at[f"{name_type} gas boiler", "efficiency"] + / costs.at[heat_system.gas_boiler_costs_name, "efficiency"] ), build_year=int(grouping_year), - lifetime=costs.at[f"{name_type} gas boiler", "lifetime"], + lifetime=costs.at[heat_system.gas_boiler_costs_name, "lifetime"], ) n.madd( "Link", nodes, - suffix=f" {name} oil boiler-{grouping_year}", + suffix=f" {heat_system} oil boiler-{grouping_year}", bus0=spatial.oil.nodes, - bus1=nodes + " " + name + " heat", + bus1=nodes + " " + heat_system.value + " heat", bus2="co2 atmosphere", - carrier=name + " oil boiler", - efficiency=costs.at["decentral oil boiler", "efficiency"], + carrier=heat_system.value + " oil boiler", + efficiency=costs.at[heat_system.oil_boiler_costs_name, "efficiency"], efficiency2=costs.at["oil", "CO2 intensity"], - capital_cost=costs.at["decentral oil boiler", "efficiency"] - * costs.at["decentral oil boiler", "fixed"], + capital_cost=costs.at[heat_system.oil_boiler_costs_name, "efficiency"] + * costs.at[heat_system.oil_boiler_costs_name, "fixed"], p_nom=( - existing_heating.loc[nodes, (name, "oil boiler")] + existing_heating.loc[nodes, (heat_system.value, "oil boiler")] * ratio - / costs.at["decentral oil boiler", "efficiency"] + / costs.at[heat_system.oil_boiler_costs_name, "efficiency"] ), build_year=int(grouping_year), - lifetime=costs.at[f"{name_type} gas boiler", "lifetime"], + lifetime=costs.at[ + f"{heat_system.central_or_decentral} gas boiler", "lifetime" + ], ) # delete links with p_nom=nan corresponding to extra nodes in country @@ -639,29 +660,22 @@ if __name__ == "__main__": ) if options["heating"]: - time_dep_hp_cop = options["time_dep_hp_cop"] - ashp_cop = ( - xr.open_dataarray(snakemake.input.cop_air_total) - .to_pandas() - .reindex(index=n.snapshots) - ) - gshp_cop = ( - xr.open_dataarray(snakemake.input.cop_soil_total) - .to_pandas() - .reindex(index=n.snapshots) - ) - default_lifetime = snakemake.params.existing_capacities[ - "default_heating_lifetime" - ] + add_heating_capacities_installed_before_baseyear( - n, - baseyear, - grouping_years_heat, - ashp_cop, - gshp_cop, - time_dep_hp_cop, - costs, - default_lifetime, + n=n, + baseyear=baseyear, + grouping_years=grouping_years_heat, + cop=xr.open_dataarray(snakemake.input.cop_profiles), + time_dep_hp_cop=options["time_dep_hp_cop"], + costs=costs, + default_lifetime=snakemake.params.existing_capacities[ + "default_heating_lifetime" + ], + existing_heating=pd.read_csv( + snakemake.input.existing_heating_distribution, + header=[0, 1], + index_col=0, + ), ) if options.get("cluster_heat_buses", False): diff --git a/scripts/build_biomass_potentials.py b/scripts/build_biomass_potentials.py index 883734eb..0a2692e8 100644 --- a/scripts/build_biomass_potentials.py +++ b/scripts/build_biomass_potentials.py @@ -13,11 +13,51 @@ import geopandas as gpd import numpy as np import pandas as pd from _helpers import configure_logging, set_scenario_config +from build_energy_totals import build_eurostat logger = logging.getLogger(__name__) AVAILABLE_BIOMASS_YEARS = [2010, 2020, 2030, 2040, 2050] +def _calc_unsustainable_potential(df, df_unsustainable, share_unsus, resource_type): + """ + Calculate the unsustainable biomass potential for a given resource type or + regex. + + Parameters + ---------- + df : pd.DataFrame + The dataframe with sustainable biomass potentials. + df_unsustainable : pd.DataFrame + The dataframe with unsustainable biomass potentials. + share_unsus : float + The share of unsustainable biomass potential retained. + resource_type : str or regex + The resource type to calculate the unsustainable potential for. + + Returns + ------- + pd.Series + The unsustainable biomass potential for the given resource type or regex. + """ + + if "|" in resource_type: + resource_potential = df_unsustainable.filter(regex=resource_type).sum(axis=1) + else: + resource_potential = df_unsustainable[resource_type] + + return ( + df.apply( + lambda c: c.sum() + / df.loc[df.index.str[:2] == c.name[:2]].sum().sum() + * resource_potential.loc[c.name[:2]], + axis=1, + ) + .mul(share_unsus) + .clip(lower=0) + ) + + def build_nuts_population_data(year=2013): pop = pd.read_csv( snakemake.input.nuts3_population, @@ -211,15 +251,104 @@ def convert_nuts2_to_regions(bio_nuts2, regions): return bio_regions +def add_unsustainable_potentials(df): + """ + Add unsustainable biomass potentials to the given dataframe. The difference + between the data of JRC and Eurostat is assumed to be unsustainable + biomass. + + Parameters + ---------- + df : pd.DataFrame + The dataframe with sustainable biomass potentials. + unsustainable_biomass : str + Path to the file with unsustainable biomass potentials. + + Returns + ------- + pd.DataFrame + The dataframe with added unsustainable biomass potentials. + """ + if "GB" in snakemake.config["countries"]: + latest_year = 2019 + else: + latest_year = 2021 + idees_rename = {"GR": "EL", "GB": "UK"} + df_unsustainable = ( + build_eurostat( + countries=snakemake.config["countries"], + input_eurostat=snakemake.input.eurostat, + nprocesses=int(snakemake.threads), + ) + .xs( + max(min(latest_year, int(snakemake.wildcards.planning_horizons)), 1990), + level=1, + ) + .xs("Primary production", level=2) + .droplevel([1, 2, 3]) + ) + + df_unsustainable.index = df_unsustainable.index.str.strip() + df_unsustainable = df_unsustainable.rename( + {v: k for k, v in idees_rename.items()}, axis=0 + ) + + bio_carriers = [ + "Primary solid biofuels", + "Biogases", + "Renewable municipal waste", + "Pure biogasoline", + "Blended biogasoline", + "Pure biodiesels", + "Blended biodiesels", + "Pure bio jet kerosene", + "Blended bio jet kerosene", + "Other liquid biofuels", + ] + + df_unsustainable = df_unsustainable[bio_carriers] + + # Phase out unsustainable biomass potentials linearly from 2020 to 2035 while phasing in sustainable potentials + share_unsus = params.get("share_unsustainable_use_retained").get(investment_year) + + df_wo_ch = df.drop(df.filter(regex="CH\d", axis=0).index) + + # Calculate unsustainable solid biomass + df_wo_ch["unsustainable solid biomass"] = _calc_unsustainable_potential( + df_wo_ch, df_unsustainable, share_unsus, "Primary solid biofuels" + ) + + # Calculate unsustainable biogas + df_wo_ch["unsustainable biogas"] = _calc_unsustainable_potential( + df_wo_ch, df_unsustainable, share_unsus, "Biogases" + ) + + # Calculate unsustainable bioliquids + df_wo_ch["unsustainable bioliquids"] = _calc_unsustainable_potential( + df_wo_ch, + df_unsustainable, + share_unsus, + resource_type="gasoline|diesel|kerosene|liquid", + ) + + share_sus = params.get("share_sustainable_potential_available").get(investment_year) + df *= share_sus + + df = df.join(df_wo_ch.filter(like="unsustainable")).fillna(0) + + return df + + if __name__ == "__main__": if "snakemake" not in globals(): + from _helpers import mock_snakemake snakemake = mock_snakemake( "build_biomass_potentials", simpl="", - clusters="5", - planning_horizons=2050, + clusters="37", + planning_horizons=2020, ) configure_logging(snakemake) @@ -269,6 +398,8 @@ if __name__ == "__main__": grouper = {v: k for k, vv in params["classes"].items() for v in vv} df = df.T.groupby(grouper).sum().T + df = add_unsustainable_potentials(df) + df *= 1e6 # TWh/a to MWh/a df.index.name = "MWh/a" diff --git a/scripts/build_cop_profiles.py b/scripts/build_cop_profiles.py deleted file mode 100644 index a6d99947..00000000 --- a/scripts/build_cop_profiles.py +++ /dev/null @@ -1,69 +0,0 @@ -# -*- coding: utf-8 -*- -# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors -# -# SPDX-License-Identifier: MIT -""" -Build coefficient of performance (COP) time series for air- or ground-sourced -heat pumps. - -The COP is approximated as a quatratic function of the temperature difference between source and -sink, based on Staffell et al. 2012. - -This rule is executed in ``build_sector.smk``. - -Relevant Settings ------------------ - -.. code:: yaml - heat_pump_sink_T: - - -Inputs: -------- -- ``resources//temp_soil_total_elec_s_.nc``: Soil temperature (total) time series. -- ``resources//temp_air_total_elec_s_.nc``: Ambient air temperature (total) time series. - -Outputs: --------- -- ``resources/cop_soil_total_elec_s_.nc``: COP (ground-sourced) time series (total). -- ``resources/cop_air_total_elec_s_.nc``: COP (air-sourced) time series (total). - - -References ----------- -[1] Staffell et al., Energy & Environmental Science 11 (2012): A review of domestic heat pumps, https://doi.org/10.1039/C2EE22653G. -""" - -import xarray as xr -from _helpers import set_scenario_config - - -def coefficient_of_performance(delta_T, source="air"): - if source == "air": - return 6.81 - 0.121 * delta_T + 0.000630 * delta_T**2 - elif source == "soil": - return 8.77 - 0.150 * delta_T + 0.000734 * delta_T**2 - else: - raise NotImplementedError("'source' must be one of ['air', 'soil']") - - -if __name__ == "__main__": - if "snakemake" not in globals(): - from _helpers import mock_snakemake - - snakemake = mock_snakemake( - "build_cop_profiles", - simpl="", - clusters=48, - ) - - set_scenario_config(snakemake) - - for source in ["air", "soil"]: - source_T = xr.open_dataarray(snakemake.input[f"temp_{source}_total"]) - - delta_T = snakemake.params.heat_pump_sink_T - source_T - - cop = coefficient_of_performance(delta_T, source) - - cop.to_netcdf(snakemake.output[f"cop_{source}_total"]) diff --git a/scripts/build_cop_profiles/BaseCopApproximator.py b/scripts/build_cop_profiles/BaseCopApproximator.py new file mode 100644 index 00000000..10019b3d --- /dev/null +++ b/scripts/build_cop_profiles/BaseCopApproximator.py @@ -0,0 +1,111 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + +from abc import ABC, abstractmethod +from typing import Union + +import numpy as np +import xarray as xr + + +class BaseCopApproximator(ABC): + """ + Abstract class for approximating the coefficient of performance (COP) of a + heat pump. + + Attributes: + ---------- + forward_temperature_celsius : Union[xr.DataArray, np.array] + The forward temperature in Celsius. + source_inlet_temperature_celsius : Union[xr.DataArray, np.array] + The source inlet temperature in Celsius. + + Methods: + ------- + __init__(self, forward_temperature_celsius, source_inlet_temperature_celsius) + Initialize CopApproximator. + approximate_cop(self) + Approximate heat pump coefficient of performance (COP). + celsius_to_kelvin(t_celsius) + Convert temperature from Celsius to Kelvin. + logarithmic_mean(t_hot, t_cold) + Calculate the logarithmic mean temperature difference. + """ + + def __init__( + self, + forward_temperature_celsius: Union[xr.DataArray, np.array], + source_inlet_temperature_celsius: Union[xr.DataArray, np.array], + ): + """ + Initialize CopApproximator. + + Parameters: + ---------- + forward_temperature_celsius : Union[xr.DataArray, np.array] + The forward temperature in Celsius. + source_inlet_temperature_celsius : Union[xr.DataArray, np.array] + The source inlet temperature in Celsius. + """ + pass + + @abstractmethod + def approximate_cop(self) -> Union[xr.DataArray, np.array]: + """ + Approximate heat pump coefficient of performance (COP). + + Returns: + ------- + Union[xr.DataArray, np.array] + The calculated COP values. + """ + pass + + @staticmethod + def celsius_to_kelvin( + t_celsius: Union[float, xr.DataArray, np.array] + ) -> Union[float, xr.DataArray, np.array]: + """ + Convert temperature from Celsius to Kelvin. + + Parameters: + ---------- + t_celsius : Union[float, xr.DataArray, np.array] + Temperature in Celsius. + + Returns: + ------- + Union[float, xr.DataArray, np.array] + Temperature in Kelvin. + """ + if (np.asarray(t_celsius) > 200).any(): + raise ValueError( + "t_celsius > 200. Are you sure you are using the right units?" + ) + return t_celsius + 273.15 + + @staticmethod + def logarithmic_mean( + t_hot: Union[float, xr.DataArray, np.ndarray], + t_cold: Union[float, xr.DataArray, np.ndarray], + ) -> Union[float, xr.DataArray, np.ndarray]: + """ + Calculate the logarithmic mean temperature difference. + + Parameters: + ---------- + t_hot : Union[float, xr.DataArray, np.ndarray] + Hot temperature. + t_cold : Union[float, xr.DataArray, np.ndarray] + Cold temperature. + + Returns: + ------- + Union[float, xr.DataArray, np.ndarray] + Logarithmic mean temperature difference. + """ + if (np.asarray(t_hot <= t_cold)).any(): + raise ValueError("t_hot must be greater than t_cold") + return (t_hot - t_cold) / np.log(t_hot / t_cold) diff --git a/scripts/build_cop_profiles/CentralHeatingCopApproximator.py b/scripts/build_cop_profiles/CentralHeatingCopApproximator.py new file mode 100644 index 00000000..08dd6a1a --- /dev/null +++ b/scripts/build_cop_profiles/CentralHeatingCopApproximator.py @@ -0,0 +1,392 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + + +from typing import Union + +import numpy as np +import xarray as xr +from BaseCopApproximator import BaseCopApproximator + + +class CentralHeatingCopApproximator(BaseCopApproximator): + """ + Approximate the coefficient of performance (COP) for a heat pump in a + central heating system (district heating). + + Uses an approximation method proposed by Jensen et al. (2018) and + default parameters from Pieper et al. (2020). The method is based on + a thermodynamic heat pump model with some hard-to-know parameters + being approximated. + + Attributes: + ---------- + forward_temperature_celsius : Union[xr.DataArray, np.array] + The forward temperature in Celsius. + return_temperature_celsius : Union[xr.DataArray, np.array] + The return temperature in Celsius. + source_inlet_temperature_celsius : Union[xr.DataArray, np.array] + The source inlet temperature in Celsius. + source_outlet_temperature_celsius : Union[xr.DataArray, np.array] + The source outlet temperature in Celsius. + delta_t_pinch_point : float, optional + The pinch point temperature difference, by default 5. + isentropic_compressor_efficiency : float, optional + The isentropic compressor efficiency, by default 0.8. + heat_loss : float, optional + The heat loss, by default 0.0. + + Methods: + ------- + __init__( + forward_temperature_celsius: Union[xr.DataArray, np.array], + source_inlet_temperature_celsius: Union[xr.DataArray, np.array], + return_temperature_celsius: Union[xr.DataArray, np.array], + source_outlet_temperature_celsius: Union[xr.DataArray, np.array], + delta_t_pinch_point: float = 5, + isentropic_compressor_efficiency: float = 0.8, + heat_loss: float = 0.0, + ) -> None: + Initializes the CentralHeatingCopApproximator object. + + approximate_cop(self) -> Union[xr.DataArray, np.array]: + Calculate the coefficient of performance (COP) for the system. + + _approximate_delta_t_refrigerant_source( + self, delta_t_source: Union[xr.DataArray, np.array] + ) -> Union[xr.DataArray, np.array]: + Approximates the temperature difference between the refrigerant and the source. + + _approximate_delta_t_refrigerant_sink( + self, + refrigerant: str = "ammonia", + a: float = {"ammonia": 0.2, "isobutane": -0.0011}, + b: float = {"ammonia": 0.2, "isobutane": 0.3}, + c: float = {"ammonia": 0.016, "isobutane": 2.4}, + ) -> Union[xr.DataArray, np.array]: + Approximates the temperature difference between the refrigerant and heat sink. + + _ratio_evaporation_compression_work_approximation( + self, + refrigerant: str = "ammonia", + a: float = {"ammonia": 0.0014, "isobutane": 0.0035}, + ) -> Union[xr.DataArray, np.array]: + Calculate the ratio of evaporation to compression work based on approximation. + + _approximate_delta_t_refrigerant_sink( + self, + refrigerant: str = "ammonia", + a: float = {"ammonia": 0.2, "isobutane": -0.0011}, + b: float = {"ammonia": 0.2, "isobutane": 0.3}, + c: float = {"ammonia": 0.016, "isobutane": 2.4}, + ) -> Union[xr.DataArray, np.array]: + Approximates the temperature difference between the refrigerant and heat sink. + + _ratio_evaporation_compression_work_approximation( + self, + refrigerant: str = "ammonia", + a: float = {"ammonia": 0.0014, "isobutane": 0.0035}, + ) -> Union[xr.DataArray, np.array]: + Calculate the ratio of evaporation to compression work based on approximation. + """ + + def __init__( + self, + forward_temperature_celsius: Union[xr.DataArray, np.array], + source_inlet_temperature_celsius: Union[xr.DataArray, np.array], + return_temperature_celsius: Union[xr.DataArray, np.array], + source_outlet_temperature_celsius: Union[xr.DataArray, np.array], + delta_t_pinch_point: float = 5, + isentropic_compressor_efficiency: float = 0.8, + heat_loss: float = 0.0, + ) -> None: + """ + Initializes the CentralHeatingCopApproximator object. + + Parameters: + ---------- + forward_temperature_celsius : Union[xr.DataArray, np.array] + The forward temperature in Celsius. + return_temperature_celsius : Union[xr.DataArray, np.array] + The return temperature in Celsius. + source_inlet_temperature_celsius : Union[xr.DataArray, np.array] + The source inlet temperature in Celsius. + source_outlet_temperature_celsius : Union[xr.DataArray, np.array] + The source outlet temperature in Celsius. + delta_t_pinch_point : float, optional + The pinch point temperature difference, by default 5. + isentropic_compressor_efficiency : float, optional + The isentropic compressor efficiency, by default 0.8. + heat_loss : float, optional + The heat loss, by default 0.0. + """ + self.t_source_in_kelvin = BaseCopApproximator.celsius_to_kelvin( + source_inlet_temperature_celsius + ) + self.t_sink_out_kelvin = BaseCopApproximator.celsius_to_kelvin( + forward_temperature_celsius + ) + + self.t_sink_in_kelvin = BaseCopApproximator.celsius_to_kelvin( + return_temperature_celsius + ) + self.t_source_out = BaseCopApproximator.celsius_to_kelvin( + source_outlet_temperature_celsius + ) + + self.isentropic_efficiency_compressor_kelvin = isentropic_compressor_efficiency + self.heat_loss = heat_loss + self.delta_t_pinch = delta_t_pinch_point + + def approximate_cop(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the coefficient of performance (COP) for the system. + + Returns: + -------- + Union[xr.DataArray, np.array]: The calculated COP values. + """ + return ( + self.ideal_lorenz_cop + * ( + ( + 1 + + (self.delta_t_refrigerant_sink + self.delta_t_pinch) + / self.t_sink_mean_kelvin + ) + / ( + 1 + + ( + self.delta_t_refrigerant_sink + + self.delta_t_refrigerant_source + + 2 * self.delta_t_pinch + ) + / self.delta_t_lift + ) + ) + * self.isentropic_efficiency_compressor_kelvin + * (1 - self.ratio_evaporation_compression_work) + + 1 + - self.isentropic_efficiency_compressor_kelvin + - self.heat_loss + ) + + @property + def t_sink_mean_kelvin(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the logarithmic mean temperature difference between the cold + and hot sinks. + + Returns + ------- + Union[xr.DataArray, np.array] + The mean temperature difference. + """ + return BaseCopApproximator.logarithmic_mean( + t_cold=self.t_sink_in_kelvin, t_hot=self.t_sink_out_kelvin + ) + + @property + def t_source_mean_kelvin(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the logarithmic mean temperature of the heat source. + + Returns + ------- + Union[xr.DataArray, np.array] + The mean temperature of the heat source. + """ + return BaseCopApproximator.logarithmic_mean( + t_hot=self.t_source_in_kelvin, t_cold=self.t_source_out + ) + + @property + def delta_t_lift(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the temperature lift as the difference between the + logarithmic sink and source temperatures. + + Returns + ------- + Union[xr.DataArray, np.array] + The temperature difference between the sink and source. + """ + return self.t_sink_mean_kelvin - self.t_source_mean_kelvin + + @property + def ideal_lorenz_cop(self) -> Union[xr.DataArray, np.array]: + """ + Ideal Lorenz coefficient of performance (COP). + + The ideal Lorenz COP is calculated as the ratio of the mean sink temperature + to the lift temperature difference. + + Returns + ------- + np.array + The ideal Lorenz COP. + """ + return self.t_sink_mean_kelvin / self.delta_t_lift + + @property + def delta_t_refrigerant_source(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the temperature difference between the refrigerant source + inlet and outlet. + + Returns + ------- + Union[xr.DataArray, np.array] + The temperature difference between the refrigerant source inlet and outlet. + """ + return self._approximate_delta_t_refrigerant_source( + delta_t_source=self.t_source_in_kelvin - self.t_source_out + ) + + @property + def delta_t_refrigerant_sink(self) -> Union[xr.DataArray, np.array]: + """ + Temperature difference between the refrigerant and the sink based on + approximation. + + Returns + ------- + Union[xr.DataArray, np.array] + The temperature difference between the refrigerant and the sink. + """ + return self._approximate_delta_t_refrigerant_sink() + + @property + def ratio_evaporation_compression_work(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the ratio of evaporation to compression work based on + approximation. + + Returns + ------- + Union[xr.DataArray, np.array] + The calculated ratio of evaporation to compression work. + """ + return self._ratio_evaporation_compression_work_approximation() + + @property + def delta_t_sink(self) -> Union[xr.DataArray, np.array]: + """ + Calculate the temperature difference at the sink. + + Returns + ------- + Union[xr.DataArray, np.array] + The temperature difference at the sink. + """ + return self.t_sink_out_kelvin - self.t_sink_in_kelvin + + def _approximate_delta_t_refrigerant_source( + self, delta_t_source: Union[xr.DataArray, np.array] + ) -> Union[xr.DataArray, np.array]: + """ + Approximates the temperature difference between the refrigerant and the + source. + + Parameters + ---------- + delta_t_source : Union[xr.DataArray, np.array] + The temperature difference for the refrigerant source. + + Returns + ------- + Union[xr.DataArray, np.array] + The approximate temperature difference between the refrigerant and heat source. + """ + return delta_t_source / 2 + + def _approximate_delta_t_refrigerant_sink( + self, + refrigerant: str = "ammonia", + a: float = {"ammonia": 0.2, "isobutane": -0.0011}, + b: float = {"ammonia": 0.2, "isobutane": 0.3}, + c: float = {"ammonia": 0.016, "isobutane": 2.4}, + ) -> Union[xr.DataArray, np.array]: + """ + Approximates the temperature difference between the refrigerant and + heat sink. + + Parameters: + ---------- + refrigerant : str, optional + The refrigerant used in the system. Either 'isobutane' or 'ammonia. Default is 'ammonia'. + a : float, optional + Coefficient for the temperature difference between the sink and source, default is 0.2. + b : float, optional + Coefficient for the temperature difference at the sink, default is 0.2. + c : float, optional + Constant term, default is 0.016. + + Returns: + ------- + Union[xr.DataArray, np.array] + The approximate temperature difference between the refrigerant and heat sink. + + Notes: + ------ + This function assumes ammonia as the refrigerant. + + The approximate temperature difference at the refrigerant sink is calculated using the following formula: + a * (t_sink_out - t_source_out + 2 * delta_t_pinch) + b * delta_t_sink + c + """ + if refrigerant not in a.keys(): + raise ValueError( + f"Invalid refrigerant '{refrigerant}'. Must be one of {a.keys()}" + ) + return ( + a[refrigerant] + * (self.t_sink_out_kelvin - self.t_source_out + 2 * self.delta_t_pinch) + + b[refrigerant] * self.delta_t_sink + + c[refrigerant] + ) + + def _ratio_evaporation_compression_work_approximation( + self, + refrigerant: str = "ammonia", + a: float = {"ammonia": 0.0014, "isobutane": 0.0035}, + b: float = {"ammonia": -0.0015, "isobutane": -0.0033}, + c: float = {"ammonia": 0.039, "isobutane": 0.053}, + ) -> Union[xr.DataArray, np.array]: + """ + Calculate the ratio of evaporation to compression work approximation. + + Parameters: + ---------- + refrigerant : str, optional + The refrigerant used in the system. Either 'isobutane' or 'ammonia. Default is 'ammonia'. + a : float, optional + Coefficient 'a' in the approximation equation. Default is 0.0014. + b : float, optional + Coefficient 'b' in the approximation equation. Default is -0.0015. + c : float, optional + Coefficient 'c' in the approximation equation. Default is 0.039. + + Returns: + ------- + Union[xr.DataArray, np.array] + The approximated ratio of evaporation to compression work. + + Notes: + ------ + This function assumes ammonia as the refrigerant. + + The approximation equation used is: + ratio = a * (t_sink_out - t_source_out + 2 * delta_t_pinch) + b * delta_t_sink + c + """ + if refrigerant not in a.keys(): + raise ValueError( + f"Invalid refrigerant '{refrigerant}'. Must be one of {a.keys()}" + ) + return ( + a[refrigerant] + * (self.t_sink_out_kelvin - self.t_source_out + 2 * self.delta_t_pinch) + + b[refrigerant] * self.delta_t_sink + + c[refrigerant] + ) diff --git a/scripts/build_cop_profiles/DecentralHeatingCopApproximator.py b/scripts/build_cop_profiles/DecentralHeatingCopApproximator.py new file mode 100644 index 00000000..e5f58af4 --- /dev/null +++ b/scripts/build_cop_profiles/DecentralHeatingCopApproximator.py @@ -0,0 +1,110 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + + +from typing import Union + +import numpy as np +import xarray as xr +from BaseCopApproximator import BaseCopApproximator + + +class DecentralHeatingCopApproximator(BaseCopApproximator): + """ + Approximate the coefficient of performance (COP) for a heat pump in a + decentral heating system (individual/household heating). + + Uses a quadratic regression on the temperature difference between the source and sink based on empirical data proposed by Staffell et al. 2012. + + Attributes + ---------- + forward_temperature_celsius : Union[xr.DataArray, np.array] + The forward temperature in Celsius. + source_inlet_temperature_celsius : Union[xr.DataArray, np.array] + The source inlet temperature in Celsius. + source_type : str + The source of the heat pump. Must be either 'air' or 'ground'. + + Methods + ------- + __init__(forward_temperature_celsius, source_inlet_temperature_celsius, source_type) + Initialize the DecentralHeatingCopApproximator object. + approximate_cop() + Compute the COP values using quadratic regression for air-/ground-source heat pumps. + _approximate_cop_air_source() + Evaluate quadratic regression for an air-sourced heat pump. + _approximate_cop_ground_source() + Evaluate quadratic regression for a ground-sourced heat pump. + + References + ---------- + [1] Staffell et al., Energy & Environmental Science 11 (2012): A review of domestic heat pumps, https://doi.org/10.1039/C2EE22653G. + """ + + def __init__( + self, + forward_temperature_celsius: Union[xr.DataArray, np.array], + source_inlet_temperature_celsius: Union[xr.DataArray, np.array], + source_type: str, + ): + """ + Initialize the DecentralHeatingCopApproximator object. + + Parameters + ---------- + forward_temperature_celsius : Union[xr.DataArray, np.array] + The forward temperature in Celsius. + source_inlet_temperature_celsius : Union[xr.DataArray, np.array] + The source inlet temperature in Celsius. + source_type : str + The source of the heat pump. Must be either 'air' or 'ground'. + """ + + self.delta_t = forward_temperature_celsius - source_inlet_temperature_celsius + if source_type not in ["air", "ground"]: + raise ValueError("'source_type' must be one of ['air', 'ground']") + else: + self.source_type = source_type + + def approximate_cop(self) -> Union[xr.DataArray, np.array]: + """ + Compute the COP values using quadratic regression for air-/ground- + source heat pumps. + + Returns + ------- + Union[xr.DataArray, np.array] + The calculated COP values. + """ + if self.source_type == "air": + return self._approximate_cop_air_source() + elif self.source_type == "ground": + return self._approximate_cop_ground_source() + + def _approximate_cop_air_source(self) -> Union[xr.DataArray, np.array]: + """ + Evaluate quadratic regression for an air-sourced heat pump. + + COP = 6.81 - 0.121 * delta_T + 0.000630 * delta_T^2 + + Returns + ------- + Union[xr.DataArray, np.array] + The calculated COP values. + """ + return 6.81 - 0.121 * self.delta_t + 0.000630 * self.delta_t**2 + + def _approximate_cop_ground_source(self) -> Union[xr.DataArray, np.array]: + """ + Evaluate quadratic regression for a ground-sourced heat pump. + + COP = 8.77 - 0.150 * delta_T + 0.000734 * delta_T^2 + + Returns + ------- + Union[xr.DataArray, np.array] + The calculated COP values. + """ + return 8.77 - 0.150 * self.delta_t + 0.000734 * self.delta_t**2 diff --git a/scripts/build_cop_profiles/run.py b/scripts/build_cop_profiles/run.py new file mode 100644 index 00000000..4d57db31 --- /dev/null +++ b/scripts/build_cop_profiles/run.py @@ -0,0 +1,95 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + +import sys + +import numpy as np +import pandas as pd +import xarray as xr +from _helpers import set_scenario_config +from CentralHeatingCopApproximator import CentralHeatingCopApproximator +from DecentralHeatingCopApproximator import DecentralHeatingCopApproximator + +from scripts.definitions.heat_system_type import HeatSystemType + +sys.path.append("..") + + +def get_cop( + heat_system_type: str, + heat_source: str, + source_inlet_temperature_celsius: xr.DataArray, +) -> xr.DataArray: + """ + Calculate the coefficient of performance (COP) for a heating system. + + Parameters + ---------- + heat_system_type : str + The type of heating system. + heat_source : str + The heat source used in the heating system. + source_inlet_temperature_celsius : xr.DataArray + The inlet temperature of the heat source in Celsius. + + Returns + ------- + xr.DataArray + The calculated coefficient of performance (COP) for the heating system. + """ + if HeatSystemType(heat_system_type).is_central: + return CentralHeatingCopApproximator( + forward_temperature_celsius=snakemake.params.forward_temperature_central_heating, + return_temperature_celsius=snakemake.params.return_temperature_central_heating, + source_inlet_temperature_celsius=source_inlet_temperature_celsius, + source_outlet_temperature_celsius=source_inlet_temperature_celsius + - snakemake.params.heat_source_cooling_central_heating, + ).approximate_cop() + + else: + return DecentralHeatingCopApproximator( + forward_temperature_celsius=snakemake.params.heat_pump_sink_T_decentral_heating, + source_inlet_temperature_celsius=source_inlet_temperature_celsius, + source_type=heat_source, + ).approximate_cop() + + +if __name__ == "__main__": + if "snakemake" not in globals(): + from _helpers import mock_snakemake + + snakemake = mock_snakemake( + "build_cop_profiles", + simpl="", + clusters=48, + ) + + set_scenario_config(snakemake) + + cop_all_system_types = [] + for heat_system_type, heat_sources in snakemake.params.heat_pump_sources.items(): + cop_this_system_type = [] + for heat_source in heat_sources: + source_inlet_temperature_celsius = xr.open_dataarray( + snakemake.input[f"temp_{heat_source.replace('ground', 'soil')}_total"] + ) + cop_da = get_cop( + heat_system_type=heat_system_type, + heat_source=heat_source, + source_inlet_temperature_celsius=source_inlet_temperature_celsius, + ) + cop_this_system_type.append(cop_da) + cop_all_system_types.append( + xr.concat( + cop_this_system_type, dim=pd.Index(heat_sources, name="heat_source") + ) + ) + + cop_dataarray = xr.concat( + cop_all_system_types, + dim=pd.Index(snakemake.params.heat_pump_sources.keys(), name="heat_system"), + ) + + cop_dataarray.to_netcdf(snakemake.output.cop_profiles) diff --git a/scripts/build_shapes.py b/scripts/build_shapes.py index 93a73858..411d56a4 100644 --- a/scripts/build_shapes.py +++ b/scripts/build_shapes.py @@ -26,7 +26,7 @@ Inputs .. image:: img/countries.png :scale: 33 % -- ``data/bundle/eez/World_EEZ_v8_2014.shp``: World `exclusive economic zones `_ (EEZ) +- ``data/eez/World_EEZ_v12_20231025_gpkg/eez_v12.gpkg ``: World `exclusive economic zones `_ (EEZ) .. image:: img/eez.png :scale: 33 % @@ -73,22 +73,16 @@ from functools import reduce from itertools import takewhile from operator import attrgetter +import country_converter as coco import geopandas as gpd import numpy as np import pandas as pd -import pycountry as pyc from _helpers import configure_logging, set_scenario_config from shapely.geometry import MultiPolygon, Polygon logger = logging.getLogger(__name__) - -def _get_country(target, **keys): - assert len(keys) == 1 - try: - return getattr(pyc.countries.get(**keys), target) - except (KeyError, AttributeError): - return np.nan +cc = coco.CountryConverter() def _simplify_polys(polys, minarea=0.1, tolerance=None, filterremote=True): @@ -135,22 +129,15 @@ def countries(naturalearth, country_list): return s -def eez(country_shapes, eez, country_list): +def eez(eez, country_list): df = gpd.read_file(eez) - df = df.loc[ - df["ISO_3digit"].isin( - [_get_country("alpha_3", alpha_2=c) for c in country_list] - ) - ] - df["name"] = df["ISO_3digit"].map(lambda c: _get_country("alpha_2", alpha_3=c)) + iso3_list = cc.convert(country_list, src="ISO2", to="ISO3") + df = df.query("ISO_TER1 in @iso3_list and POL_TYPE == '200NM'").copy() + df["name"] = cc.convert(df.ISO_TER1, src="ISO3", to="ISO2") s = df.set_index("name").geometry.map( lambda s: _simplify_polys(s, filterremote=False) ) - s = gpd.GeoSeries( - {k: v for k, v in s.items() if v.distance(country_shapes[k]) < 1e-3}, - crs=df.crs, - ) - s = s.to_frame("geometry") + s = s.to_frame("geometry").set_crs(df.crs) s.index.name = "name" return s @@ -262,9 +249,7 @@ if __name__ == "__main__": country_shapes = countries(snakemake.input.naturalearth, snakemake.params.countries) country_shapes.reset_index().to_file(snakemake.output.country_shapes) - offshore_shapes = eez( - country_shapes, snakemake.input.eez, snakemake.params.countries - ) + offshore_shapes = eez(snakemake.input.eez, snakemake.params.countries) offshore_shapes.reset_index().to_file(snakemake.output.offshore_shapes) europe_shape = gpd.GeoDataFrame( diff --git a/scripts/definitions/heat_sector.py b/scripts/definitions/heat_sector.py new file mode 100644 index 00000000..03bcaffd --- /dev/null +++ b/scripts/definitions/heat_sector.py @@ -0,0 +1,28 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + +from enum import Enum + + +class HeatSector(Enum): + """ + Enumeration class representing different heat sectors. + + Attributes: + RESIDENTIAL (str): Represents the residential heat sector. + SERVICES (str): Represents the services heat sector. + """ + + RESIDENTIAL = "residential" + SERVICES = "services" + + def __str__(self) -> str: + """ + Returns the string representation of the heat sector. + + Returns: + str: The string representation of the heat sector. + """ + return self.value diff --git a/scripts/definitions/heat_system.py b/scripts/definitions/heat_system.py new file mode 100644 index 00000000..b907b0fe --- /dev/null +++ b/scripts/definitions/heat_system.py @@ -0,0 +1,267 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + +from enum import Enum + +from scripts.definitions.heat_sector import HeatSector +from scripts.definitions.heat_system_type import HeatSystemType + + +class HeatSystem(Enum): + """ + Enumeration representing different heat systems. + + Attributes + ---------- + RESIDENTIAL_RURAL : str + Heat system for residential areas in rural locations. + SERVICES_RURAL : str + Heat system for service areas in rural locations. + RESIDENTIAL_URBAN_DECENTRAL : str + Heat system for residential areas in urban decentralized locations. + SERVICES_URBAN_DECENTRAL : str + Heat system for service areas in urban decentralized locations. + URBAN_CENTRAL : str + Heat system for urban central areas. + + Methods + ------- + __str__() + Returns the string representation of the heat system. + central_or_decentral() + Returns whether the heat system is central or decentralized. + system_type() + Returns the type of the heat system. + sector() + Returns the sector of the heat system. + rural() + Returns whether the heat system is for rural areas. + urban_decentral() + Returns whether the heat system is for urban decentralized areas. + urban() + Returns whether the heat system is for urban areas. + heat_demand_weighting(urban_fraction=None, dist_fraction=None) + Calculates the heat demand weighting based on urban fraction and distribution fraction. + heat_pump_costs_name(heat_source) + Generates the name for the heat pump costs based on the heat source. + """ + + RESIDENTIAL_RURAL = "residential rural" + SERVICES_RURAL = "services rural" + RESIDENTIAL_URBAN_DECENTRAL = "residential urban decentral" + SERVICES_URBAN_DECENTRAL = "services urban decentral" + URBAN_CENTRAL = "urban central" + + def __init__(self, *args): + super().__init__(*args) + + def __str__(self) -> str: + """ + Returns the string representation of the heat system. + + Returns + ------- + str + The string representation of the heat system. + """ + return self.value + + @property + def central_or_decentral(self) -> str: + """ + Returns whether the heat system is central or decentralized. + + Returns + ------- + str + "central" if the heat system is central, "decentral" otherwise. + """ + if self == HeatSystem.URBAN_CENTRAL: + return "central" + else: + return "decentral" + + @property + def system_type(self) -> HeatSystemType: + """ + Returns the type of the heat system. + + Returns + ------- + str + The type of the heat system. + + Raises + ------ + RuntimeError + If the heat system is invalid. + """ + if self == HeatSystem.URBAN_CENTRAL: + return HeatSystemType.URBAN_CENTRAL + elif ( + self == HeatSystem.RESIDENTIAL_URBAN_DECENTRAL + or self == HeatSystem.SERVICES_URBAN_DECENTRAL + ): + return HeatSystemType.URBAN_DECENTRAL + elif self == HeatSystem.RESIDENTIAL_RURAL or self == HeatSystem.SERVICES_RURAL: + return HeatSystemType.RURAL + else: + raise RuntimeError(f"Invalid heat system: {self}") + + @property + def sector(self) -> HeatSector: + """ + Returns the sector of the heat system. + + Returns + ------- + HeatSector + The sector of the heat system. + """ + if ( + self == HeatSystem.RESIDENTIAL_RURAL + or self == HeatSystem.RESIDENTIAL_URBAN_DECENTRAL + ): + return HeatSector.RESIDENTIAL + elif ( + self == HeatSystem.SERVICES_RURAL + or self == HeatSystem.SERVICES_URBAN_DECENTRAL + ): + return HeatSector.SERVICES + else: + "tot" + + @property + def is_rural(self) -> bool: + """ + Returns whether the heat system is for rural areas. + + Returns + ------- + bool + True if the heat system is for rural areas, False otherwise. + """ + if self == HeatSystem.RESIDENTIAL_RURAL or self == HeatSystem.SERVICES_RURAL: + return True + else: + return False + + @property + def is_urban_decentral(self) -> bool: + """ + Returns whether the heat system is for urban decentralized areas. + + Returns + ------- + bool + True if the heat system is for urban decentralized areas, False otherwise. + """ + if ( + self == HeatSystem.RESIDENTIAL_URBAN_DECENTRAL + or self == HeatSystem.SERVICES_URBAN_DECENTRAL + ): + return True + else: + return False + + @property + def is_urban(self) -> bool: + """ + Returns whether the heat system is for urban areas. + + Returns + ------- + bool True if the heat system is for urban areas, False otherwise. + """ + return not self.is_rural + + def heat_demand_weighting(self, urban_fraction=None, dist_fraction=None) -> float: + """ + Calculates the heat demand weighting based on urban fraction and + distribution fraction. + + Parameters + ---------- + urban_fraction : float, optional + The fraction of urban heat demand. + dist_fraction : float, optional + The fraction of distributed heat demand. + + Returns + ------- + float + The heat demand weighting. + + Raises + ------ + RuntimeError + If the heat system is invalid. + """ + if "rural" in self.value: + return 1 - urban_fraction + elif "urban central" in self.value: + return dist_fraction + elif "urban decentral" in self.value: + return urban_fraction - dist_fraction + else: + raise RuntimeError(f"Invalid heat system: {self}") + + def heat_pump_costs_name(self, heat_source: str) -> str: + """ + Generates the name for the heat pump costs based on the heat source and + system. + Used to retrieve data from `technology-data `. + + Parameters + ---------- + heat_source : str + The heat source. + + Returns + ------- + str + The name for the heat pump costs. + """ + return f"{self.central_or_decentral} {heat_source}-sourced heat pump" + + @property + def resistive_heater_costs_name(self) -> str: + """ + Generates the name for the resistive heater costs based on the heat + system. + Used to retrieve data from `technology-data `. + + Returns + ------- + str + The name for the heater costs. + """ + return f"{self.central_or_decentral} resistive heater" + + @property + def gas_boiler_costs_name(self) -> str: + """ + Generates the name for the gas boiler costs based on the heat system. + Used to retrieve data from `technology-data `. + + Returns + ------- + str + The name for the gas boiler costs. + """ + return f"{self.central_or_decentral} gas boiler" + + @property + def oil_boiler_costs_name(self) -> str: + """ + Generates the name for the oil boiler costs based on the heat system. + Used to retrieve data from `technology-data `. + + Returns + ------- + str + The name for the oil boiler costs. + """ + return "decentral oil boiler" diff --git a/scripts/definitions/heat_system_type.py b/scripts/definitions/heat_system_type.py new file mode 100644 index 00000000..5bd1bee3 --- /dev/null +++ b/scripts/definitions/heat_system_type.py @@ -0,0 +1,35 @@ +# -*- coding: utf-8 -*- +# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors +# +# SPDX-License-Identifier: MIT + +from enum import Enum + + +class HeatSystemType(Enum): + """ + Enumeration representing different types of heat systems. + """ + + URBAN_CENTRAL = "urban central" + URBAN_DECENTRAL = "urban decentral" + RURAL = "rural" + + def __str__(self) -> str: + """ + Returns the string representation of the heat system type. + + Returns: + str: The string representation of the heat system type. + """ + return self.value + + @property + def is_central(self) -> bool: + """ + Returns whether the heat system type is central. + + Returns: + bool: True if the heat system type is central, False otherwise. + """ + return self == HeatSystemType.URBAN_CENTRAL diff --git a/scripts/prepare_links_p_nom.py b/scripts/prepare_links_p_nom.py deleted file mode 100644 index 7c1ed211..00000000 --- a/scripts/prepare_links_p_nom.py +++ /dev/null @@ -1,95 +0,0 @@ -#!/usr/bin/env python -# -*- coding: utf-8 -*- - -# SPDX-FileCopyrightText: : 2017-2024 The PyPSA-Eur Authors -# -# SPDX-License-Identifier: MIT -""" -Extracts capacities of HVDC links from `Wikipedia. - -`_. - -Relevant Settings ------------------ - -.. code:: yaml - - enable: - prepare_links_p_nom: - -.. seealso:: - Documentation of the configuration file ``config/config.yaml`` at - :ref:`toplevel_cf` - -Inputs ------- - -*None* - -Outputs -------- - -- ``data/links_p_nom.csv``: A plain download of https://en.wikipedia.org/wiki/List_of_HVDC_projects#Europe plus extracted coordinates. - -Description ------------ - -*None* -""" - -import logging - -import pandas as pd -from _helpers import configure_logging, set_scenario_config - -logger = logging.getLogger(__name__) - - -def multiply(s): - return s.str[0].astype(float) * s.str[1].astype(float) - - -def extract_coordinates(s): - regex = ( - r"(\d{1,2})°(\d{1,2})′(\d{1,2})″(N|S) " r"(\d{1,2})°(\d{1,2})′(\d{1,2})″(E|W)" - ) - e = s.str.extract(regex, expand=True) - lat = ( - e[0].astype(float) + (e[1].astype(float) + e[2].astype(float) / 60.0) / 60.0 - ) * e[3].map({"N": +1.0, "S": -1.0}) - lon = ( - e[4].astype(float) + (e[5].astype(float) + e[6].astype(float) / 60.0) / 60.0 - ) * e[7].map({"E": +1.0, "W": -1.0}) - return lon, lat - - -if __name__ == "__main__": - if "snakemake" not in globals(): - from _helpers import mock_snakemake # rule must be enabled in config - - snakemake = mock_snakemake("prepare_links_p_nom", simpl="") - configure_logging(snakemake) - set_scenario_config(snakemake) - - links_p_nom = pd.read_html( - "https://en.wikipedia.org/wiki/List_of_HVDC_projects", header=0, match="SwePol" - )[0] - - mw = "Power (MW)" - m_b = links_p_nom[mw].str.contains("x").fillna(False) - - links_p_nom.loc[m_b, mw] = links_p_nom.loc[m_b, mw].str.split("x").pipe(multiply) - links_p_nom[mw] = ( - links_p_nom[mw].str.extract("[-/]?([\d.]+)", expand=False).astype(float) - ) - - links_p_nom["x1"], links_p_nom["y1"] = extract_coordinates( - links_p_nom["Converterstation 1"] - ) - links_p_nom["x2"], links_p_nom["y2"] = extract_coordinates( - links_p_nom["Converterstation 2"] - ) - - links_p_nom.dropna(subset=["x1", "y1", "x2", "y2"]).to_csv( - snakemake.output[0], index=False - ) diff --git a/scripts/prepare_sector_network.py b/scripts/prepare_sector_network.py index 5d61821c..bb850c3a 100644 --- a/scripts/prepare_sector_network.py +++ b/scripts/prepare_sector_network.py @@ -37,6 +37,10 @@ from pypsa.geo import haversine_pts from pypsa.io import import_components_from_dataframe from scipy.stats import beta +from scripts.definitions.heat_sector import HeatSector +from scripts.definitions.heat_system import HeatSystem +from scripts.definitions.heat_system_type import HeatSystemType + spatial = SimpleNamespace() logger = logging.getLogger(__name__) @@ -56,19 +60,27 @@ def define_spatial(nodes, options): # biomass spatial.biomass = SimpleNamespace() + spatial.msw = SimpleNamespace() if options.get("biomass_spatial", options["biomass_transport"]): spatial.biomass.nodes = nodes + " solid biomass" + spatial.biomass.bioliquids = nodes + " bioliquids" spatial.biomass.locations = nodes spatial.biomass.industry = nodes + " solid biomass for industry" spatial.biomass.industry_cc = nodes + " solid biomass for industry CC" + spatial.msw.nodes = nodes + " municipal solid waste" + spatial.msw.locations = nodes else: spatial.biomass.nodes = ["EU solid biomass"] + spatial.biomass.bioliquids = ["EU unsustainable bioliquids"] spatial.biomass.locations = ["EU"] spatial.biomass.industry = ["solid biomass for industry"] spatial.biomass.industry_cc = ["solid biomass for industry CC"] + spatial.msw.nodes = ["EU municipal solid waste"] + spatial.msw.locations = ["EU"] spatial.biomass.df = pd.DataFrame(vars(spatial.biomass), index=nodes) + spatial.msw.df = pd.DataFrame(vars(spatial.msw), index=nodes) # co2 @@ -2102,7 +2114,7 @@ def build_heat_demand(n): .unstack(level=1) ) - sectors = ["residential", "services"] + sectors = [sector.value for sector in HeatSector] uses = ["water", "space"] heat_demand = {} @@ -2130,10 +2142,21 @@ def build_heat_demand(n): return heat_demand -def add_heat(n, costs): +def add_heat(n: pypsa.Network, costs: pd.DataFrame, cop: xr.DataArray): + """ + Add heat sector to the network. + + Parameters: + n (pypsa.Network): The PyPSA network object. + costs (pd.DataFrame): DataFrame containing cost information. + cop (xr.DataArray): DataArray containing coefficient of performance (COP) values. + + Returns: + None + """ logger.info("Add heat sector") - sectors = ["residential", "services"] + sectors = [sector.value for sector in HeatSector] heat_demand = build_heat_demand(n) @@ -2152,23 +2175,6 @@ def add_heat(n, costs): for sector in sectors: heat_demand[sector + " space"] = (1 - dE) * heat_demand[sector + " space"] - heat_systems = [ - "residential rural", - "services rural", - "residential urban decentral", - "services urban decentral", - "urban central", - ] - - cop = { - "air": xr.open_dataarray(snakemake.input.cop_air_total) - .to_pandas() - .reindex(index=n.snapshots), - "ground": xr.open_dataarray(snakemake.input.cop_soil_total) - .to_pandas() - .reindex(index=n.snapshots), - } - if options["solar_thermal"]: solar_thermal = ( xr.open_dataarray(snakemake.input.solar_thermal_total) @@ -2178,31 +2184,34 @@ def add_heat(n, costs): # 1e3 converts from W/m^2 to MW/(1000m^2) = kW/m^2 solar_thermal = options["solar_cf_correction"] * solar_thermal / 1e3 - for name in heat_systems: - name_type = "central" if name == "urban central" else "decentral" + for ( + heat_system + ) in ( + HeatSystem + ): # this loops through all heat systems defined in _entities.HeatSystem - if name == "urban central": + if heat_system == HeatSystem.URBAN_CENTRAL: nodes = dist_fraction.index[dist_fraction > 0] else: nodes = pop_layout.index - n.add("Carrier", name + " heat") + n.add("Carrier", f"{heat_system} heat") n.madd( "Bus", - nodes + f" {name} heat", + nodes + f" {heat_system.value} heat", location=nodes, - carrier=name + " heat", + carrier=f"{heat_system.value} heat", unit="MWh_th", ) - if name == "urban central" and options.get("central_heat_vent"): + if heat_system == HeatSystem.URBAN_CENTRAL and options.get("central_heat_vent"): n.madd( "Generator", - nodes + f" {name} heat vent", - bus=nodes + f" {name} heat", + nodes + f" {heat_system} heat vent", + bus=nodes + f" {heat_system} heat", location=nodes, - carrier=name + " heat vent", + carrier=f"{heat_system} heat vent", p_nom_extendable=True, p_max_pu=0, p_min_pu=-1, @@ -2210,30 +2219,24 @@ def add_heat(n, costs): ) ## Add heat load + factor = heat_system.heat_demand_weighting( + urban_fraction=urban_fraction[nodes], dist_fraction=dist_fraction[nodes] + ) + if not heat_system == HeatSystem.URBAN_CENTRAL: + heat_load = ( + heat_demand[ + [ + heat_system.sector.value + " water", + heat_system.sector.value + " space", + ] + ] + .T.groupby(level=1) + .sum() + .T[nodes] + .multiply(factor) + ) - for sector in sectors: - # heat demand weighting - if "rural" in name: - factor = 1 - urban_fraction[nodes] - elif "urban central" in name: - factor = dist_fraction[nodes] - elif "urban decentral" in name: - factor = urban_fraction[nodes] - dist_fraction[nodes] - else: - raise NotImplementedError( - f" {name} not in " f"heat systems: {heat_systems}" - ) - - if sector in name: - heat_load = ( - heat_demand[[sector + " water", sector + " space"]] - .T.groupby(level=1) - .sum() - .T[nodes] - .multiply(factor) - ) - - if name == "urban central": + if heat_system == HeatSystem.URBAN_CENTRAL: heat_load = ( heat_demand.T.groupby(level=1) .sum() @@ -2246,20 +2249,25 @@ def add_heat(n, costs): n.madd( "Load", nodes, - suffix=f" {name} heat", - bus=nodes + f" {name} heat", - carrier=name + " heat", + suffix=f" {heat_system} heat", + bus=nodes + f" {heat_system} heat", + carrier=f"{heat_system} heat", p_set=heat_load, ) ## Add heat pumps - - heat_pump_types = ["air"] if "urban" in name else ["ground", "air"] - - for heat_pump_type in heat_pump_types: - costs_name = f"{name_type} {heat_pump_type}-sourced heat pump" + for heat_source in snakemake.params.heat_pump_sources[ + heat_system.system_type.value + ]: + costs_name = heat_system.heat_pump_costs_name(heat_source) efficiency = ( - cop[heat_pump_type][nodes] + cop.sel( + heat_system=heat_system.system_type.value, + heat_source=heat_source, + name=nodes, + ) + .to_pandas() + .reindex(index=n.snapshots) if options["time_dep_hp_cop"] else costs.at[costs_name, "efficiency"] ) @@ -2267,10 +2275,10 @@ def add_heat(n, costs): n.madd( "Link", nodes, - suffix=f" {name} {heat_pump_type} heat pump", + suffix=f" {heat_system} {heat_source} heat pump", bus0=nodes, - bus1=nodes + f" {name} heat", - carrier=f"{name} {heat_pump_type} heat pump", + bus1=nodes + f" {heat_system} heat", + carrier=f"{heat_system} {heat_source} heat pump", efficiency=efficiency, capital_cost=costs.at[costs_name, "efficiency"] * costs.at[costs_name, "fixed"] @@ -2280,59 +2288,65 @@ def add_heat(n, costs): ) if options["tes"]: - n.add("Carrier", name + " water tanks") + n.add("Carrier", f"{heat_system} water tanks") n.madd( "Bus", - nodes + f" {name} water tanks", + nodes + f" {heat_system} water tanks", location=nodes, - carrier=name + " water tanks", + carrier=f"{heat_system} water tanks", unit="MWh_th", ) n.madd( "Link", - nodes + f" {name} water tanks charger", - bus0=nodes + f" {name} heat", - bus1=nodes + f" {name} water tanks", + nodes + f" {heat_system} water tanks charger", + bus0=nodes + f" {heat_system} heat", + bus1=nodes + f" {heat_system} water tanks", efficiency=costs.at["water tank charger", "efficiency"], - carrier=name + " water tanks charger", + carrier=f"{heat_system} water tanks charger", p_nom_extendable=True, ) n.madd( "Link", - nodes + f" {name} water tanks discharger", - bus0=nodes + f" {name} water tanks", - bus1=nodes + f" {name} heat", - carrier=name + " water tanks discharger", + nodes + f" {heat_system} water tanks discharger", + bus0=nodes + f" {heat_system} water tanks", + bus1=nodes + f" {heat_system} heat", + carrier=f"{heat_system} water tanks discharger", efficiency=costs.at["water tank discharger", "efficiency"], p_nom_extendable=True, ) - tes_time_constant_days = options["tes_tau"][name_type] + tes_time_constant_days = options["tes_tau"][ + heat_system.central_or_decentral + ] n.madd( "Store", - nodes + f" {name} water tanks", - bus=nodes + f" {name} water tanks", + nodes + f" {heat_system} water tanks", + bus=nodes + f" {heat_system} water tanks", e_cyclic=True, e_nom_extendable=True, - carrier=name + " water tanks", + carrier=f"{heat_system} water tanks", standing_loss=1 - np.exp(-1 / 24 / tes_time_constant_days), - capital_cost=costs.at[name_type + " water tank storage", "fixed"], - lifetime=costs.at[name_type + " water tank storage", "lifetime"], + capital_cost=costs.at[ + heat_system.central_or_decentral + " water tank storage", "fixed" + ], + lifetime=costs.at[ + heat_system.central_or_decentral + " water tank storage", "lifetime" + ], ) if options["resistive_heaters"]: - key = f"{name_type} resistive heater" + key = f"{heat_system.central_or_decentral} resistive heater" n.madd( "Link", - nodes + f" {name} resistive heater", + nodes + f" {heat_system} resistive heater", bus0=nodes, - bus1=nodes + f" {name} heat", - carrier=name + " resistive heater", + bus1=nodes + f" {heat_system} heat", + carrier=f"{heat_system} resistive heater", efficiency=costs.at[key, "efficiency"], capital_cost=costs.at[key, "efficiency"] * costs.at[key, "fixed"] @@ -2342,16 +2356,16 @@ def add_heat(n, costs): ) if options["boilers"]: - key = f"{name_type} gas boiler" + key = f"{heat_system.central_or_decentral} gas boiler" n.madd( "Link", - nodes + f" {name} gas boiler", + nodes + f" {heat_system} gas boiler", p_nom_extendable=True, bus0=spatial.gas.df.loc[nodes, "nodes"].values, - bus1=nodes + f" {name} heat", + bus1=nodes + f" {heat_system} heat", bus2="co2 atmosphere", - carrier=name + " gas boiler", + carrier=f"{heat_system} gas boiler", efficiency=costs.at[key, "efficiency"], efficiency2=costs.at["gas", "CO2 intensity"], capital_cost=costs.at[key, "efficiency"] @@ -2361,22 +2375,26 @@ def add_heat(n, costs): ) if options["solar_thermal"]: - n.add("Carrier", name + " solar thermal") + n.add("Carrier", f"{heat_system} solar thermal") n.madd( "Generator", nodes, - suffix=f" {name} solar thermal collector", - bus=nodes + f" {name} heat", - carrier=name + " solar thermal", + suffix=f" {heat_system} solar thermal collector", + bus=nodes + f" {heat_system} heat", + carrier=f"{heat_system} solar thermal", p_nom_extendable=True, - capital_cost=costs.at[name_type + " solar thermal", "fixed"] + capital_cost=costs.at[ + heat_system.central_or_decentral + " solar thermal", "fixed" + ] * overdim_factor, p_max_pu=solar_thermal[nodes], - lifetime=costs.at[name_type + " solar thermal", "lifetime"], + lifetime=costs.at[ + heat_system.central_or_decentral + " solar thermal", "lifetime" + ], ) - if options["chp"] and name == "urban central": + if options["chp"] and heat_system == HeatSystem.URBAN_CENTRAL: # add gas CHP; biomass CHP is added in biomass section n.madd( "Link", @@ -2433,16 +2451,20 @@ def add_heat(n, costs): lifetime=costs.at["central gas CHP", "lifetime"], ) - if options["chp"] and options["micro_chp"] and name != "urban central": + if ( + options["chp"] + and options["micro_chp"] + and heat_system.value != "urban central" + ): n.madd( "Link", - nodes + f" {name} micro gas CHP", + nodes + f" {heat_system} micro gas CHP", p_nom_extendable=True, bus0=spatial.gas.df.loc[nodes, "nodes"].values, bus1=nodes, - bus2=nodes + f" {name} heat", + bus2=nodes + f" {heat_system} heat", bus3="co2 atmosphere", - carrier=name + " micro gas CHP", + carrier=heat_system.value + " micro gas CHP", efficiency=costs.at["micro CHP", "efficiency"], efficiency2=costs.at["micro CHP", "efficiency-heat"], efficiency3=costs.at["gas", "CO2 intensity"], @@ -2478,7 +2500,7 @@ def add_heat(n, costs): ) / heat_demand.T.groupby(level=[1]).sum().T for name in n.loads[ - n.loads.carrier.isin([x + " heat" for x in heat_systems]) + n.loads.carrier.isin([x + " heat" for x in HeatSystem]) ].index: node = n.buses.loc[name, "location"] ct = pop_layout.loc[node, "ct"] @@ -2632,19 +2654,83 @@ def add_biomass(n, costs): biogas_potentials_spatial = biomass_potentials["biogas"].rename( index=lambda x: x + " biogas" ) + unsustainable_biogas_potentials_spatial = biomass_potentials[ + "unsustainable biogas" + ].rename(index=lambda x: x + " biogas") else: biogas_potentials_spatial = biomass_potentials["biogas"].sum() + unsustainable_biogas_potentials_spatial = biomass_potentials[ + "unsustainable biogas" + ].sum() if options.get("biomass_spatial", options["biomass_transport"]): solid_biomass_potentials_spatial = biomass_potentials["solid biomass"].rename( index=lambda x: x + " solid biomass" ) + msw_biomass_potentials_spatial = biomass_potentials[ + "municipal solid waste" + ].rename(index=lambda x: x + " municipal solid waste") + unsustainable_solid_biomass_potentials_spatial = biomass_potentials[ + "unsustainable solid biomass" + ].rename(index=lambda x: x + " solid biomass") + else: solid_biomass_potentials_spatial = biomass_potentials["solid biomass"].sum() + msw_biomass_potentials_spatial = biomass_potentials[ + "municipal solid waste" + ].sum() + unsustainable_solid_biomass_potentials_spatial = biomass_potentials[ + "unsustainable solid biomass" + ].sum() + + if options["regional_oil_demand"]: + unsustainable_liquid_biofuel_potentials_spatial = biomass_potentials[ + "unsustainable bioliquids" + ].rename(index=lambda x: x + " bioliquids") + else: + unsustainable_liquid_biofuel_potentials_spatial = biomass_potentials[ + "unsustainable bioliquids" + ].sum() n.add("Carrier", "biogas") n.add("Carrier", "solid biomass") + if ( + options["municipal_solid_waste"] + and not options["industry"] + and cf_industry["waste_to_energy"] + or cf_industry["waste_to_energy_cc"] + ): + logger.warning( + "Flag municipal_solid_waste can be only used with industry " + "sector waste to energy." + "Setting municipal_solid_waste=False." + ) + options["municipal_solid_waste"] = False + + if options["municipal_solid_waste"]: + + n.add("Carrier", "municipal solid waste") + + n.madd( + "Bus", + spatial.msw.nodes, + location=spatial.msw.locations, + carrier="municipal solid waste", + ) + + e_max_pu = pd.Series([1] * (len(n.snapshots) - 1) + [0], index=n.snapshots) + n.madd( + "Store", + spatial.msw.nodes, + bus=spatial.msw.nodes, + carrier="municipal solid waste", + e_nom=msw_biomass_potentials_spatial, + marginal_cost=0, # costs.at["municipal solid waste", "fuel"], + e_max_pu=e_max_pu, + e_initial=msw_biomass_potentials_spatial, + ) + n.madd( "Bus", spatial.gas.biogas, @@ -2729,6 +2815,81 @@ def add_biomass(n, costs): p_nom_extendable=True, ) + if biomass_potentials.filter(like="unsustainable").sum().sum() > 0: + + # Create timeseries to force usage of unsustainable potentials + e_max_pu = pd.DataFrame(1, index=n.snapshots, columns=spatial.gas.biogas) + e_max_pu.iloc[-1] = 0 + + n.madd( + "Store", + spatial.gas.biogas, + suffix=" unsustainable", + bus=spatial.gas.biogas, + carrier="unsustainable biogas", + e_nom=unsustainable_biogas_potentials_spatial, + marginal_cost=costs.at["biogas", "fuel"], + e_initial=unsustainable_biogas_potentials_spatial, + e_nom_extendable=False, + e_max_pu=e_max_pu, + ) + + e_max_pu = pd.DataFrame(1, index=n.snapshots, columns=spatial.biomass.nodes) + e_max_pu.iloc[-1] = 0 + + n.madd( + "Store", + spatial.biomass.nodes, + suffix=" unsustainable", + bus=spatial.biomass.nodes, + carrier="unsustainable solid biomass", + e_nom=unsustainable_solid_biomass_potentials_spatial, + marginal_cost=costs.at["fuelwood", "fuel"], + e_initial=unsustainable_solid_biomass_potentials_spatial, + e_nom_extendable=False, + e_max_pu=e_max_pu, + ) + + n.madd( + "Bus", + spatial.biomass.bioliquids, + location=spatial.biomass.locations, + carrier="unsustainable bioliquids", + unit="MWh_LHV", + ) + + e_max_pu = pd.DataFrame( + 1, index=n.snapshots, columns=spatial.biomass.bioliquids + ) + e_max_pu.iloc[-1] = 0 + + n.madd( + "Store", + spatial.biomass.bioliquids, + suffix=" unsustainable", + bus=spatial.biomass.bioliquids, + carrier="unsustainable bioliquids", + e_nom=unsustainable_liquid_biofuel_potentials_spatial, + marginal_cost=costs.at["biodiesel crops", "fuel"], + e_initial=unsustainable_liquid_biofuel_potentials_spatial, + e_nom_extendable=False, + e_max_pu=e_max_pu, + ) + + n.madd( + "Link", + spatial.biomass.bioliquids, + bus0=spatial.biomass.bioliquids, + bus1=spatial.oil.nodes, + bus2="co2 atmosphere", + carrier="unsustainable bioliquids", + efficiency=1, + efficiency2=-costs.at["solid biomass", "CO2 intensity"] + + costs.at["BtL", "CO2 stored"], + p_nom=unsustainable_liquid_biofuel_potentials_spatial, + marginal_cost=costs.at["BtL", "VOM"], + ) + n.madd( "Link", spatial.gas.biogas_to_gas, @@ -2800,6 +2961,19 @@ def add_biomass(n, costs): carrier="solid biomass transport", ) + if options["municipal_solid_waste"]: + n.madd( + "Link", + biomass_transport.index, + bus0=biomass_transport.bus0 + " municipal solid waste", + bus1=biomass_transport.bus1 + " municipal solid waste", + p_nom_extendable=False, + p_nom=5e4, + length=biomass_transport.length.values, + marginal_cost=biomass_transport.costs * biomass_transport.length.values, + carrier="municipal solid waste transport", + ) + elif options["biomass_spatial"]: # add artificial biomass generators at nodes which include transport costs transport_costs = pd.read_csv( @@ -2829,6 +3003,26 @@ def add_biomass(n, costs): type="operational_limit", ) + if options["municipal_solid_waste"]: + # Add municipal solid waste + n.madd( + "Generator", + spatial.msw.nodes, + bus=spatial.msw.nodes, + carrier="municipal solid waste", + p_nom=10000, + marginal_cost=0 # costs.at["municipal solid waste", "fuel"] + + bus_transport_costs * average_distance, + ) + n.add( + "GlobalConstraint", + "msw limit", + carrier_attribute="municipal solid waste", + sense="<=", + constant=biomass_potentials["municipal solid waste"].sum(), + type="operational_limit", + ) + # AC buses with district heating urban_central = n.buses.index[n.buses.carrier == "urban central heat"] if not urban_central.empty and options["chp"]: @@ -3420,27 +3614,23 @@ def add_industry(n, costs): if options["oil_boilers"]: nodes = pop_layout.index - for name in [ - "residential rural", - "services rural", - "residential urban decentral", - "services urban decentral", - ]: - n.madd( - "Link", - nodes + f" {name} oil boiler", - p_nom_extendable=True, - bus0=spatial.oil.nodes, - bus1=nodes + f" {name} heat", - bus2="co2 atmosphere", - carrier=f"{name} oil boiler", - efficiency=costs.at["decentral oil boiler", "efficiency"], - efficiency2=costs.at["oil", "CO2 intensity"], - capital_cost=costs.at["decentral oil boiler", "efficiency"] - * costs.at["decentral oil boiler", "fixed"] - * options["overdimension_individual_heating"], - lifetime=costs.at["decentral oil boiler", "lifetime"], - ) + for heat_system in HeatSystem: + if not heat_system == HeatSystem.URBAN_CENTRAL: + n.madd( + "Link", + nodes + f" {heat_system} oil boiler", + p_nom_extendable=True, + bus0=spatial.oil.nodes, + bus1=nodes + f" {heat_system} heat", + bus2="co2 atmosphere", + carrier=f"{heat_system} oil boiler", + efficiency=costs.at["decentral oil boiler", "efficiency"], + efficiency2=costs.at["oil", "CO2 intensity"], + capital_cost=costs.at["decentral oil boiler", "efficiency"] + * costs.at["decentral oil boiler", "fixed"] + * options["overdimension_individual_heating"], + lifetime=costs.at["decentral oil boiler", "lifetime"], + ) n.madd( "Link", @@ -3535,6 +3725,17 @@ def add_industry(n, costs): efficiency3=process_co2_per_naphtha, ) + if options.get("biomass", True) and options["municipal_solid_waste"]: + n.madd( + "Link", + spatial.msw.locations, + bus0=spatial.msw.nodes, + bus1=non_sequestered_hvc_locations, + carrier="municipal solid waste", + p_nom_extendable=True, + efficiency=1.0, + ) + n.madd( "Link", spatial.oil.demand_locations, @@ -4427,13 +4628,14 @@ def add_enhanced_geothermal(n, egs_potentials, egs_overlap, costs): # %% if __name__ == "__main__": if "snakemake" not in globals(): + from _helpers import mock_snakemake snakemake = mock_snakemake( "prepare_sector_network", simpl="", opts="", - clusters="1", + clusters="37", ll="vopt", sector_opts="", planning_horizons="2050", @@ -4492,7 +4694,7 @@ if __name__ == "__main__": add_land_transport(n, costs) if options["heating"]: - add_heat(n, costs) + add_heat(n=n, costs=costs, cop=xr.open_dataarray(snakemake.input.cop_profiles)) if options["biomass"]: add_biomass(n, costs)