def prepare_network(config: dict, costs: pd.DataFrame, paths: dict) -> pypsa.Network:
"""Prepares/makes the network object for myopic mode according to config &
at 1 node per region/province
Args:
config (dict): the snakemake config
costs (pd.DataFrame): the costs dataframe (anualised capex and marginal costs)
paths (dict): dictionary of paths to input data
Returns:
pypsa.Network: the pypsa network object
"""
# derive the config
config["add_gas"] = (
True if [tech for tech in config["Techs"]["conv_techs"] if "gas" in tech] else False
)
config["add_coal"] = (
True if [tech for tech in config["Techs"]["conv_techs"] if "coal" in tech] else False
)
planning_horizons = snakemake.wildcards["planning_horizons"]
# empty network object
network = pypsa.Network()
# load graph
nodes = pd.Index(PROV_NAMES)
# set times
# make snapshots (drop leap days)
snapshot_cfg = config["snapshots"]
snapshots = make_periodic_snapshots(
year=snakemake.wildcards.planning_horizons,
freq=snapshot_cfg["freq"],
start_day_hour=snapshot_cfg["start"],
end_day_hour=snapshot_cfg["end"],
bounds=snapshot_cfg["bounds"],
# naive local tz
tz=None,
end_year=(
None
if not snapshot_cfg["end_year_plus1"]
else snakemake.wildcards.planning_horizons + 1
),
)
network.set_snapshots(snapshots.values)
network.snapshot_weightings[:] = config["snapshots"]["frequency"]
represented_hours = network.snapshot_weightings.sum()[0]
# TODO: what about leap years?
n_years = represented_hours / YEAR_HRS
prov_shapes = read_province_shapes(snakemake.input.province_shape)
prov_centroids = prov_shapes.to_crs("+proj=cea").centroid.to_crs(CRS)
# TODO split by carrier, make transparent
# add buses
for suffix in config["bus_suffix"]:
carrier = config["bus_carrier"][suffix]
add_buses(network, nodes, suffix, carrier, prov_centroids)
add_carriers(network, config, costs)
# ===== add load demand data =======
demand_path = snakemake.input.elec_load.replace("{planning_horizons}", cost_year)
with pd.HDFStore(demand_path, mode="r") as store:
load = store["load"].loc[network.snapshots] # MWh !!
load.columns = PROV_NAMES
network.add("Load", nodes, bus=nodes, p_set=load[nodes])
ws_carriers = [c for c in config["Techs"]["vre_techs"] if c.find("wind") >= 0 or c == "solar"]
add_wind_and_solar(network, ws_carriers, paths, planning_horizons, costs)
if config["heat_coupling"]:
central_fraction = pd.read_hdf(snakemake.input.central_fraction)
with pd.HDFStore(snakemake.input.heat_demand_profile, mode="r") as store:
heat_demand = store["heat_demand_profiles"]
# TODO fix this possilby not working
heat_demand.index = heat_demand.index.tz_localize(None)
heat_demand = heat_demand.loc[network.snapshots]
network.add(
"Load",
nodes,
suffix=" decentral heat",
bus=nodes + " decentral heat",
p_set=heat_demand[nodes].multiply(1 - central_fraction),
)
network.add(
"Load",
nodes,
suffix=" central heat",
bus=nodes + " central heat",
p_set=heat_demand[nodes].multiply(central_fraction),
)
# ====== add gas techs ======
if [tech for tech in config["Techs"]["conv_techs"] if "gas" in tech]:
# add converter from fuel source
network.add(
"Generator",
nodes,
suffix=" gas fuel",
bus=nodes + " gas",
carrier="gas",
p_nom_extendable=False,
p_nom=1e8,
marginal_cost=costs.at["OCGT", "fuel"],
)
network.add(
"Store",
nodes + " gas Store",
bus=nodes + " gas",
e_nom_extendable=False,
e_nom=1e8,
e_cyclic=True,
carrier="gas",
)
if "OCGT gas" in config["Techs"]["conv_techs"]:
network.add(
"Link",
nodes,
suffix=" OCGT",
bus0=nodes + " gas",
bus1=nodes,
carrier="OCGT gas",
marginal_cost=costs.at["OCGT", "efficiency"]
* costs.at["OCGT", "VOM"], # NB: VOM is per MWel
capital_cost=costs.at["OCGT", "efficiency"]
* costs.at["OCGT", "capital_cost"], # NB: capital cost is per MWel
p_nom_extendable=True,
efficiency=costs.at["OCGT", "efficiency"],
lifetime=costs.at["OCGT", "lifetime"],
)
if "gas boiler" in config["Techs"]["conv_techs"] and config["heat_coupling"]:
for cat in [" decentral "]:
network.add(
"Link",
nodes + cat + "gas boiler",
p_nom_extendable=True,
bus0=nodes + " gas",
bus1=nodes + cat + "heat",
carrier="gas boiler",
efficiency=costs.at[cat.lstrip() + "gas boiler", "efficiency"],
marginal_cost=costs.at[cat.lstrip() + "gas boiler", "efficiency"]
* costs.at[cat.lstrip() + "gas boiler", "VOM"],
capital_cost=costs.at[cat.lstrip() + "gas boiler", "efficiency"]
* costs.at[cat.lstrip() + "gas boiler", "capital_cost"],
lifetime=costs.at[cat.lstrip() + "gas boiler", "lifetime"],
)
# TODO missing second bus?
if "CHP gas" in config["Techs"]["conv_techs"]:
network.add(
"Link",
nodes,
suffix=" central CHP gas generator",
bus0=nodes + " gas",
bus1=nodes,
carrier="CHP gas",
p_nom_extendable=True,
marginal_cost=costs.at["central gas CHP", "efficiency"]
* costs.at["central gas CHP", "VOM"], # NB: VOM is per MWel
capital_cost=costs.at["central gas CHP", "efficiency"]
* costs.at["central gas CHP", "capital_cost"], # NB: capital cost is per MWel
efficiency=costs.at["central gas CHP", "efficiency"],
p_nom_ratio=1.0,
c_b=costs.at["central gas CHP", "c_b"],
lifetime=costs.at["central gas CHP", "lifetime"],
)
network.add(
"Link",
nodes,
suffix=" central CHP gas boiler",
bus0=nodes + " gas",
bus1=nodes + " central heat",
carrier="CHP gas",
p_nom_extendable=True,
marginal_cost=costs.at["central gas CHP", "efficiency"]
* costs.at["central gas CHP", "VOM"], # NB: VOM is per MWel
efficiency=costs.at["central gas CHP", "efficiency"]
/ costs.at["central gas CHP", "c_v"],
lifetime=costs.at["central gas CHP", "lifetime"],
)
# TODO separate retrofit in config from coal power plant
if "coal power plant" in config["Techs"]["conv_techs"] and config["Techs"]["coal_ccs_retrofit"]:
network.add(
"Generator",
nodes,
suffix=" coal cc",
bus=nodes,
carrier="coal cc",
p_nom_extendable=True,
efficiency=costs.at["coal", "efficiency"],
marginal_cost=costs.at["coal", "marginal_cost"],
capital_cost=costs.at["coal", "capital_cost"]
+ costs.at["retrofit", "capital_cost"], # NB: capital cost is per MWel
lifetime=costs.at["coal", "lifetime"],
)
# TODO FIXME harcoded
for year in range(int(planning_horizons) - 25, 2021, 5):
network.add(
"Generator",
nodes,
suffix=" coal-" + str(year) + "-retrofit",
bus=nodes,
carrier="coal cc",
p_nom_extendable=True,
capital_cost=costs.at["coal", "capital_cost"]
+ costs.at["retrofit", "capital_cost"]
+ 2021
- year,
efficiency=costs.at["coal", "efficiency"],
lifetime=costs.at["coal", "lifetime"],
build_year=year,
marginal_cost=costs.at["coal", "marginal_cost"],
)
# ===== add coal techs =====
if [tech for tech in config["Techs"]["conv_techs"] if "coal" in tech]:
# TODO check if this is needed (added Ivan), add for gas too, also why is it node resolved?
# network.add(
# "Bus",
# nodes,
# suffix=" coal fuel",
# x=prov_centroids.x,
# y=prov_centroids.y,
# carrier="coal",
# )
network.add(
"Generator",
nodes + " coal fuel",
bus=nodes + " coal",
carrier="coal",
p_nom_extendable=False,
p_nom=1e8,
marginal_cost=costs.at["coal", "marginal_cost"],
)
if "coal boiler" in config["Techs"]["conv_techs"]:
for cat in [" decentral ", " central "]:
network.add(
"Link",
nodes + cat + "coal boiler",
p_nom_extendable=True,
bus0=nodes + " coal",
bus1=nodes + cat + "heat",
carrier="coal boiler",
efficiency=costs.at[cat.lstrip() + "coal boiler", "efficiency"],
marginal_cost=costs.at[cat.lstrip() + "coal boiler", "efficiency"]
* costs.at[cat.lstrip() + "coal boiler", "VOM"],
capital_cost=costs.at[cat.lstrip() + "coal boiler", "efficiency"]
* costs.at[cat.lstrip() + "coal boiler", "capital_cost"],
lifetime=costs.at[cat.lstrip() + "coal boiler", "lifetime"],
)
# TODO missing second bus?
if "CHP coal" in config["Techs"]["conv_techs"]:
network.add(
"Link",
nodes,
suffix=" central CHP coal generator",
bus0=nodes + " coal",
bus1=nodes,
carrier="CHP coal",
p_nom_extendable=True,
marginal_cost=costs.at["central coal CHP", "efficiency"]
* costs.at["central coal CHP", "VOM"], # NB: VOM is per MWel
capital_cost=costs.at["central coal CHP", "efficiency"]
* costs.at["central coal CHP", "capital_cost"], # NB: capital cost is per MWel
efficiency=costs.at["central coal CHP", "efficiency"],
p_nom_ratio=1.0,
c_b=costs.at["central coal CHP", "c_b"],
lifetime=costs.at["central coal CHP", "lifetime"],
)
network.add(
"Link",
nodes,
suffix=" central CHP coal boiler",
bus0=nodes + " coal",
bus1=nodes + " central heat",
carrier="CHP coal",
p_nom_extendable=True,
marginal_cost=costs.at["central coal CHP", "efficiency"]
* costs.at["central coal CHP", "VOM"], # NB: VOM is per MWel
efficiency=costs.at["central coal CHP", "efficiency"]
/ costs.at["central coal CHP", "c_v"],
lifetime=costs.at["central coal CHP", "lifetime"],
)
if config["add_biomass"]:
network.add(
"Bus",
nodes,
suffix=" biomass",
x=prov_centroids.x,
y=prov_centroids.y,
carrier="biomass",
)
biomass_potential = pd.read_hdf(snakemake.input.biomass_potential)
biomass_potential.index = biomass_potential.index + " biomass"
network.add(
"Store",
nodes + " biomass",
bus=nodes + " biomass",
e_nom_extendable=False,
e_nom=biomass_potential,
e_initial=biomass_potential,
carrier="biomass",
)
network.add("Carrier", "CO2", co2_emissions=0)
network.add(
"Bus",
nodes,
suffix=" CO2",
x=prov_centroids.x,
y=prov_centroids.y,
carrier="CO2",
)
network.add("Store", nodes + " CO2", bus=nodes + " CO2", carrier="CO2")
# normally taking away from carrier generates CO2, but here we are
# adding CO2 stored, so the emissions will point the other way ?
network.add("Carrier", "CO2 capture", co2_emissions=1)
network.add(
"Bus",
nodes,
suffix=" CO2 capture",
x=prov_centroids.x,
y=prov_centroids.y,
carrier="CO2 capture",
)
network.add(
"Store",
nodes + " CO2 capture",
bus=nodes + " CO2 capture",
e_nom_extendable=True,
carrier="CO2 capture",
)
network.add(
"Link",
nodes + " central biomass CHP capture",
bus0=nodes + " CO2",
bus1=nodes + " CO2 capture",
bus2=nodes,
p_nom_extendable=True,
carrier="CO2 capture",
efficiency=costs.at["biomass CHP capture", "capture_rate"],
efficiency2=-1
* costs.at["biomass CHP capture", "capture_rate"]
* costs.at["biomass CHP capture", "electricity-input"],
capital_cost=costs.at["biomass CHP capture", "capture_rate"]
* costs.at["biomass CHP capture", "capital_cost"],
lifetime=costs.at["biomass CHP capture", "lifetime"],
)
# TODO rmemoe hard coded
network.add(
"Link",
nodes + " central biomass CHP",
bus0=nodes + " biomass",
bus1=nodes,
bus2=nodes + " central heat",
bus3=nodes + " CO2",
p_nom_extendable=True,
carrier="biomass",
efficiency=costs.at["biomass CHP", "efficiency"],
efficiency2=costs.at["biomass CHP", "efficiency-heat"],
# 4187.0095385594495TWh equates to 0.79*(5.24/3.04) Gt CO2 # tCO2/MWh
# TODO centralise
efficiency3=0.32522269504651985,
capital_cost=costs.at["biomass CHP", "efficiency"]
* costs.at["biomass CHP", "capital_cost"],
marginal_cost=costs.at["biomass CHP", "efficiency"]
* costs.at["biomass CHP", "marginal_cost"]
+ costs.at["solid biomass", "fuel"],
lifetime=costs.at["biomass CHP", "lifetime"],
)
network.add(
"Link",
nodes + " decentral biomass boiler",
bus0=nodes + " biomass",
bus1=nodes + " decentral heat",
p_nom_extendable=True,
carrier="biomass",
efficiency=costs.at["biomass boiler", "efficiency"],
capital_cost=costs.at["biomass boiler", "efficiency"]
* costs.at["biomass boiler", "capital_cost"],
marginal_cost=costs.at["biomass boiler", "efficiency"]
* costs.at["biomass boiler", "marginal_cost"]
+ costs.at["biomass boiler", "pelletizing cost"]
+ costs.at["solid biomass", "fuel"],
lifetime=costs.at["biomass boiler", "lifetime"],
)
if config["add_hydro"]:
# load dams
df = pd.read_csv(config["hydro_dams"]["dams_path"], index_col=0)
points = df.apply(lambda row: Point(row.Lon, row.Lat), axis=1)
dams = gpd.GeoDataFrame(df, geometry=points, crs=CRS)
hourly_rng = pd.date_range(
config["hydro_dams"]["inflow_date_start"],
config["hydro_dams"]["inflow_date_end"],
freq="1h",
inclusive="left",
)
inflow = pd.read_pickle(config["hydro_dams"]["inflow_path"]).reindex(
hourly_rng, fill_value=0
)
inflow.columns = dams.index
# convert to naive local timezone
inflow.index = inflow.index.tz_localize("UTC").tz_convert(TIMEZONE).tz_localize(None)
inflow = inflow.loc[str(INFLOW_DATA_YR)]
inflow = shift_profile_to_planning_year(inflow, planning_horizons)
inflow = inflow.loc[network.snapshots]
water_consumption_factor = (
dams.loc[:, "Water_consumption_factor_avg"] * 1e3
) # m^3/KWh -> m^3/MWh
###
# # Add hydro stations as buses
network.add(
"Bus",
dams.index,
suffix=" station",
carrier="stations",
x=dams["geometry"].to_crs("+proj=cea").centroid.to_crs(CRS).x,
y=dams["geometry"].to_crs("+proj=cea").centroid.to_crs(CRS).y,
location=dams["Province"],
)
dam_buses = network.buses[network.buses.carrier == "stations"]
# ===== add hydro reservoirs as stores ======
initial_capacity = pd.read_pickle(config["hydro_dams"]["reservoir_initial_capacity_path"])
effective_capacity = pd.read_pickle(
config["hydro_dams"]["reservoir_effective_capacity_path"]
)
initial_capacity.index = dams.index
effective_capacity.index = dams.index
initial_capacity = initial_capacity / water_consumption_factor
effective_capacity = effective_capacity / water_consumption_factor
network.add(
"Store",
dams.index,
suffix=" reservoir",
bus=dam_buses.index,
e_nom=effective_capacity,
e_initial=initial_capacity,
e_cyclic=True,
marginal_cost=config["hyro"]["marginal_cost"]["reservoir"] # EUR/MWh"
)
# add hydro turbines to link stations to provinces
network.add(
"Link",
dams.index,
suffix=" turbines",
bus0=dam_buses.index,
bus1=dams["Province"],
carrier="hydroelectricity",
p_nom=10 * dams["installed_capacity_10MW"],
capital_cost=costs.at["hydro", "capital_cost"],
efficiency=1,
)
# add rivers to link station to station
bus0s = [
0,
21,
11,
19,
22,
29,
8,
40,
25,
1,
7,
4,
10,
15,
12,
20,
26,
6,
3,
39,
]
bus1s = [
5,
11,
19,
22,
32,
8,
40,
25,
35,
2,
4,
10,
9,
12,
20,
23,
6,
17,
14,
16,
]
for bus0, bus2 in list(zip(dams.index[bus0s], dam_buses.iloc[bus1s].index)):
# normal flow
network.links.at[bus0 + " turbines", "bus2"] = bus2
network.links.at[bus0 + " turbines", "efficiency2"] = 1.0
# spillage
for bus0, bus1 in list(zip(dam_buses.iloc[bus0s].index, dam_buses.iloc[bus1s].index)):
network.add(
"Link",
f"{bus0}-{bus1}" + " spillage",
bus0=bus0,
bus1=bus1,
p_nom_extendable=True,
)
dam_ends = [
dam
for dam in range(len(dams.index))
if (dam in bus1s and dam not in bus0s) or (dam not in bus0s + bus1s)
]
for bus0 in dam_buses.iloc[dam_ends].index:
network.add(
"Link",
bus0 + " spillage",
bus0=bus0,
bus1="Tibet",
p_nom_extendable=True,
efficiency=0.0,
)
# == add inflow as generators
# only feed into hydro stations which are the first of a cascade
inflow_stations = [dam for dam in range(len(dams.index)) if dam not in bus1s]
for inflow_station in inflow_stations:
# p_nom = 1 and p_max_pu & p_min_pu = p_pu, compulsory inflow
p_nom = (inflow / water_consumption_factor).iloc[:, inflow_station].max()
p_pu = (inflow / water_consumption_factor).iloc[:, inflow_station] / p_nom
network.add(
"Generator",
dams.index[inflow_station] + " inflow",
bus=dam_buses.iloc[inflow_station].name,
carrier="hydro_inflow",
p_max_pu=p_pu.clip(1.0e-6),
p_min_pu=p_pu.clip(1.0e-6),
p_nom=p_nom,
)
# p_nom*p_pu = XXX m^3 then use turbines efficiency to convert to power
# ======= add other existing hydro power
hydro_p_nom = pd.read_hdf(config["hydro_dams"]["p_nom_path"])
hydro_p_max_pu = pd.read_hdf(
config["hydro_dams"]["p_max_pu_path"],
key=config["hydro_dams"]["p_max_pu_key"],
).tz_localize(None)
hydro_p_max_pu = shift_profile_to_planning_year(hydro_p_max_pu, planning_horizons)
# sort buses (columns) otherwise stuff will break
hydro_p_max_pu.sort_index(axis=1, inplace=True)
hydro_p_max_pu = hydro_p_max_pu.loc[snapshots]
hydro_p_max_pu.index = network.snapshots
network.add(
"Generator",
nodes,
suffix=" hydroelectricity",
bus=nodes,
carrier="hydroelectricity",
p_nom=hydro_p_nom,
p_nom_min=hydro_p_nom,
p_nom_extendable=False,
capital_cost=costs.at["hydro", "capital_cost"],
p_max_pu=hydro_p_max_pu,
)
if config["add_H2"]:
network.add(
"Bus",
nodes,
suffix=" H2",
x=prov_centroids.x,
y=prov_centroids.y,
carrier="H2",
)
network.add(
"Link",
nodes + " H2 Electrolysis",
bus0=nodes,
bus1=nodes + " H2",
bus2=nodes + " central heat",
p_nom_extendable=True,
carrier="H2",
efficiency=costs.at["electrolysis", "efficiency"],
efficiency2=costs.at["electrolysis", "efficiency-heat"],
capital_cost=costs.at["electrolysis", "capital_cost"],
lifetime=costs.at["electrolysis", "lifetime"],
)
network.add(
"Link",
nodes + " central H2 CHP",
bus0=nodes + " H2",
bus1=nodes,
bus2=nodes + " central heat",
p_nom_extendable=True,
carrier="H2 CHP",
efficiency=costs.at["central hydrogen CHP", "efficiency"],
efficiency2=costs.at["central hydrogen CHP", "efficiency"]
/ costs.at["central hydrogen CHP", "c_b"],
capital_cost=costs.at["central hydrogen CHP", "efficiency"]
* costs.at["central hydrogen CHP", "capital_cost"],
lifetime=costs.at["central hydrogen CHP", "lifetime"],
)
# TODO fix hard coded
H2_under_nodes = pd.Index(
[
"Sichuan",
"Chongqing",
"Hubei",
"Jiangxi",
"Anhui",
"Jiangsu",
"Shandong",
"Guangdong",
]
)
H2_type1_nodes = nodes.difference(H2_under_nodes)
network.add(
"Store",
H2_under_nodes + " H2 Store",
bus=H2_under_nodes + " H2",
e_nom_extendable=True,
e_cyclic=True,
capital_cost=costs.at["hydrogen storage underground", "capital_cost"],
lifetime=costs.at["hydrogen storage underground", "lifetime"],
)
network.add(
"Store",
H2_type1_nodes + " H2 Store",
bus=H2_type1_nodes + " H2",
e_nom_extendable=True,
e_cyclic=True,
capital_cost=costs.at[
"hydrogen storage tank type 1 including compressor", "capital_cost"
],
lifetime=costs.at["hydrogen storage tank type 1 including compressor", "lifetime"],
)
if config["add_methanation"]:
network.add(
"Link",
nodes + " Sabatier",
bus0=nodes + " H2",
bus1=nodes + " gas",
p_nom_extendable=True,
carrier="Sabatier",
efficiency=costs.at["methanation", "efficiency"],
capital_cost=costs.at["methanation", "efficiency"]
* costs.at["methanation", "capital_cost"]
+ costs.at["direct air capture", "capital_cost"]
* costs.at["gas", "co2_emissions"]
* costs.at["methanation", "efficiency"],
# TODO fix hardcoded
marginal_cost=(400 - 5 * (int(cost_year) - 2020))
* costs.at["gas", "co2_emissions"]
* costs.at["methanation", "efficiency"],
lifetime=costs.at["methanation", "lifetime"],
)
if "nuclear" in config["Techs"]["vre_techs"]:
nuclear_nodes = pd.Index(NUCLEAR_EXTENDABLE)
network.add(
"Generator",
nuclear_nodes,
suffix=" nuclear",
p_nom_extendable=True,
p_min_pu=0.7,
bus=nuclear_nodes,
carrier="nuclear",
efficiency=costs.at["nuclear", "efficiency"],
capital_cost=costs.at["nuclear", "capital_cost"], # NB: capital cost is per MWel
marginal_cost=costs.at["nuclear", "marginal_cost"],
lifetime=costs.at["nuclear", "lifetime"],
)
if "heat pump" in config["Techs"]["vre_techs"]:
with pd.HDFStore(snakemake.input.cop_name, mode="r") as store:
ashp_cop = store["ashp_cop_profiles"]
ashp_cop.index = ashp_cop.index.tz_localize(None)
ashp_cop = shift_profile_to_planning_year(ashp_cop, planning_horizons)
gshp_cop = store["gshp_cop_profiles"]
gshp_cop.index = gshp_cop.index.tz_localize(None)
gshp_cop = shift_profile_to_planning_year(gshp_cop, planning_horizons)
ashp_cop = ashp_cop.loc[snapshots]
gshp_cop = gshp_cop.loc[snapshots]
for cat in [" decentral ", " central "]:
network.add(
"Link",
nodes,
suffix=cat + "heat pump",
bus0=nodes,
bus1=nodes + cat + "heat",
carrier="heat pump",
efficiency=(
ashp_cop[nodes]
if config["time_dep_hp_cop"]
else costs.at[cat.lstrip() + "air-sourced heat pump", "efficiency"]
),
capital_cost=costs.at[cat.lstrip() + "air-sourced heat pump", "efficiency"]
* costs.at[cat.lstrip() + "air-sourced heat pump", "capital_cost"],
marginal_cost=costs.at[cat.lstrip() + "air-sourced heat pump", "efficiency"]
* costs.at[cat.lstrip() + "air-sourced heat pump", "marginal_cost"],
p_nom_extendable=True,
lifetime=costs.at[cat.lstrip() + "air-sourced heat pump", "lifetime"],
)
# TODO not valid for decentral
network.add(
"Link",
nodes,
suffix=cat + " ground heat pump",
bus0=nodes,
bus1=nodes + cat + "heat",
carrier="heat pump",
efficiency=(
gshp_cop[nodes]
if config["time_dep_hp_cop"]
else costs.at["decentral ground-sourced heat pump", "efficiency"]
),
capital_cost=costs.at[cat.lstrip() + "ground-sourced heat pump", "efficiency"]
* costs.at["decentral ground-sourced heat pump", "capital_cost"],
marginal_cost=costs.at[cat.lstrip() + "ground-sourced heat pump", "efficiency"]
* costs.at[cat.lstrip() + "ground-sourced heat pump", "marginal_cost"],
p_nom_extendable=True,
lifetime=costs.at["decentral ground-sourced heat pump", "lifetime"],
)
if "resistive heater" in config["Techs"]["vre_techs"]:
for cat in [" decentral ", " central "]:
network.add(
"Link",
nodes + cat + "resistive heater",
bus0=nodes,
bus1=nodes + cat + "heat",
carrier="resistive heater",
efficiency=costs.at[cat.lstrip() + "resistive heater", "efficiency"],
capital_cost=costs.at[cat.lstrip() + "resistive heater", "efficiency"]
* costs.at[cat.lstrip() + "resistive heater", "capital_cost"],
marginal_cost=costs.at[cat.lstrip() + "resistive heater", "efficiency"]
* costs.at[cat.lstrip() + "resistive heater", "marginal_cost"],
p_nom_extendable=True,
lifetime=costs.at[cat.lstrip() + "resistive heater", "lifetime"],
)
if "solar thermal" in config["Techs"]["vre_techs"]:
# this is the amount of heat collected in W per m^2, accounting
# for efficiency
with pd.HDFStore(snakemake.input.solar_thermal_name, mode="r") as store:
# 1e3 converts from W/m^2 to MW/(1000m^2) = kW/m^2
solar_thermal = config["solar_cf_correction"] * store["solar_thermal_profiles"] / 1e3
solar_thermal.index = solar_thermal.index.tz_localize(None)
solar_thermal = shift_profile_to_planning_year(solar_thermal, planning_horizons)
solar_thermal = solar_thermal.loc[snapshots]
for cat in [" decentral ", " central "]:
network.add(
"Generator",
nodes,
suffix=cat + "solar thermal",
bus=nodes + cat + "heat",
carrier="solar thermal",
p_nom_extendable=True,
capital_cost=costs.at[cat.lstrip() + "solar thermal", "capital_cost"],
p_max_pu=solar_thermal[nodes].clip(1.0e-4),
lifetime=costs.at[cat.lstrip() + "solar thermal", "lifetime"],
)
if "water tanks" in config["Techs"]["store_techs"]:
for cat in [" decentral ", " central "]:
network.add(
"Bus",
nodes,
suffix=cat + "water tanks",
x=prov_centroids.x,
y=prov_centroids.y,
carrier="water tanks",
)
network.add(
"Link",
nodes + cat + "water tanks charger",
bus0=nodes + cat + "heat",
bus1=nodes + cat + "water tanks",
carrier="water tanks",
efficiency=costs.at["water tank charger", "efficiency"],
p_nom_extendable=True,
)
network.add(
"Link",
nodes + cat + "water tanks discharger",
bus0=nodes + cat + "water tanks",
bus1=nodes + cat + "heat",
carrier="water tanks",
efficiency=costs.at["water tank discharger", "efficiency"],
p_nom_extendable=True,
)
# [HP] 180 day time constant for centralised, 3 day for decentralised
tes_tau = config["water_tanks"]["tes_tau"][cat.strip()]
network.add(
"Store",
nodes + cat + "water tank",
bus=nodes + cat + "water tanks",
carrier="water tanks",
e_cyclic=True,
e_nom_extendable=True,
standing_loss=1 - np.exp(-1 / (24.0 * tes_tau)),
capital_cost=costs.at[cat.lstrip() + "water tank storage", "capital_cost"],
lifetime=costs.at[cat.lstrip() + "water tank storage", "lifetime"],
)
if "battery" in config["Techs"]["store_techs"]:
network.add(
"Bus",
nodes,
suffix=" battery",
x=prov_centroids.x,
y=prov_centroids.y,
carrier="battery",
)
network.add(
"Store",
nodes + " battery",
bus=nodes + " battery",
e_cyclic=True,
e_nom_extendable=True,
capital_cost=costs.at["battery storage", "capital_cost"],
lifetime=costs.at["battery storage", "lifetime"],
)
network.add(
"Link",
nodes + " battery charger",
bus0=nodes,
bus1=nodes + " battery",
efficiency=costs.at["battery inverter", "efficiency"] ** 0.5,
capital_cost=costs.at["battery inverter", "capital_cost"],
p_nom_extendable=True,
carrier="battery",
lifetime=costs.at["battery inverter", "lifetime"],
)
network.add(
"Link",
nodes + " battery discharger",
bus0=nodes + " battery",
bus1=nodes,
efficiency=costs.at["battery inverter", "efficiency"] ** 0.5,
marginal_cost=0.0,
carrier="battery discharger",
p_nom_extendable=True,
)
if "PHS" in config["Techs"]["store_techs"]:
# pure pumped hydro storage, fixed, 6h energy by default, no inflow
hydrocapa_df = pd.read_csv("resources/data/hydro/PHS_p_nom.csv", index_col=0)
phss = hydrocapa_df.index[hydrocapa_df["MW"] > 0].intersection(nodes)
if config["hydro"]["hydro_capital_cost"]:
cc = costs.at["PHS", "capital_cost"]
else:
cc = 0.0
network.add(
"StorageUnit",
phss,
suffix=" PHS",
bus=phss,
carrier="PHS",
p_nom_extendable=False,
p_nom=hydrocapa_df.loc[phss]["MW"],
p_nom_min=hydrocapa_df.loc[phss]["MW"],
max_hours=config["hydro"]["PHS_max_hours"],
efficiency_store=np.sqrt(costs.at["PHS", "efficiency"]),
efficiency_dispatch=np.sqrt(costs.at["PHS", "efficiency"]),
cyclic_state_of_charge=True,
capital_cost=cc,
marginal_cost=0.0,
)
# ============= add lines =========
# The lines are implemented according to the transport model (no KVL) with losses.
# This requires two directions
# see Neumann et al 10.1016/j.apenergy.2022.118859
# TODO make lossless optional (speed up)
if not config["no_lines"]:
edges = pd.read_csv(snakemake.input.edges, header=None)
lengths = NON_LIN_PATH_SCALING * np.array(
[
haversine(
[network.buses.at[name0, "x"], network.buses.at[name0, "y"]],
[network.buses.at[name1, "x"], network.buses.at[name1, "y"]],
)
for name0, name1 in edges[[0, 1]].values
]
)
cc = (
(config["line_cost_factor"] * lengths * [HVAC_cost_curve(len_) for len_ in lengths])
* LINE_SECURITY_MARGIN
* FOM_LINES
* n_years
* annuity(ECON_LIFETIME_LINES, config["costs"]["discountrate"])
)
network.add(
"Link",
edges[0] + "-" + edges[1],
bus0=edges[0].values,
bus1=edges[1].values,
suffix=" positive",
p_nom_extendable=True,
p_min_pu=0,
efficiency=config["transmission_efficiency"]["DC"]["efficiency_static"]
* config["transmission_efficiency"]["DC"]["efficiency_per_1000km"] ** (lengths / 1000),
length=lengths,
capital_cost=cc,
)
network.add(
"Link",
edges[1] + "-" + edges[0],
bus0=edges[1].values,
bus1=edges[0].values,
suffix=" reversed",
p_nom_extendable=True,
p_min_pu=0,
efficiency=config["transmission_efficiency"]["DC"]["efficiency_static"]
* config["transmission_efficiency"]["DC"]["efficiency_per_1000km"] ** (lengths / 1000),
length=lengths,
capital_cost=0,
)
if config["Techs"]["hydrogen_lines"]:
edges = pd.read_csv(snakemake.input.edges, header=None)
lengths = NON_LIN_PATH_SCALING * np.array(
[
haversine(
[network.buses.at[name0, "x"], network.buses.at[name0, "y"]],
[network.buses.at[name1, "x"], network.buses.at[name1, "y"]],
)
for name0, name1 in edges[[0, 1]].values
]
)
cc = costs.at["H2 (g) pipeline", "capital_cost"] * lengths
network.add(
"Link",
edges[0] + "-" + edges[1] + " H2 pipeline",
suffix=" positive",
bus0=edges[0].values + " H2",
bus1=edges[1].values + " H2",
bus2=edges[0].values,
p_nom_extendable=True,
p_nom=0,
p_nom_min=0,
p_min_pu=0,
efficiency=config["transmission_efficiency"]["H2 pipeline"]["efficiency_static"]
* config["transmission_efficiency"]["H2 pipeline"]["efficiency_per_1000km"]
** (lengths / 1000),
efficiency2=-config["transmission_efficiency"]["H2 pipeline"]["compression_per_1000km"]
* lengths
/ 1e3,
length=lengths,
lifetime=costs.at["H2 (g) pipeline", "lifetime"],
capital_cost=cc,
)
network.add(
"Link",
edges[1] + "-" + edges[0] + " H2 pipeline",
suffix=" reversed",
bus0=edges[1].values + " H2",
bus1=edges[0].values + " H2",
bus2=edges[1].values,
p_nom_extendable=True,
p_nom=0,
p_nom_min=0,
p_min_pu=0,
efficiency=config["transmission_efficiency"]["H2 pipeline"]["efficiency_static"]
* config["transmission_efficiency"]["H2 pipeline"]["efficiency_per_1000km"]
** (lengths / 1000),
efficiency2=-config["transmission_efficiency"]["H2 pipeline"]["compression_per_1000km"]
* lengths
/ 1e3,
length=lengths,
lifetime=costs.at["H2 (g) pipeline", "lifetime"],
capital_cost=0,
)
return network