# -*- coding: utf-8 -*-
# SPDX-FileCopyrightText: : 2020-2024 The PyPSA-Eur Authors
#
# SPDX-License-Identifier: MIT
"""
Prepares brownfield data from previous planning horizon.
"""
import logging
import numpy as np
import pandas as pd
import pypsa
import xarray as xr
from _helpers import (
configure_logging,
get_snapshots,
set_scenario_config,
update_config_from_wildcards,
)
from add_existing_baseyear import add_build_year_to_new_assets
from pypsa.clustering.spatial import normed_or_uniform
logger = logging.getLogger(__name__)
idx = pd.IndexSlice
def add_brownfield(n, n_p, year):
logger.info(f"Preparing brownfield for the year {year}")
# electric transmission grid set optimised capacities of previous as minimum
n.lines.s_nom_min = n_p.lines.s_nom_opt
dc_i = n.links[n.links.carrier == "DC"].index
n.links.loc[dc_i, "p_nom_min"] = n_p.links.loc[dc_i, "p_nom_opt"]
for c in n_p.iterate_components(["Link", "Generator", "Store"]):
attr = "e" if c.name == "Store" else "p"
# first, remove generators, links and stores that track
# CO2 or global EU values since these are already in n
n_p.mremove(c.name, c.df.index[c.df.lifetime == np.inf])
# remove assets whose build_year + lifetime <= year
n_p.mremove(c.name, c.df.index[c.df.build_year + c.df.lifetime <= year])
# remove assets if their optimized nominal capacity is lower than a threshold
# since CHP heat Link is proportional to CHP electric Link, make sure threshold is compatible
chp_heat = c.df.index[
(c.df[f"{attr}_nom_extendable"] & c.df.index.str.contains("urban central"))
& c.df.index.str.contains("CHP")
& c.df.index.str.contains("heat")
]
threshold = snakemake.params.threshold_capacity
if not chp_heat.empty:
threshold_chp_heat = (
threshold
* c.df.efficiency[chp_heat.str.replace("heat", "electric")].values
* c.df.p_nom_ratio[chp_heat.str.replace("heat", "electric")].values
/ c.df.efficiency[chp_heat].values
)
n_p.mremove(
c.name,
chp_heat[c.df.loc[chp_heat, f"{attr}_nom_opt"] < threshold_chp_heat],
)
n_p.mremove(
c.name,
c.df.index[
(c.df[f"{attr}_nom_extendable"] & ~c.df.index.isin(chp_heat))
& (c.df[f"{attr}_nom_opt"] < threshold)
],
)
# copy over assets but fix their capacity
c.df[f"{attr}_nom"] = c.df[f"{attr}_nom_opt"]
c.df[f"{attr}_nom_extendable"] = False
n.import_components_from_dataframe(c.df, c.name)
# copy time-dependent
selection = n.component_attrs[c.name].type.str.contains(
"series"
) & n.component_attrs[c.name].status.str.contains("Input")
for tattr in n.component_attrs[c.name].index[selection]:
n.import_series_from_dataframe(c.pnl[tattr], c.name, tattr)
# deal with gas network
pipe_carrier = ["gas pipeline"]
if snakemake.params.H2_retrofit:
# subtract the already retrofitted from today's gas grid capacity
h2_retrofitted_fixed_i = n.links[
(n.links.carrier == "H2 pipeline retrofitted")
& (n.links.build_year != year)
].index
gas_pipes_i = n.links[n.links.carrier.isin(pipe_carrier)].index
CH4_per_H2 = 1 / snakemake.params.H2_retrofit_capacity_per_CH4
fr = "H2 pipeline retrofitted"
to = "gas pipeline"
# today's pipe capacity
pipe_capacity = n.links.loc[gas_pipes_i, "p_nom"]
# already retrofitted capacity from gas -> H2
already_retrofitted = (
n.links.loc[h2_retrofitted_fixed_i, "p_nom"]
.rename(lambda x: x.split("-2")[0].replace(fr, to) + f"-{year}")
.groupby(level=0)
.sum()
)
remaining_capacity = pipe_capacity - CH4_per_H2 * already_retrofitted.reindex(
index=pipe_capacity.index
).fillna(0)
n.links.loc[gas_pipes_i, "p_nom"] = remaining_capacity
def disable_grid_expansion_if_limit_hit(n):
"""
Check if transmission expansion limit is already reached; then turn off.
In particular, this function checks if the total transmission
capital cost or volume implied by s_nom_min and p_nom_min are
numerically close to the respective global limit set in
n.global_constraints. If so, the nominal capacities are set to the
minimum and extendable is turned off; the corresponding global
constraint is then dropped.
"""
types = {"expansion_cost": "capital_cost", "volume_expansion": "length"}
for limit_type in types:
glcs = n.global_constraints.query(f"type == 'transmission_{limit_type}_limit'")
for name, glc in glcs.iterrows():
total_expansion = (
(
n.lines.query("s_nom_extendable")
.eval(f"s_nom_min * {types[limit_type]}")
.sum()
)
+ (
n.links.query("carrier == 'DC' and p_nom_extendable")
.eval(f"p_nom_min * {types[limit_type]}")
.sum()
)
).sum()
# Allow small numerical differences
if np.abs(glc.constant - total_expansion) / glc.constant < 1e-6:
logger.info(
f"Transmission expansion {limit_type} is already reached, disabling expansion and limit"
)
extendable_acs = n.lines.query("s_nom_extendable").index
n.lines.loc[extendable_acs, "s_nom_extendable"] = False
n.lines.loc[extendable_acs, "s_nom"] = n.lines.loc[
extendable_acs, "s_nom_min"
]
extendable_dcs = n.links.query(
"carrier == 'DC' and p_nom_extendable"
).index
n.links.loc[extendable_dcs, "p_nom_extendable"] = False
n.links.loc[extendable_dcs, "p_nom"] = n.links.loc[
extendable_dcs, "p_nom_min"
]
n.global_constraints.drop(name, inplace=True)
def adjust_renewable_profiles(n, input_profiles, params, year):
"""
Adjusts renewable profiles according to the renewable technology specified,
using the latest year below or equal to the selected year.
"""
# spatial clustering
cluster_busmap = pd.read_csv(snakemake.input.cluster_busmap, index_col=0).squeeze()
simplify_busmap = pd.read_csv(
snakemake.input.simplify_busmap, index_col=0
).squeeze()
clustermaps = simplify_busmap.map(cluster_busmap)
clustermaps.index = clustermaps.index.astype(str)
# temporal clustering
dr = get_snapshots(params["snapshots"], params["drop_leap_day"])
snapshotmaps = (
pd.Series(dr, index=dr).where(lambda x: x.isin(n.snapshots), pd.NA).ffill()
)
for carrier in params["carriers"]:
if carrier == "hydro":
continue
with xr.open_dataset(getattr(input_profiles, "profile_" + carrier)) as ds:
if ds.indexes["bus"].empty or "year" not in ds.indexes:
continue
closest_year = max(
(y for y in ds.year.values if y <= year), default=min(ds.year.values)
)
p_max_pu = (
ds["profile"]
.sel(year=closest_year)
.transpose("time", "bus")
.to_pandas()
)
# spatial clustering
weight = ds["weight"].sel(year=closest_year).to_pandas()
weight = weight.groupby(clustermaps).transform(normed_or_uniform)
p_max_pu = (p_max_pu * weight).T.groupby(clustermaps).sum().T
p_max_pu.columns = p_max_pu.columns + f" {carrier}"
# temporal_clustering
p_max_pu = p_max_pu.groupby(snapshotmaps).mean()
# replace renewable time series
n.generators_t.p_max_pu.loc[:, p_max_pu.columns] = p_max_pu
if __name__ == "__main__":
if "snakemake" not in globals():
from _helpers import mock_snakemake
snakemake = mock_snakemake(
"add_brownfield",
simpl="",
clusters="37",
opts="",
ll="v1.0",
sector_opts="168H-T-H-B-I-dist1",
planning_horizons=2030,
)
configure_logging(snakemake)
set_scenario_config(snakemake)
update_config_from_wildcards(snakemake.config, snakemake.wildcards)
logger.info(f"Preparing brownfield from the file {snakemake.input.network_p}")
year = int(snakemake.wildcards.planning_horizons)
n = pypsa.Network(snakemake.input.network)
adjust_renewable_profiles(n, snakemake.input, snakemake.params, year)
add_build_year_to_new_assets(n, year)
n_p = pypsa.Network(snakemake.input.network_p)
add_brownfield(n, n_p, year)
disable_grid_expansion_if_limit_hit(n)
n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
n.export_to_netcdf(snakemake.output[0])