Build renewable profiles

Adapted from pypsa-EUR by the pypsa China-PIK authors

Calculates for each clustered region the (i) installable capacity (based on land-use from :mod:determine_availability_matrix) (ii) the available generation time series (based on weather data) (iii) the average distanc from the node for onshore wind, AC-connected offshore wind, DC-connected offshore wind and solar PV generators.

Outputs

• resources/profile_{technology}.nc with the following structure

  =================== ==================== ===================================================== Field Dimensions Description =================== ==================== ===================================================== profile year, bus, bin, time the per unit hourly availability factors for each bus

  ---

  p_nom_max bus, bin maximal installable capacity at the bus (in MW)

  ---

  average_distance bus, bin average distance of units in the region to the grid bus for onshore techs and to the shoreline for offshore technologies (in km) =================== ==================== =====================================================

Description

This script functions at two main spatial resolutions: the resolution of the clustered network regions, and the resolution of the cutout grid cells for the weather data. Typically the weather data grid is finer than the network regions, so we have to work out the distribution of generators across the grid cells within each region. This is done by taking account of a combination of the available land at each grid cell (computed in :mod:determine_availability_matrix) and the capacity factor there.

Based on the availability matrix, the script first computes how much of the technology can be installed at each cutout grid cell. To compute the layout of generators in each clustered region, the installable potential in each grid cell is multiplied with the capacity factor at each grid cell. This is done since we assume more generators are installed at cells with a higher capacity factor.

Based on the average capacity factor, the potentials are further divided into a configurable number of resource classes (bins).

This layout is then used to compute the generation availability time series from the weather data cutout from atlite.

The maximal installable potential for the node (p_nom_max) is computed by adding up the installable potentials of the individual grid cells.

build_resource_classes(cutout, nbins, regions, capacity_factor, params)

Bin resources based on their capacity factor The number of bins can be dynamically reduced based on a min delta cf

Parameters:

Name	Type	Description	Default
cutout	Cutout	the atlite cutout	required
nbins	int	the number of bins	required
regions	GeoSeries	the regions	required
capacity_factor	(DataArray,)	the capacity factor	required
params	dict	the config for VREs	required

Returns:

Type	Description
DataArray	xr.DataArray: the mask for the resource classes
GeoSeries	gpd.GeoSeries: multi-indexed series [bus, bin]: geometry

Source code in workflow/scripts/build_renewable_profiles.py

