Merge branch 'master' into complete-carrier-set

.pre-commit-config.yaml

@@ -39,7 +39,7 @@ repos:
 
   # Make docstrings PEP 257 compliant
 - repo: https://github.com/PyCQA/docformatter
-  rev: v1.6.4
+  rev: v1.6.5
   hooks:
   - id: docformatter
     args: ["--in-place", "--make-summary-multi-line", "--pre-summary-newline"]

.readthedocs.yml

@@ -4,8 +4,14 @@
 
 version: 2
 
+build:
+  os: ubuntu-22.04
+  tools:
+    python: "3.11"
+  apt_packages:
+  - graphviz
+
 python:
-  version: 3.8
   install:
   - requirements: doc/requirements.txt
-  system_packages: true
+  system_packages: false

config/config.default.yaml

@@ -623,9 +623,9 @@ clustering:
 solving:
   #tmpdir: "path/to/tmp"
   options:
-    formulation: kirchhoff
     clip_p_max_pu: 1.e-2
     load_shedding: false
+    transmission_losses: 0
     noisy_costs: true
     skip_iterations: true
     track_iterations: false

doc/conf.py

@@ -36,6 +36,7 @@ sys.path.insert(0, os.path.abspath("../scripts"))
 extensions = [
     #'sphinx.ext.autodoc',
     #'sphinx.ext.autosummary',
+    "myst_parser",
     "sphinx.ext.autosectionlabel",
     "sphinx.ext.intersphinx",
     "sphinx.ext.todo",

doc/configtables/solving.csv

@@ -1,7 +1,7 @@
 ,Unit,Values,Description
 options,,,
--- formulation,--,"Any of {'angles', 'kirchhoff', 'cycles', 'ptdf'}","Specifies which variant of linearized power flow formulations to use in the optimisation problem. Recommended is 'kirchhoff'. Explained in `this article <https://arxiv.org/abs/1704.01881>`_."
 -- load_shedding,bool/float,"{'true','false', float}","Add generators with very high marginal cost to simulate load shedding and avoid problem infeasibilities. If load shedding is a float, it denotes the marginal cost in EUR/kWh."
+-- transmission_losses,int,"[0-9]","Add piecewise linear approximation of transmission losses based on n tangents. Defaults to 0, which means losses are ignored."
 -- noisy_costs,bool,"{'true','false'}","Add random noise to marginal cost of generators by :math:`\mathcal{U}(0.009,0,011)` and capital cost of lines and links by :math:`\mathcal{U}(0.09,0,11)`."
 -- min_iterations,--,int,"Minimum number of solving iterations in between which resistance and reactence (``x/r``) are updated for branches according to ``s_nom_opt`` of the previous run."
 -- max_iterations,--,int,"Maximum number of solving iterations in between which resistance and reactence (``x/r``) are updated for branches according to ``s_nom_opt`` of the previous run."

doc/release_notes.rst

@@ -17,6 +17,8 @@ Upcoming Release
 
 * Renamed script file from PyPSA-EUR ``build_load_data`` to ``build_electricity_demand`` and ``retrieve_load_data`` to ``retrieve_electricity_demand``.
 
+* Fix docs readthedocs built
+
 * Add plain hydrogen turbine as additional re-electrification option besides
   hydrogen fuel cell. Add switches for both re-electrification options under
   ``sector: hydrogen_turbine:`` and ``sector: hydrogen_fuel_cell:``.
@@ -25,6 +27,13 @@ Upcoming Release
 
 * Additional tech_color are added to include previously unlisted carriers.
 
+* Remove ``vresutils`` dependency.
+
+* Add option to include a piecewise linear approximation of transmission losses,
+  e.g. by setting ``solving: options: transmission_losses: 2`` for an
+  approximation with two tangents.
+  
+  
 PyPSA-Eur 0.8.0 (18th March 2023)
 =================================
 

doc/requirements.txt

@@ -2,12 +2,13 @@
 #
 # SPDX-License-Identifier: CC0-1.0
 
+setuptools
 sphinx
 sphinx_book_theme
 sphinxcontrib-bibtex
+myst-parser  # recommark is deprecated, https://stackoverflow.com/a/71660856/13573820
 
 pypsa
-vresutils>=0.3.1
 powerplantmatching>=0.5.5
 atlite>=0.2.9
 dask[distributed]

envs/environment.fixed.yaml

@@ -226,7 +226,7 @@ dependencies:
 - nspr=4.35
 - nss=3.88
 - numexpr=2.8.3
-- numpy=1.23.5
+- numpy=1.24
 - openjdk=17.0.3
 - openjpeg=2.5.0
 - openpyxl=3.1.0
@@ -378,4 +378,3 @@ dependencies:
   - highspy==1.5.0.dev0
   - pybind11==2.10.3
   - tsam==2.2.2
-  - vresutils==0.3.1

envs/environment.yaml

@@ -10,7 +10,7 @@ dependencies:
 - python>=3.8
 - pip
 
-- pypsa>=0.21.3
+- pypsa>=0.23
 - atlite>=0.2.9
 - dask
 
@@ -25,7 +25,7 @@ dependencies:
 - pytables
 - lxml
 - powerplantmatching>=0.5.5
-- numpy<1.24
+- numpy
 - pandas>=1.4
 - geopandas>=0.11.0
 - xarray
@@ -55,5 +55,4 @@ dependencies:
 - rasterio!=1.2.10
 
 - pip:
-  - vresutils>=0.3.1
   - tsam>=1.1.0

rules/solve_electricity.smk

@@ -14,8 +14,6 @@ rule solve_network:
         ),
         python=LOGS
         + "solve_network/elec_s{simpl}_{clusters}_ec_l{ll}_{opts}_python.log",
