Merge branch 'year-specific-techs' of github.com:p-glaum/pypsa-eur into p-glaum-year-specific-techs

config/config.default.yaml

@@ -164,11 +164,14 @@ atlite:
 
 # docs in https://pypsa-eur.readthedocs.io/en/latest/configuration.html#renewable
 renewable:
+  year: 2020
   onwind:
     cutout: europe-2013-era5
     resource:
       method: wind
-      turbine: Vestas_V112_3MW
+      turbine:
+        2020: Vestas_V112_3MW
+        2030: NREL_ReferenceTurbine_2020ATB_5.5MW
       add_cutout_windspeed: true
     capacity_per_sqkm: 3
     # correction_factor: 0.93
@@ -187,7 +190,9 @@ renewable:
     cutout: europe-2013-era5
     resource:
       method: wind
-      turbine: NREL_ReferenceTurbine_2020ATB_5.5MW
+      turbine:
+        2020: NREL_ReferenceTurbine_5MW_offshore.yaml
+        2030: NREL_ReferenceTurbine_2020ATB_15MW_offshore
       add_cutout_windspeed: true
     capacity_per_sqkm: 2
     correction_factor: 0.8855
@@ -203,7 +208,10 @@ renewable:
     cutout: europe-2013-era5
     resource:
       method: wind
-      turbine: NREL_ReferenceTurbine_2020ATB_5.5MW
+      turbine:
+        2020: Vestas_V164_7MW_offshore
+        2025: NREL_ReferenceTurbine_2020ATB_15MW_offshore
+        2030: NREL_ReferenceTurbine_2020ATB_18MW_offshore
       add_cutout_windspeed: true
     capacity_per_sqkm: 2
     correction_factor: 0.8855
@@ -219,7 +227,8 @@ renewable:
     cutout: europe-2013-sarah
     resource:
       method: pv
-      panel: CSi
+      panel:
+        2020: CSi
       orientation:
         slope: 35.
         azimuth: 180.

rules/build_electricity.smk

@@ -261,6 +261,7 @@ rule build_renewable_profiles:
     params:
         snapshots={k: config["snapshots"][k] for k in ["start", "end", "inclusive"]},
         renewable=config["renewable"],
+        foresight=config["foresight"],
     input:
         **opt,
         base_network=RESOURCES + "networks/base.nc",

rules/solve_myopic.smk

@@ -86,6 +86,12 @@ rule add_brownfield:
         H2_retrofit_capacity_per_CH4=config["sector"]["H2_retrofit_capacity_per_CH4"],
         threshold_capacity=config["existing_capacities"]["threshold_capacity"],
     input:
+        **{
+            f"profile_{tech}": RESOURCES + f"profile_{tech}.nc"
+            for tech in config["electricity"]["renewable_carriers"]
+        },
+        simplify_busmap=RESOURCES + "busmap_elec_s{simpl}.csv",
+        cluster_busmap=RESOURCES + "busmap_elec_s{simpl}_{clusters}.csv",
         network=RESULTS
         + "prenetworks/elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}.nc",
         network_p=solved_previous_horizon,  #solved network at previous time step

scripts/add_brownfield.py

@@ -11,8 +11,10 @@ import logging
 import numpy as np
 import pandas as pd
 import pypsa
+import xarray as xr
 from _helpers import update_config_with_sector_opts
 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
@@ -143,6 +145,82 @@ def disable_grid_expansion_if_LV_limit_hit(n):
         n.global_constraints.drop("lv_limit", inplace=True)
 
 