def build_resource_classes(
    cutout: Cutout,
    nbins: int,
    regions: gpd.GeoSeries,
    capacity_factor: xr.DataArray,
    params: dict,
) -> tuple[xr.DataArray, gpd.GeoSeries]:
    """Bin resources based on their capacity factor
    The number of bins can be dynamically reduced based on a min delta cf

    Args:
        cutout (Cutout): the atlite cutout
        nbins (int): the number of bins
        regions (gpd.GeoSeries): the regions
        capacity_factor (xr.DataArray,): the capacity factor
        params (dict): the config for VREs

    Returns:
        xr.DataArray: the mask for the resource classes
        gpd.GeoSeries: multi-indexed series [bus, bin]: geometry
    """
    resource_classes = params.get("resource_classes", {})
    nbins = resource_classes.get("n", 1)
    min_cf_delta = resource_classes.get("min_cf_delta", 0.0)
    buses = regions.index

    # indicator matrix for which cells touch which regions
    IndMat = np.ceil(cutout.availabilitymatrix(regions, ExclusionContainer()))
    cf_by_bus = capacity_factor * IndMat.where(IndMat > 0)

    epsilon = 1e-3
    cf_min, cf_max = (
        cf_by_bus.min(dim=["x", "y"]) - epsilon,
        cf_by_bus.max(dim=["x", "y"]) + epsilon,
    )

    # avoid binning resources that are very similar
    nbins_per_bus = [int(min(nbins, x)) for x in (cf_max - cf_min) // min_cf_delta]
    normed_bins = xr.DataArray(
        np.vstack(
            [np.hstack([[0] * (nbins - n), np.linspace(0, 1, n + 1)]) for n in nbins_per_bus]
        ),
        dims=["bus", "bin"],
        coords={"bus": regions.index},
    )
    bins = cf_min + (cf_max - cf_min) * normed_bins

    cf_by_bus_bin = cf_by_bus.expand_dims(bin=range(nbins))
    lower_edges = bins[:, :-1]
    upper_edges = bins[:, 1:]
    class_masks = (cf_by_bus_bin >= lower_edges) & (cf_by_bus_bin < upper_edges)

    if nbins == 1:
        bus_bin_mi = pd.MultiIndex.from_product([regions.index, [0]], names=["bus", "bin"])
        class_regions = regions.set_axis(bus_bin_mi)
        class_regions["cf"] = bins.to_series()
    else:
        grid = cutout.grid.set_index(["y", "x"])
        class_regions = {}
        for bus, bin_id in product(buses, range(nbins)):
            bus_bin_mask = (
                class_masks.sel(bus=bus, bin=bin_id).stack(spatial=["y", "x"]).to_pandas()
            )
            grid_cells = grid.loc[bus_bin_mask]
            geometry = grid_cells.intersection(regions.loc[bus, "geometry"]).union_all().buffer(0)
            class_regions[(bus, bin_id)] = geometry

        class_regions = gpd.GeoDataFrame(
            {"geometry": class_regions.values()},
            index=pd.MultiIndex.from_tuples(class_regions.keys(), names=["bus", "bin"]),
        )
        class_regions["cf"] = bins.to_series()

    return class_masks, class_regions

localize_cutout_time(cutout, drop_leap=True)

Localize the time to the local timezone

Parameters:

Name	Type	Description	Default
cutout	Cutout	the atlite cutout object	required
drop_leap	bool	drop 29th Feb. Defaults to True.	True

Returns:

Name	Type	Description
Cutout	Cutout	the updated cutout

Source code in workflow/scripts/build_renewable_profiles.py

def localize_cutout_time(cutout: Cutout, drop_leap=True) -> Cutout:
    """Localize the time to the local timezone

    Args:
        cutout (Cutout): the atlite cutout object
        drop_leap (bool, optional): drop 29th Feb. Defaults to True.

    Returns:
        Cutout: the updated cutout
    """

    data = cutout.data

    timestamps = pd.DatetimeIndex(data.time)
    # go from ECMWF/atlite UTC to local time
    ts_naive = timestamps.tz_localize("UTC").tz_convert(TIMEZONE).tz_localize(None)
    cutout.data = cutout.data.assign_coords(time=ts_naive)

    if drop_leap:
        data = cutout.data
        cutout.data = data.sel(time=~((data.time.dt.month == 2) & (data.time.dt.day == 29)))

    return cutout

prepare_resource_config(params, nprocesses, noprogress=True)

Parse the resource config (atlite calc config)

Parameters:

Name	Type	Description	Default
params	dict	the renewable options	required
nprocesses	int	the number or processes	required
noprogress	bool	whether to show progress bars	True

Returns:

Type	Description
(dict, dict)	the resource config for the atlite calcs, the turbine/panel models

Source code in workflow/scripts/build_renewable_profiles.py

def prepare_resource_config(params: dict, nprocesses: int, noprogress=True) -> tuple[dict]:
    """Parse the resource config (atlite calc config)

    Args:
        params (dict): the renewable options
        nprocesses (int): the number or processes
        noprogress (bool): whether to show progress bars

    Returns:
        (dict, dict): the resource config for the atlite calcs, the turbine/panel models
    """

    resource = params["resource"]  # pv panel params / wind turbine params
    resource["show_progress"] = not noprogress
    tech = "panel" if "panel" in resource else "turbine"

    # in case of multiple years
    models = resource[tech]
    if not isinstance(models, dict):
        models = {0: models}

    if nprocesses > 1:
        client = Client(n_workers=nprocesses, threads_per_worker=1)
        resource["dask_kwargs"] = {"scheduler": client}

    return resource, models