-        memory=LOGS
-        + "solve_network/elec_s{simpl}_{clusters}_ec_l{ll}_{opts}_memory.log",
     benchmark:
         BENCHMARKS + "solve_network/elec_s{simpl}_{clusters}_ec_l{ll}_{opts}"
     threads: 4
@@ -41,8 +39,6 @@ rule solve_operations_network:
         ),
         python=LOGS
         + "solve_operations_network/elec_s{simpl}_{clusters}_ec_l{ll}_{opts}_op_python.log",
-        memory=LOGS
-        + "solve_operations_network/elec_s{simpl}_{clusters}_ec_l{ll}_{opts}_op_memory.log",
     benchmark:
         (
             BENCHMARKS

rules/solve_myopic.smk

@@ -90,8 +90,6 @@ rule solve_sector_network_myopic:
         + "elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}_solver.log",
         python=LOGS
         + "elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}_python.log",
-        memory=LOGS
-        + "elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}_memory.log",
     threads: 4
     resources:
         mem_mb=config["solving"]["mem"],

rules/solve_overnight.smk

@@ -21,8 +21,6 @@ rule solve_sector_network:
         + "elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}_solver.log",
         python=LOGS
         + "elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}_python.log",
-        memory=LOGS
-        + "elec_s{simpl}_{clusters}_l{ll}_{opts}_{sector_opts}_{planning_horizons}_memory.log",
     threads: config["solving"]["solver"].get("threads", 4)
     resources:
         mem_mb=config["solving"]["mem"],

scripts/add_electricity.py

@@ -85,16 +85,18 @@ It further adds extendable ``generators`` with **zero** capacity for
 """
 
 import logging
+from itertools import product
 
 import geopandas as gpd
 import numpy as np
 import pandas as pd
 import powerplantmatching as pm
 import pypsa
+import scipy.sparse as sparse
 import xarray as xr
 from _helpers import configure_logging, update_p_nom_max
 from powerplantmatching.export import map_country_bus
-from vresutils import transfer as vtransfer
+from shapely.prepared import prep
 
 idx = pd.IndexSlice
 
@@ -257,6 +259,21 @@ def load_powerplants(ppl_fn):
     )
 
 
+def shapes_to_shapes(orig, dest):
+    """
+    Adopted from vresutils.transfer.Shapes2Shapes()
+    """
+    orig_prepped = list(map(prep, orig))
+    transfer = sparse.lil_matrix((len(dest), len(orig)), dtype=float)
+
+    for i, j in product(range(len(dest)), range(len(orig))):
+        if orig_prepped[j].intersects(dest[i]):
+            area = orig[j].intersection(dest[i]).area
+            transfer[i, j] = area / dest[i].area
+
+    return transfer
+
+
 def attach_load(n, regions, load, nuts3_shapes, countries, scaling=1.0):
     substation_lv_i = n.buses.index[n.buses["substation_lv"]]
     regions = gpd.read_file(regions).set_index("name").reindex(substation_lv_i)
@@ -273,9 +290,7 @@ def attach_load(n, regions, load, nuts3_shapes, countries, scaling=1.0):
             return pd.DataFrame({group.index[0]: l})
         else:
             nuts3_cntry = nuts3.loc[nuts3.country == cntry]
-            transfer = vtransfer.Shapes2Shapes(
+            transfer = shapes_to_shapes(group, nuts3_cntry.geometry).T.tocsr()
-                group, nuts3_cntry.geometry, normed=False
-            ).T.tocsr()
             gdp_n = pd.Series(
                 transfer.dot(nuts3_cntry["gdp"].fillna(1.0).values), index=group.index
             )
@@ -444,7 +459,9 @@ def attach_conventional_generators(
             if f"conventional_{carrier}_{attr}" in conventional_inputs:
                 # Values affecting generators of technology k country-specific
                 # First map generator buses to countries; then map countries to p_max_pu
-                values = pd.read_csv(values, index_col=0).iloc[:, 0]
+                values = pd.read_csv(
+                    snakemake.input[f"conventional_{carrier}_{attr}"], index_col=0
+                ).iloc[:, 0]
                 bus_values = n.buses.country.map(values)
                 n.generators[attr].update(
                     n.generators.loc[idx].bus.map(bus_values).dropna()

scripts/prepare_sector_network.py

@@ -22,14 +22,13 @@ from _helpers import (
     override_component_attrs,
     update_config_with_sector_opts,
 )
-from add_electricity import sanitize_carriers
+from add_electricity import sanitize_carriers, calculate_annuity
 from build_energy_totals import build_co2_totals, build_eea_co2, build_eurostat_co2
 from networkx.algorithms import complement
 from networkx.algorithms.connectivity.edge_augmentation import k_edge_augmentation
 from pypsa.geo import haversine_pts
 from pypsa.io import import_components_from_dataframe
 from scipy.stats import beta
-from vresutils.costdata import annuity
 
 logger = logging.getLogger(__name__)
 
@@ -743,7 +742,7 @@ def prepare_costs(cost_file, config, nyears):
     costs = costs.fillna(config["fill_values"])
 
     def annuity_factor(v):
-        return annuity(v["lifetime"], v["discount rate"]) + v["FOM"] / 100
+        return calculate_annuity(v["lifetime"], v["discount rate"]) + v["FOM"] / 100
 
     costs["fixed"] = [
         annuity_factor(v) * v["investment"] * nyears for i, v in costs.iterrows()
@@ -852,7 +851,7 @@ def add_wave(n, wave_cost_factor):
     capacity = pd.Series({"Attenuator": 750, "F2HB": 1000, "MultiPA": 600})
 
     # in EUR/MW
-    annuity_factor = annuity(25, 0.07) + 0.03
+    annuity_factor = calculate_annuity(25, 0.07) + 0.03
     costs = (
         1e6
         * wave_cost_factor

scripts/retrieve_databundle.py

@@ -58,9 +58,8 @@ if __name__ == "__main__":
     else:
         url = "https://zenodo.org/record/3517935/files/pypsa-eur-data-bundle.tar.xz"
 
-    # Save locations
     tarball_fn = Path(f"{rootpath}/bundle.tar.xz")
-    to_fn = Path(f"{rootpath}/data")
+    to_fn = Path(rootpath) / Path(snakemake.output[0]).parent.parent
 
     logger.info(f"Downloading databundle from '{url}'.")
     disable_progress = snakemake.config["run"].get("disable_progressbar", False)

scripts/retrieve_gas_infrastructure_data.py

@@ -29,7 +29,7 @@ if __name__ == "__main__":
 
     # Save locations
     zip_fn = Path(f"{rootpath}/IGGIELGN.zip")
-    to_fn = Path(f"{rootpath}/data/gas_network/scigrid-gas")
+    to_fn = Path(rootpath) / Path(snakemake.output[0]).parent.parent
 
     logger.info(f"Downloading databundle from '{url}'.")
     disable_progress = snakemake.config["run"].get("disable_progressbar", False)

scripts/retrieve_sector_databundle.py

@@ -10,23 +10,25 @@ import logging
 
 logger = logging.getLogger(__name__)
 
-import os
-import sys
 import tarfile
 from pathlib import Path
 
-# Add pypsa-eur scripts to path for import of _helpers
-sys.path.insert(0, os.getcwd() + "/../pypsa-eur/scripts")
-
 from _helpers import configure_logging, progress_retrieve
 
 if __name__ == "__main__":
+    if "snakemake" not in globals():
+        from _helpers import mock_snakemake
+
+        snakemake = mock_snakemake("retrieve_databundle")
+        rootpath = ".."
+    else:
+        rootpath = "."
     configure_logging(snakemake)
 
     url = "https://zenodo.org/record/5824485/files/pypsa-eur-sec-data-bundle.tar.gz"
 
-    tarball_fn = Path("sector-bundle.tar.gz")
+    tarball_fn = Path(f"{rootpath}/sector-bundle.tar.gz")
-    to_fn = Path("data")
+    to_fn = Path(rootpath) / Path(snakemake.output[0]).parent.parent
 
     logger.info(f"Downloading databundle from '{url}'.")
     disable_progress = snakemake.config["run"].get("disable_progressbar", False)

scripts/solve_network.py

@@ -38,7 +38,6 @@ from _helpers import (
     override_component_attrs,
     update_config_with_sector_opts,
 )
-from vresutils.benchmark import memory_logger
 
 logger = logging.getLogger(__name__)
 pypsa.pf.logger.setLevel(logging.WARNING)
@@ -601,6 +600,7 @@ def solve_network(n, config, opts="", **kwargs):
     track_iterations = cf_solving.get("track_iterations", False)
     min_iterations = cf_solving.get("min_iterations", 4)
     max_iterations = cf_solving.get("max_iterations", 6)
+    transmission_losses = cf_solving.get("transmission_losses", 0)
 
     # add to network for extra_functionality
     n.config = config
@@ -614,6 +614,7 @@ def solve_network(n, config, opts="", **kwargs):
     if skip_iterations:
         status, condition = n.optimize(
             solver_name=solver_name,
+            transmission_losses=transmission_losses,
             extra_functionality=extra_functionality,
             **solver_options,
             **kwargs,
@@ -624,6 +625,7 @@ def solve_network(n, config, opts="", **kwargs):
             track_iterations=track_iterations,
             min_iterations=min_iterations,
             max_iterations=max_iterations,
+            transmission_losses=transmission_losses,
             extra_functionality=extra_functionality,
             **solver_options,
             **kwargs,
@@ -667,23 +669,17 @@ if __name__ == "__main__":
 
     np.random.seed(solve_opts.get("seed", 123))
 
-    fn = getattr(snakemake.log, "memory", None)
+    if "overrides" in snakemake.input.keys():
-    with memory_logger(filename=fn, interval=30.0) as mem:
+        overrides = override_component_attrs(snakemake.input.overrides)
-        if "overrides" in snakemake.input.keys():
+        n = pypsa.Network(snakemake.input.network, override_component_attrs=overrides)
-            overrides = override_component_attrs(snakemake.input.overrides)
+    else:
-            n = pypsa.Network(
+        n = pypsa.Network(snakemake.input.network)
-                snakemake.input.network, override_component_attrs=overrides
-            )
-        else:
-            n = pypsa.Network(snakemake.input.network)
 
-        n = prepare_network(n, solve_opts, config=snakemake.config)
+    n = prepare_network(n, solve_opts, config=snakemake.config)
 
-        n = solve_network(
+    n = solve_network(
-            n, config=snakemake.config, opts=opts, log_fn=snakemake.log.solver
+        n, config=snakemake.config, opts=opts, log_fn=snakemake.log.solver
-        )
+    )
 
-        n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
+    n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
-        n.export_to_netcdf(snakemake.output[0])
+    n.export_to_netcdf(snakemake.output[0])
-
-    logger.info("Maximum memory usage: {}".format(mem.mem_usage))

scripts/solve_operations_network.py

@@ -17,7 +17,6 @@ from _helpers import (
     update_config_with_sector_opts,
 )
 from solve_network import prepare_network, solve_network
-from vresutils.benchmark import memory_logger
 
 logger = logging.getLogger(__name__)
 
@@ -46,23 +45,17 @@ if __name__ == "__main__":
 
     np.random.seed(solve_opts.get("seed", 123))
 
-    fn = getattr(snakemake.log, "memory", None)
+    if "overrides" in snakemake.input:
-    with memory_logger(filename=fn, interval=30.0) as mem:
+        overrides = override_component_attrs(snakemake.input.overrides)
-        if "overrides" in snakemake.input:
+        n = pypsa.Network(snakemake.input.network, override_component_attrs=overrides)
-            overrides = override_component_attrs(snakemake.input.overrides)
+    else:
-            n = pypsa.Network(
+        n = pypsa.Network(snakemake.input.network)
-                snakemake.input.network, override_component_attrs=overrides
-            )
-        else:
-            n = pypsa.Network(snakemake.input.network)
 
-        n.optimize.fix_optimal_capacities()
+    n.optimize.fix_optimal_capacities()
-        n = prepare_network(n, solve_opts, config=snakemake.config)
+    n = prepare_network(n, solve_opts, config=snakemake.config)
-        n = solve_network(
+    n = solve_network(
-            n, config=snakemake.config, opts=opts, log_fn=snakemake.log.solver
+        n, config=snakemake.config, opts=opts, log_fn=snakemake.log.solver
-        )
+    )
 
-        n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
+    n.meta = dict(snakemake.config, **dict(wildcards=dict(snakemake.wildcards)))
-        n.export_to_netcdf(snakemake.output[0])
+    n.export_to_netcdf(snakemake.output[0])
-
-    logger.info("Maximum memory usage: {}".format(mem.mem_usage))