+def adjust_renewable_profiles(n, input_profiles, config, year):
+    """
+    Adjusts renewable profiles according to the renewable technology specified.
+
+    If the planning horizon is not available, the closest year is used
+    instead.
+    """
+
+    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)
+    dr = pd.date_range(**config["snapshots"], freq="H")
+    snapshotmaps = (
+        pd.Series(dr, index=dr).where(lambda x: x.isin(n.snapshots), pd.NA).ffill()
+    )
+
+    for carrier in config["electricity"]["renewable_carriers"]:
+        if carrier == "hydro":
+            continue
+
+    clustermaps.index = clustermaps.index.astype(str)
+    dr = pd.date_range(**config["snapshots"], freq="H")
+    snapshotmaps = (
+        pd.Series(dr, index=dr).where(lambda x: x.isin(n.snapshots), pd.NA).ffill()
+    )
+    for carrier in config["electricity"]["renewable_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
+            if year in ds.indexes["year"]:
+                p_max_pu = (
+                    ds["year_profiles"]
+                    .sel(year=year)
+                    .transpose("time", "bus")
+                    .to_pandas()
+                )
+            else:
+                available_previous_years = [
+                    available_year
+                    for available_year in ds.indexes["year"]
+                    if available_year < year
+                ]
+                available_following_years = [
+                    available_year
+                    for available_year in ds.indexes["year"]
+                    if available_year > year
+                ]
+                if available_previous_years:
+                    closest_year = max(available_previous_years)
+                if available_following_years:
+                    closest_year = min(available_following_years)
+                logging.warning(
+                    f"Planning horizon {year} not in {carrier} profiles. Using closest year {closest_year} instead."
+                )
+                p_max_pu = (
+                    ds["year_profiles"]
+                    .sel(year=closest_year)
+                    .transpose("time", "bus")
+                    .to_pandas()
+                )
+            # spatial clustering
+            weight = ds["weight"].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
@@ -167,6 +245,8 @@ if __name__ == "__main__":
 
     n = pypsa.Network(snakemake.input.network)
 
+    adjust_renewable_profiles(n, snakemake.input, snakemake.config, year)
+
     add_build_year_to_new_assets(n, year)
 
     n_p = pypsa.Network(snakemake.input.network_p)

scripts/build_renewable_profiles.py

@@ -200,14 +200,25 @@ if __name__ == "__main__":
     if "snakemake" not in globals():
         from _helpers import mock_snakemake
 
-        snakemake = mock_snakemake("build_renewable_profiles", technology="solar")
+        snakemake = mock_snakemake("build_renewable_profiles", technology="onwind")
     configure_logging(snakemake)
 
     nprocesses = int(snakemake.threads)
     noprogress = snakemake.config["run"].get("disable_progressbar", True)
     noprogress = noprogress or not snakemake.config["atlite"]["show_progress"]
+    year = snakemake.params.renewable["year"]
+    foresight = snakemake.params.foresight
     params = snakemake.params.renewable[snakemake.wildcards.technology]
     resource = params["resource"]  # pv panel params / wind turbine params
+
+    year_dependent_techs = {
+        k: resource.get(k)
+        for k in ["panel", "turbine"]
+        if isinstance(resource.get(k), dict)
+    }
+    for key, techs in year_dependent_techs.items():
+        resource[key] = resource[key][year]
+
     correction_factor = params.get("correction_factor", 1.0)
     capacity_per_sqkm = params["capacity_per_sqkm"]
     snapshots = snakemake.params.snapshots
@@ -335,6 +346,29 @@ if __name__ == "__main__":
         **resource,
     )
 
+    if year_dependent_techs and foresight != "overnight":
+        for key, techs in year_dependent_techs.items():
+            year_profiles = list()
+            tech_profiles = dict()
+            tech_profiles[resource[key]] = profile
+            for year, tech in techs.items():
+                resource[key] = tech
+                if tech not in tech_profiles:
+                    tech_profiles[tech] = func(
+                        matrix=availability.stack(spatial=["y", "x"]),
+                        layout=layout,
+                        index=buses,
+                        per_unit=True,
+                        return_capacity=False,
+                        **resource,
+                    )
+                year_profile = tech_profiles[tech]
+                year_profile = year_profile.expand_dims({"year": [year]}).rename(
+                    "year_profiles"
+                )
+                year_profiles.append(year_profile)
+        year_profiles = xr.merge(year_profiles)
+
     duration = time.time() - start
     logger.info(
         f"Completed weighted capacity factor time series calculation ({duration:2.2f}s)"
@@ -373,6 +407,9 @@ if __name__ == "__main__":
         ]
     )
 
+    if year_dependent_techs:
+        ds = xr.merge([ds, year_profiles * correction_factor])
+
     if snakemake.wildcards.technology.startswith("offwind"):
         logger.info("Calculate underwater fraction of connections.")
         offshore_shape = gpd.read_file(snakemake.input["offshore_shapes"]).unary_union
@@ -396,6 +433,9 @@ if __name__ == "__main__":
     if "clip_p_max_pu" in params:
         min_p_max_pu = params["clip_p_max_pu"]
         ds["profile"] = ds["profile"].where(ds["profile"] >= min_p_max_pu, 0)
+        ds["year_profiles"] = ds["year_profiles"].where(
+            ds["year_profiles"] >= min_p_max_pu, 0
+        )
 
     ds.to_netcdf(snakemake.output.profile)