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 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
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 = "{} 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,
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.25.1. 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.components: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.components:The following lines have zero r, which could break the linear load flow:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
INFO:linopy.model: Solve problem using Glpk solver
INFO:linopy.io: Writing time: 0.24s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 490 primals, 1104 duals
Objective: 1.47e+07
Solver model: not available
Solver message: optimal
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.
GLPSOL--GLPK LP/MIP Solver 5.0
Parameter(s) specified in the command line:
--lp /tmp/linopy-problem-zqcucctl.lp --output /tmp/linopy-solve-20ukpja1.sol
Reading problem data from '/tmp/linopy-problem-zqcucctl.lp'...
1104 rows, 490 columns, 2359 non-zeros
6379 lines were read
GLPK Simplex Optimizer 5.0
1104 rows, 490 columns, 2359 non-zeros
Preprocessing...
671 rows, 462 columns, 1895 non-zeros
Scaling...
A: min|aij| = 7.921e-05 max|aij| = 6.000e+00 ratio = 7.575e+04
GM: min|aij| = 1.287e-01 max|aij| = 7.771e+00 ratio = 6.038e+01
EQ: min|aij| = 1.656e-02 max|aij| = 1.000e+00 ratio = 6.038e+01
Constructing initial basis...
Size of triangular part is 671
0: obj = 4.758946148e+07 inf = 8.727e+06 (102)
280: obj = 3.145991071e+07 inf = 0.000e+00 (0) 2
* 536: obj = 1.467050883e+07 inf = 2.482e-12 (0) 2
OPTIMAL LP SOLUTION FOUND
Time used: 0.0 secs
Memory used: 1.3 Mb (1324551 bytes)
Writing basic solution to '/tmp/linopy-solve-20ukpja1.sol'...
[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.25.1. 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.components: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.components:The following lines have zero r, which could break the linear load flow:
Index(['0', '1', '2', '3', '4', '5'], dtype='object', name='Line')
INFO:linopy.model: Solve problem using Glpk solver
INFO:linopy.io: Writing time: 0.28s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 504 primals, 1144 duals
Objective: 1.47e+07
Solver model: not available
Solver message: optimal
GLPSOL--GLPK LP/MIP Solver 5.0
Parameter(s) specified in the command line:
--lp /tmp/linopy-problem-jaorgnst.lp --output /tmp/linopy-solve-flmqlikw.sol
Reading problem data from '/tmp/linopy-problem-jaorgnst.lp'...
1144 rows, 504 columns, 2413 non-zeros
6561 lines were read
GLPK Simplex Optimizer 5.0
1144 rows, 504 columns, 2413 non-zeros
Preprocessing...
685 rows, 476 columns, 1923 non-zeros
Scaling...
A: min|aij| = 7.921e-05 max|aij| = 6.000e+00 ratio = 7.575e+04
GM: min|aij| = 1.287e-01 max|aij| = 7.771e+00 ratio = 6.038e+01
EQ: min|aij| = 1.656e-02 max|aij| = 1.000e+00 ratio = 6.038e+01
Constructing initial basis...
Size of triangular part is 685
0: obj = 4.758946148e+07 inf = 8.736e+06 (102)
275: obj = 3.517083992e+07 inf = 3.839e-13 (0) 2
* 547: obj = 1.467050883e+07 inf = 1.808e-11 (0) 2
OPTIMAL LP SOLUTION FOUND
Time used: 0.0 secs
Memory used: 1.3 Mb (1360043 bytes)
Writing basic solution to '/tmp/linopy-solve-flmqlikw.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.