Note
You can download this example as a Jupyter notebook or start it in interactive mode.
Replace StorageUnits with fundamental Links and Stores#
This notebook demonstrates how storage units can be replaced by more fundamental components, and how their parameters map to each other.
[1]:
import pandas as pd
from numpy.testing import assert_almost_equal, assert_array_almost_equal
import pypsa
We define two functions we use in the following.
[2]:
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 = f"{su_to_replace} inflow"
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,
bus0=su["bus"],
bus1=bus_name,
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.add_constraints(
model["Store-e_nom"][store_name]
- model["Link-p_nom"][link_1_name] * ratio1
== 0,
name="store_fix_1",
)
model.add_constraints(
model["Store-e_nom"][store_name]
- model["Link-p_nom"][link_2_name] * ratio2
== 0,
name="store_fix_2",
)
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.optimize()
WARNING:pypsa.io:Importing network from PyPSA version v0.17.1 while current version is v0.31.0. Read the release notes at https://pypsa.readthedocs.io/en/latest/release_notes.html to prepare your network for import.
INFO:pypsa.io:Imported network storage-hvdc.nc has buses, carriers, generators, global_constraints, lines, links, loads, storage_units
WARNING:pypsa.consistency:The following links have carriers which are not defined:
Index(['TL 0', 'TL 1'], dtype='object', name='Link')
WARNING:pypsa.consistency:The following lines have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following lines have zero r, which could break the linear load flow:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following buses have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Bus')
WARNING:pypsa.consistency:The following links have carriers which are not defined:
Index(['TL 0', 'TL 1'], dtype='object', name='Link')
WARNING:pypsa.consistency:The following lines have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following lines have zero r, which could break the linear load flow:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following buses have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Bus')
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/envs/v0.31.0/lib/python3.12/site-packages/linopy/common.py:147: UserWarning: coords for dimension(s) ['Generator'] is not aligned with the pandas object. Previously, the indexes of the pandas were ignored and overwritten in these cases. Now, the pandas object's coordinates are taken considered for alignment.
warn(
INFO:linopy.model: Solve problem using Highs solver
INFO:linopy.io: Writing time: 0.19s
INFO:linopy.solvers:Log file at /tmp/highs.log
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 490 primals, 1116 duals
Objective: 1.47e+07
Solver model: available
Solver message: optimal
Running HiGHS 1.7.2 (git hash: 184e327): Copyright (c) 2024 HiGHS under MIT licence terms
Coefficient ranges:
Matrix [8e-05, 6e+00]
Cost [3e-03, 3e+03]
Bound [2e+01, 1e+06]
RHS [2e+00, 1e+03]
Presolving model
673 rows, 487 cols, 1922 nonzeros 0s
585 rows, 401 cols, 2258 nonzeros 0s
Presolve : Reductions: rows 585(-531); columns 401(-89); elements 2258(-245)
Solving the presolved LP
Using EKK dual simplex solver - serial
Iteration Objective Infeasibilities num(sum)
0 -1.4501656367e+05 Ph1: 385(1.3432e+06); Du: 125(265.685) 0s
413 1.4670508825e+07 Pr: 0(0); Du: 0(6.71685e-15) 0s
Solving the original LP from the solution after postsolve
Model status : Optimal
Simplex iterations: 413
Objective value : 1.4670508825e+07
HiGHS run time : 0.01
Writing the solution to /tmp/linopy-solve-iwa1u5yf.sol
INFO:pypsa.optimization.optimize:The shadow-prices of the constraints Generator-fix-p-lower, Generator-fix-p-upper, Generator-ext-p-lower, Generator-ext-p-upper, Line-fix-s-lower, Line-fix-s-upper, Line-ext-s-lower, Line-ext-s-upper, Link-fix-p-lower, Link-fix-p-upper, Link-ext-p-lower, Link-ext-p-upper, StorageUnit-fix-p_dispatch-lower, StorageUnit-fix-p_dispatch-upper, StorageUnit-ext-p_dispatch-lower, StorageUnit-ext-p_dispatch-upper, StorageUnit-fix-p_store-lower, StorageUnit-fix-p_store-upper, StorageUnit-ext-p_store-lower, StorageUnit-ext-p_store-upper, StorageUnit-fix-state_of_charge-lower, StorageUnit-fix-state_of_charge-upper, StorageUnit-ext-state_of_charge-lower, StorageUnit-ext-state_of_charge-upper, StorageUnit-state_of_charge_set, Kirchhoff-Voltage-Law, StorageUnit-energy_balance were not assigned to the network.
[3]:
('ok', 'optimal')
[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.optimize(extra_functionality=extra_functionality)
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"],
)
WARNING:pypsa.io:Importing network from PyPSA version v0.17.1 while current version is v0.31.0. Read the release notes at https://pypsa.readthedocs.io/en/latest/release_notes.html to prepare your network for import.
INFO:pypsa.io:Imported network storage-hvdc.nc has buses, carriers, generators, global_constraints, lines, links, loads, storage_units
WARNING:pypsa.consistency:The following links have carriers which are not defined:
Index(['TL 0', 'TL 1', 'Storage 0 converter AC to battery'], dtype='object', name='Link')
WARNING:pypsa.consistency:The following lines have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following lines have zero r, which could break the linear load flow:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following buses have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Bus')
WARNING:pypsa.consistency:The following links have carriers which are not defined:
Index(['TL 0', 'TL 1', 'Storage 0 converter AC to battery'], dtype='object', name='Link')
WARNING:pypsa.consistency:The following lines have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following lines have zero r, which could break the linear load flow:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following buses have carriers which are not defined:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Bus')
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/envs/v0.31.0/lib/python3.12/site-packages/linopy/common.py:147: UserWarning: coords for dimension(s) ['Generator'] is not aligned with the pandas object. Previously, the indexes of the pandas were ignored and overwritten in these cases. Now, the pandas object's coordinates are taken considered for alignment.
warn(
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/envs/v0.31.0/lib/python3.12/site-packages/linopy/variables.py:192: FutureWarning: Accessing a single value with `Variable[...]` and return type ScalarVariable is deprecated. In future, this will return a Variable.To get a ScalarVariable use `Variable.at[...]` instead.
warn(
INFO:linopy.model: Solve problem using Highs solver
INFO:linopy.io: Writing time: 0.21s
INFO:linopy.solvers:Log file at /tmp/highs.log
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 504 primals, 1157 duals
Objective: 1.47e+07
Solver model: available
Solver message: optimal
Running HiGHS 1.7.2 (git hash: 184e327): Copyright (c) 2024 HiGHS under MIT licence terms
Coefficient ranges:
Matrix [8e-05, 6e+00]
Cost [3e-03, 3e+03]
Bound [1e+06, 1e+06]
RHS [2e+00, 1e+03]
Presolving model
682 rows, 498 cols, 1942 nonzeros 0s
585 rows, 401 cols, 2254 nonzeros 0s
585 rows, 401 cols, 2254 nonzeros 0s
Presolve : Reductions: rows 585(-572); columns 401(-103); elements 2254(-315)
Solving the presolved LP
Using EKK dual simplex solver - serial
Iteration Objective Infeasibilities num(sum)
0 -1.4523089525e+05 Ph1: 390(1.38972e+06); Du: 128(270.398) 0s
426 1.4670508825e+07 Pr: 0(0); Du: 0(1.20459e-14) 0s
Solving the original LP from the solution after postsolve
Model status : Optimal
Simplex iterations: 426
Objective value : 1.4670508825e+07
HiGHS run time : 0.01
Writing the solution to /tmp/linopy-solve-852oioxj.sol
INFO:pypsa.optimization.optimize:The shadow-prices of the constraints Generator-fix-p-lower, Generator-fix-p-upper, Generator-ext-p-lower, Generator-ext-p-upper, Line-fix-s-lower, Line-fix-s-upper, Line-ext-s-lower, Line-ext-s-upper, Link-fix-p-lower, Link-fix-p-upper, Link-ext-p-lower, Link-ext-p-upper, Store-ext-e-lower, Store-ext-e-upper, StorageUnit-fix-p_dispatch-lower, StorageUnit-fix-p_dispatch-upper, StorageUnit-ext-p_dispatch-lower, StorageUnit-ext-p_dispatch-upper, StorageUnit-fix-p_store-lower, StorageUnit-fix-p_store-upper, StorageUnit-ext-p_store-lower, StorageUnit-ext-p_store-upper, StorageUnit-fix-state_of_charge-lower, StorageUnit-fix-state_of_charge-upper, StorageUnit-ext-state_of_charge-lower, StorageUnit-ext-state_of_charge-upper, StorageUnit-state_of_charge_set, Kirchhoff-Voltage-Law, StorageUnit-energy_balance, Store-energy_balance, store_fix_1, store_fix_2 were not assigned to the network.