Note
You can download this example as a Jupyter notebook or start it in interactive mode.
Replace components with fundamental Links and Stores
This notebook demonstrates how generators and storage units can be replaced by more fundamental components, and how their parameters map to each other.
[1]:
import pypsa, os
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from numpy.testing import assert_almost_equal, assert_array_almost_equal
from pyomo.environ import Constraint
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
Cell In [1], line 1
----> 1 import pypsa, os
2 import numpy as np
3 import pandas as pd
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/pypsa/__init__.py:10
1 # -*- coding: utf-8 -*-
4 """
5 Python for Power Systems Analysis (PyPSA)
6
7 Grid calculation library.
8 """
---> 10 from pypsa import (
11 components,
12 contingency,
13 descriptors,
14 examples,
15 geo,
16 io,
17 linopf,
18 linopt,
19 networkclustering,
20 opf,
21 opt,
22 optimization,
23 pf,
24 plot,
25 )
26 from pypsa.components import Network, SubNetwork
28 __version__ = "0.21.2"
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/pypsa/components.py:50
37 from pypsa.io import (
38 export_to_csv_folder,
39 export_to_hdf5,
(...)
47 import_series_from_dataframe,
48 )
49 from pypsa.opf import network_lopf, network_opf
---> 50 from pypsa.optimization.optimize import OptimizationAccessor
51 from pypsa.pf import (
52 calculate_B_H,
53 calculate_dependent_values,
(...)
62 sub_network_pf,
63 )
64 from pypsa.plot import iplot, plot
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/pypsa/optimization/__init__.py:7
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 """
4 Build optimisation problems from PyPSA networks with Linopy.
5 """
----> 7 from pypsa.optimization import abstract, constraints, optimize, variables
8 from pypsa.optimization.optimize import create_model
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/pypsa/optimization/constraints.py:9
6 import logging
8 import pandas as pd
----> 9 from linopy.expressions import LinearExpression, merge
10 from numpy import arange, cumsum, inf, nan, roll
11 from scipy import sparse
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/linopy/__init__.py:9
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 """
4 Created on Wed Mar 10 11:03:06 2021.
5
6 @author: fabulous
7 """
----> 9 from linopy import model, remote
10 from linopy.expressions import merge
11 from linopy.io import read_netcdf
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/linopy/model.py:22
20 from linopy import solvers
21 from linopy.common import best_int, replace_by_map
---> 22 from linopy.constraints import (
23 AnonymousConstraint,
24 AnonymousScalarConstraint,
25 Constraints,
26 )
27 from linopy.eval import Expr
28 from linopy.expressions import LinearExpression, ScalarLinearExpression
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/linopy/constraints.py:21
18 from scipy.sparse import coo_matrix
19 from xarray import DataArray, Dataset
---> 21 from linopy import expressions, variables
22 from linopy.common import (
23 _merge_inplace,
24 has_assigned_model,
(...)
27 replace_by_map,
28 )
31 class Constraint(DataArray):
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/linopy/expressions.py:23
20 from xarray.core.dataarray import DataArrayCoordinates
21 from xarray.core.groupby import _maybe_reorder, peek_at
---> 23 from linopy import constraints, variables
24 from linopy.common import as_dataarray
27 def exprwrap(method, *default_args, **new_default_kwargs):
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/site-packages/linopy/variables.py:398
393 roll = varwrap(DataArray.roll)
395 rolling = varwrap(DataArray.rolling)
--> 398 @dataclass(repr=False)
399 class Variables:
400 """
401 A variables container used for storing multiple variable arrays.
402 """
404 labels: Dataset = Dataset()
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/dataclasses.py:1211, in dataclass.<locals>.wrap(cls)
1210 def wrap(cls):
-> 1211 return _process_class(cls, init, repr, eq, order, unsafe_hash,
1212 frozen, match_args, kw_only, slots,
1213 weakref_slot)
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/dataclasses.py:959, in _process_class(cls, init, repr, eq, order, unsafe_hash, frozen, match_args, kw_only, slots, weakref_slot)
956 kw_only = True
957 else:
958 # Otherwise it's a field of some type.
--> 959 cls_fields.append(_get_field(cls, name, type, kw_only))
961 for f in cls_fields:
962 fields[f.name] = f
File ~/checkouts/readthedocs.org/user_builds/pypsa/conda/v0.21.2/lib/python3.11/dataclasses.py:816, in _get_field(cls, a_name, a_type, default_kw_only)
812 # For real fields, disallow mutable defaults. Use unhashable as a proxy
813 # indicator for mutability. Read the __hash__ attribute from the class,
814 # not the instance.
815 if f._field_type is _FIELD and f.default.__class__.__hash__ is None:
--> 816 raise ValueError(f'mutable default {type(f.default)} for field '
817 f'{f.name} is not allowed: use default_factory')
819 return f
ValueError: mutable default <class 'xarray.core.dataset.Dataset'> for field labels is not allowed: use default_factory
We define two functions we use in the following.
[2]:
def replace_gen(network, gen_to_replace):
"""Replace the generator gen_to_replace with a bus for the energy
carrier, a link for the conversion from the energy carrier to electricity
and a store to keep track of the depletion of the energy carrier and its
CO2 emissions."""
gen = network.generators.loc[gen_to_replace]
bus_name = "{} {}".format(gen["bus"], gen["carrier"])
link_name = "{} converter {} to AC".format(gen_to_replace, gen["carrier"])
store_name = "{} store {}".format(gen_to_replace, gen["carrier"])
network.add("Bus", bus_name, carrier=gen["carrier"])
network.add(
"Link",
link_name,
bus0=bus_name,
bus1=gen["bus"],
capital_cost=gen["capital_cost"] * gen["efficiency"],
p_nom=gen["p_nom"] / gen["efficiency"],
p_nom_extendable=gen["p_nom_extendable"],
p_nom_max=gen["p_nom_max"] / gen["efficiency"],
p_nom_min=gen["p_nom_min"] / gen["efficiency"],
p_max_pu=network.generators_t.p_max_pu.loc[:, gen_to_replace]
if gen_to_replace in network.generators_t.p_max_pu.columns
else gen["p_max_pu"],
p_min_pu=network.generators_t.p_min_pu.loc[:, gen_to_replace]
if gen_to_replace in network.generators_t.p_min_pu.columns
else gen["p_min_pu"],
marginal_cost=gen["marginal_cost"] * gen["efficiency"],
efficiency=gen["efficiency"],
)
network.add(
"Store",
store_name,
bus=bus_name,
e_nom_min=-float("inf"),
e_nom_max=0,
e_nom_extendable=True,
e_min_pu=1.0,
e_max_pu=0.0,
)
network.remove("Generator", gen_to_replace)
return bus_name, link_name, store_name
def replace_su(network, su_to_replace):
"""Replace the storage unit su_to_replace with a bus for the energy
carrier, two links for the conversion of the energy carrier to and from electricity,
a store to keep track of the depletion of the energy carrier and its
CO2 emissions, and a variable generator for the storage inflow.
Because the energy size and power size are linked in the storage unit by the max_hours,
extra functionality must be added to the LOPF to implement this constraint."""
su = network.storage_units.loc[su_to_replace]
bus_name = "{} {}".format(su["bus"], su["carrier"])
link_1_name = "{} converter {} to AC".format(su_to_replace, su["carrier"])
link_2_name = "{} converter AC to {}".format(su_to_replace, su["carrier"])
store_name = "{} store {}".format(su_to_replace, su["carrier"])
gen_name = "{} inflow".format(su_to_replace)
network.add("Bus", bus_name, carrier=su["carrier"])
# dispatch link
network.add(
"Link",
link_1_name,
bus0=bus_name,
bus1=su["bus"],
capital_cost=su["capital_cost"] * su["efficiency_dispatch"],
p_nom=su["p_nom"] / su["efficiency_dispatch"],
p_nom_extendable=su["p_nom_extendable"],
p_nom_max=su["p_nom_max"] / su["efficiency_dispatch"],
p_nom_min=su["p_nom_min"] / su["efficiency_dispatch"],
p_max_pu=su["p_max_pu"],
marginal_cost=su["marginal_cost"] * su["efficiency_dispatch"],
efficiency=su["efficiency_dispatch"],
)
# store link
network.add(
"Link",
link_2_name,
bus1=bus_name,
bus0=su["bus"],
p_nom=su["p_nom"],
p_nom_extendable=su["p_nom_extendable"],
p_nom_max=su["p_nom_max"],
p_nom_min=su["p_nom_min"],
p_max_pu=-su["p_min_pu"],
efficiency=su["efficiency_store"],
)
if (
su_to_replace in network.storage_units_t.state_of_charge_set.columns
and (
~pd.isnull(network.storage_units_t.state_of_charge_set[su_to_replace])
).any()
):
e_max_pu = pd.Series(data=1.0, index=network.snapshots)
e_min_pu = pd.Series(data=0.0, index=network.snapshots)
non_null = ~pd.isnull(
network.storage_units_t.state_of_charge_set[su_to_replace]
)
e_max_pu[non_null] = network.storage_units_t.state_of_charge_set[su_to_replace][
non_null
]
e_min_pu[non_null] = network.storage_units_t.state_of_charge_set[su_to_replace][
non_null
]
else:
e_max_pu = 1.0
e_min_pu = 0.0
network.add(
"Store",
store_name,
bus=bus_name,
e_nom=su["p_nom"] * su["max_hours"],
e_nom_min=su["p_nom_min"] / su["efficiency_dispatch"] * su["max_hours"],
e_nom_max=su["p_nom_max"] / su["efficiency_dispatch"] * su["max_hours"],
e_nom_extendable=su["p_nom_extendable"],
e_max_pu=e_max_pu,
e_min_pu=e_min_pu,
standing_loss=su["standing_loss"],
e_cyclic=su["cyclic_state_of_charge"],
e_initial=su["state_of_charge_initial"],
)
network.add("Carrier", "rain", co2_emissions=0.0)
# inflow from a variable generator, which can be curtailed (i.e. spilled)
inflow_max = network.storage_units_t.inflow[su_to_replace].max()
if inflow_max == 0.0:
inflow_pu = 0.0
else:
inflow_pu = network.storage_units_t.inflow[su_to_replace] / inflow_max
network.add(
"Generator",
gen_name,
bus=bus_name,
carrier="rain",
p_nom=inflow_max,
p_max_pu=inflow_pu,
)
if su["p_nom_extendable"]:
ratio2 = su["max_hours"]
ratio1 = ratio2 * su["efficiency_dispatch"]
def extra_functionality(network, snapshots):
model = network.model
model.store_fix_1 = Constraint(
rule=lambda model: model.store_e_nom[store_name]
== model.link_p_nom[link_1_name] * ratio1
)
model.store_fix_2 = Constraint(
rule=lambda model: model.store_e_nom[store_name]
== model.link_p_nom[link_2_name] * ratio2
)
else:
extra_functionality = None
network.remove("StorageUnit", su_to_replace)
return bus_name, link_1_name, link_2_name, store_name, gen_name, extra_functionality
Now, take an example from the git repo which has already been solved
[3]:
network_r = pypsa.examples.storage_hvdc(from_master=True)
network_r.lopf()
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In [3], line 1
----> 1 network_r = pypsa.examples.storage_hvdc(from_master=True)
2 network_r.lopf()
NameError: name 'pypsa' is not defined
[4]:
network = pypsa.examples.storage_hvdc(from_master=True)
su_to_replace = "Storage 0"
(
bus_name,
link_1_name,
link_2_name,
store_name,
gen_name,
extra_functionality,
) = replace_su(network, su_to_replace)
network.lopf(
network.snapshots, extra_functionality=extra_functionality, formulation="kirchhoff"
)
assert_almost_equal(network_r.objective, network.objective, decimal=2)
assert_array_almost_equal(
network_r.storage_units_t.state_of_charge[su_to_replace],
network.stores_t.e[store_name],
)
assert_array_almost_equal(
network_r.storage_units_t.p[su_to_replace],
-network.links_t.p1[link_1_name] - network.links_t.p0[link_2_name],
)
# check optimised size
assert_array_almost_equal(
network_r.storage_units.at[su_to_replace, "p_nom_opt"],
network.links.at[link_2_name, "p_nom_opt"],
)
assert_array_almost_equal(
network_r.storage_units.at[su_to_replace, "p_nom_opt"],
network.links.at[link_1_name, "p_nom_opt"]
* network_r.storage_units.at[su_to_replace, "efficiency_dispatch"],
)
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In [4], line 1
----> 1 network = pypsa.examples.storage_hvdc(from_master=True)
3 su_to_replace = "Storage 0"
5 (
6 bus_name,
7 link_1_name,
(...)
11 extra_functionality,
12 ) = replace_su(network, su_to_replace)
NameError: name 'pypsa' is not defined
[5]:
network = pypsa.examples.storage_hvdc(from_master=True)
gen_to_replace = "Gas 0"
bus_name, link_name, store_name = replace_gen(network, gen_to_replace)
network.lopf(network.snapshots)
assert_almost_equal(network_r.objective, network.objective, decimal=2)
# check dispatch
assert_array_almost_equal(
-network.links_t.p1[link_name], network_r.generators_t.p[gen_to_replace]
)
# check optimised size
assert_array_almost_equal(
network_r.generators.at[gen_to_replace, "p_nom_opt"],
network.links.at[link_name, "p_nom_opt"]
* network.links.at[link_name, "efficiency"],
)
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In [5], line 1
----> 1 network = pypsa.examples.storage_hvdc(from_master=True)
3 gen_to_replace = "Gas 0"
5 bus_name, link_name, store_name = replace_gen(network, gen_to_replace)
NameError: name 'pypsa' is not defined