Security-Constrained LOPF¶
In this example, the dispatch of generators is optimised using the security-constrained linear optimal power flow (SCLOPF) functionality of PyPSA, to guarantee that no branches are overloaded in the event of certain branch outages.
In [2]:
Copied!
import pypsa
n = pypsa.examples.scigrid_de()
# correct some infeasibilties in the network
for line_name in ["316", "527", "602"]:
n.lines.loc[line_name, "s_nom"] = 1200
now = n.snapshots[0]
n.plot(bus_size=0);
import pypsa
n = pypsa.examples.scigrid_de()
# correct some infeasibilties in the network
for line_name in ["316", "527", "602"]:
n.lines.loc[line_name, "s_nom"] = 1200
now = n.snapshots[0]
n.plot(bus_size=0);
/home/runner/work/PyPSA/PyPSA/pypsa/network/io.py:2082: FutureWarning: pandas infers the `str` dtype for string data since its version 3.0. PyPSA still converts it back to numpy object dtype on import, but will keep it from PyPSA 2.0 on. Set `pypsa.options.api.legacy_string_dtype` explicitly to suppress this warning. new_static = _coerce_string_dtypes(new_static)
INFO:pypsa.network.io:Imported network 'SciGrid-DE' has buses, carriers, generators, lines, loads, storage_units, transformers
First, let's run the network without any $N-1$ security constraints and see how much it costs to operate the system.:
In [3]:
Copied!
n0 = n.copy()
n0.optimize(snapshots=now)
n0.statistics.opex().sum()
n0 = n.copy()
n0.optimize(snapshots=now)
n0.statistics.opex().sum()
/tmp/ipykernel_7805/344066763.py:2: FutureWarning: The default value of `include_objective_constant` will change from True to False in version 2.0. Set `include_objective_constant` explicitly to suppress this warning. Using False improves LP numerical conditioning by not including the objective constant as a variable.
n0.optimize(snapshots=now)
WARNING:pypsa.consistency:The following transformers have zero r, which could break the linear load flow:
Index(['2', '5', '10', '12', '13', '15', '18', '20', '22', '24', '26', '30',
'32', '37', '42', '46', '52', '56', '61', '68', '69', '74', '78', '86',
'87', '94', '95', '96', '99', '100', '104', '105', '106', '107', '117',
'120', '123', '124', '125', '128', '129', '138', '143', '156', '157',
'159', '160', '165', '184', '191', '195', '201', '220', '231', '232',
'233', '236', '247', '248', '250', '251', '252', '261', '263', '264',
'267', '272', '279', '281', '282', '292', '303', '307', '308', '312',
'315', '317', '322', '332', '334', '336', '338', '351', '353', '360',
'362', '382', '384', '385', '391', '403', '404', '413', '421', '450',
'458'],
dtype='object', name='name')
INFO:linopy.model: Solve problem using Gurobi solver
INFO:linopy.model:Solver options: - log_to_console: False
INFO:linopy.io: Writing time: 0.08s
Set parameter WLSAccessID
Set parameter WLSSecret
Set parameter LicenseID to value 2537914
Academic license 2537914 - for non-commercial use only - registered to l.___@tu-berlin.de
Read LP format model from file /tmp/linopy-problem-uhah44f5.lp
Reading time = 0.01 seconds
obj: 5957 rows, 2485 columns, 10948 nonzeros
Set parameter LogToConsole to value 0
INFO:linopy.constants: Optimization successful: Status: ok Termination condition: optimal Solution: 2485 primals, 5957 duals Objective: 3.32e+05 Solver: gurobi Runtime: 0.03s Dual bound: 3.32e+05 Solver model: available Solver message: 2
INFO:pypsa.optimization.optimize:The shadow-prices of the constraints Generator-fix-p-lower, Generator-fix-p-upper, Line-fix-s-lower, Line-fix-s-upper, Transformer-fix-s-lower, Transformer-fix-s-upper, StorageUnit-fix-p_dispatch-lower, StorageUnit-fix-p_dispatch-upper, StorageUnit-fix-p_store-lower, StorageUnit-fix-p_store-upper, StorageUnit-fix-state_of_charge-lower, StorageUnit-fix-state_of_charge-upper, Kirchhoff-Voltage-Law, StorageUnit-energy_balance were not assigned to the network.
Out[3]:
np.float64(332383.51399999997)
The security-constrained linear optimal power flow (SCLOPF) is executed using the n.optimize.optimize_security_constrained() method.
This method takes a list of snapshots (here: now) and a list of branches that are to be considered as outages (here, the 30 lines with the highest loading).
In [4]:
Copied!
branch_outages = (n0.lines_t.p0.loc[now] / n0.lines.s_nom).nlargest(30).index
n.optimize.optimize_security_constrained(
now, branch_outages=branch_outages, log_to_console=False
)
n.statistics.opex().sum()
branch_outages = (n0.lines_t.p0.loc[now] / n0.lines.s_nom).nlargest(30).index
n.optimize.optimize_security_constrained(
now, branch_outages=branch_outages, log_to_console=False
)
n.statistics.opex().sum()
WARNING:pypsa.consistency:The following transformers have zero r, which could break the linear load flow:
Index(['2', '5', '10', '12', '13', '15', '18', '20', '22', '24', '26', '30',
'32', '37', '42', '46', '52', '56', '61', '68', '69', '74', '78', '86',
'87', '94', '95', '96', '99', '100', '104', '105', '106', '107', '117',
'120', '123', '124', '125', '128', '129', '138', '143', '156', '157',
'159', '160', '165', '184', '191', '195', '201', '220', '231', '232',
'233', '236', '247', '248', '250', '251', '252', '261', '263', '264',
'267', '272', '279', '281', '282', '292', '303', '307', '308', '312',
'315', '317', '322', '332', '334', '336', '338', '351', '353', '360',
'362', '382', '384', '385', '391', '403', '404', '413', '421', '450',
'458'],
dtype='object', name='name')
/home/runner/work/PyPSA/PyPSA/pypsa/optimization/abstract.py:437: FutureWarning: The default value of `include_objective_constant` will change from True to False in version 2.0. Set `include_objective_constant` explicitly to suppress this warning. Using False improves LP numerical conditioning by not including the objective constant as a variable. m = n.optimize.create_model(
INFO:linopy.model: Solve problem using Gurobi solver
INFO:linopy.model:Solver options: - log_to_console: False
INFO:linopy.io: Writing time: 0.12s
Set parameter WLSAccessID
Set parameter WLSSecret
Set parameter LicenseID to value 2537914
Academic license 2537914 - for non-commercial use only - registered to l.___@tu-berlin.de
Read LP format model from file /tmp/linopy-problem-7h0fu67b.lp
Reading time = 0.07 seconds
obj: 62837 rows, 2485 columns, 102832 nonzeros
Set parameter LogToConsole to value 0
INFO:linopy.constants: Optimization successful: Status: ok Termination condition: optimal Solution: 2485 primals, 62837 duals Objective: 4.25e+05 Solver: gurobi Runtime: 0.26s Dual bound: 4.25e+05 Solver model: available Solver message: 2
INFO:pypsa.optimization.optimize:The shadow-prices of the constraints Generator-fix-p-lower, Generator-fix-p-upper, Line-fix-s-lower, Line-fix-s-upper, Transformer-fix-s-lower, Transformer-fix-s-upper, StorageUnit-fix-p_dispatch-lower, StorageUnit-fix-p_dispatch-upper, StorageUnit-fix-p_store-lower, StorageUnit-fix-p_store-upper, StorageUnit-fix-state_of_charge-lower, StorageUnit-fix-state_of_charge-upper, Kirchhoff-Voltage-Law, StorageUnit-energy_balance, Line-fix-s-lower-security-for-Line-outage-in-sub-network-0, Line-fix-s-upper-security-for-Line-outage-in-sub-network-0, Transformer-fix-s-lower-security-for-Line-outage-in-sub-network-0, Transformer-fix-s-upper-security-for-Line-outage-in-sub-network-0 were not assigned to the network.
Out[4]:
np.float64(424597.83679)
You can see that the cost of operating the system in the given hour rises from 332 k€ to 427 k€ when the security constraints are applied.
The maps below indicate the difference in the dispatch patterns, line flows and marginal prices between the $N-0$ (first) and $N-1$ (second) cases.
In [5]:
Copied!
def plot_network(n, snapshot):
bus_size = (
n.statistics.supply(groupby="bus", components=["Generator", "StorageUnit"])
.groupby("bus")
.sum()
)
line_flows = n.lines_t.p0.loc[snapshot]
bus_color = n.buses_t.marginal_price.loc[snapshot]
line_loading = n.lines_t.p0.abs().loc[snapshot] / n.lines.s_nom
n.plot(
bus_size=bus_size / 30000,
bus_color=bus_color,
bus_cmap="Reds",
line_color=line_loading,
line_flow=line_flows / 50,
)
def plot_network(n, snapshot):
bus_size = (
n.statistics.supply(groupby="bus", components=["Generator", "StorageUnit"])
.groupby("bus")
.sum()
)
line_flows = n.lines_t.p0.loc[snapshot]
bus_color = n.buses_t.marginal_price.loc[snapshot]
line_loading = n.lines_t.p0.abs().loc[snapshot] / n.lines.s_nom
n.plot(
bus_size=bus_size / 30000,
bus_color=bus_color,
bus_cmap="Reds",
line_color=line_loading,
line_flow=line_flows / 50,
)
In [6]:
Copied!
plot_network(n0, now)
plot_network(n0, now)
In [7]:
Copied!
plot_network(n, now)
plot_network(n, now)
We can also look at where the nodal dispatch is ramped up (red) or down (blue) in the $N-1$ case compared to the $N-0$ case. Mostly this means ramp down upstream of the potential outages, and ramp up downstream of the potential outages.
In [8]:
Copied!
bus_size0 = (
(
n0.statistics.supply(groupby="bus", components=["Generator", "StorageUnit"])
.groupby("bus")
.sum()
)
.reindex(index=n.buses.index)
.fillna(0)
)
bus_size1 = (
(
n.statistics.supply(groupby="bus", components=["Generator", "StorageUnit"])
.groupby("bus")
.sum()
)
.reindex(index=n.buses.index)
.fillna(0)
)
bus_size = bus_size1 - bus_size0
n.plot(
bus_size=bus_size.abs() / 30000,
bus_color=bus_size.map(lambda x: "b" if x < 0 else "r"),
);
bus_size0 = (
(
n0.statistics.supply(groupby="bus", components=["Generator", "StorageUnit"])
.groupby("bus")
.sum()
)
.reindex(index=n.buses.index)
.fillna(0)
)
bus_size1 = (
(
n.statistics.supply(groupby="bus", components=["Generator", "StorageUnit"])
.groupby("bus")
.sum()
)
.reindex(index=n.buses.index)
.fillna(0)
)
bus_size = bus_size1 - bus_size0
n.plot(
bus_size=bus_size.abs() / 30000,
bus_color=bus_size.map(lambda x: "b" if x < 0 else "r"),
);
We can also double-check that the $N-1$ constraints are satisfied and no lines are overloaded in the event of the outages considered:
In [9]:
Copied!
n.optimize.fix_optimal_dispatch()
p0_test = n.lpf_contingency(now, branch_outages=branch_outages)
n.optimize.fix_optimal_dispatch()
p0_test = n.lpf_contingency(now, branch_outages=branch_outages)
INFO:pypsa.network.power_flow:Performing linear load-flow on AC sub-network <pypsa.SubNetwork object at 0x7f0b8e73ad70> for snapshot(s) DatetimeIndex(['2011-01-01'], dtype='datetime64[ns]', name='snapshot', freq=None)
WARNING:pypsa.network.power_flow:No type given for 183, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 214, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 350, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 389, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 448, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 586, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 707, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 809, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 818, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 390, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 587, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 124, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 608, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 682, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 683, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 29, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 218, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 714, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 595, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 596, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 371, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 17, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 514, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 511, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 510, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 58, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 254, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 792, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 436, assuming it is a line
WARNING:pypsa.network.power_flow:No type given for 11, assuming it is a line
Check the maximum loading as per unit of s_nom in each contingency:
In [10]:
Copied!
max_loading = abs(p0_test.divide(n.passive_branches().s_nom, axis=0)).max()
max_loading
max_loading = abs(p0_test.divide(n.passive_branches().s_nom, axis=0)).max()
max_loading
Out[10]:
base 1.0 (Line, 183) 1.0 (Line, 214) 1.0 (Line, 350) 1.0 (Line, 389) 1.0 (Line, 448) 1.0 (Line, 586) 1.0 (Line, 707) 1.0 (Line, 809) 1.0 (Line, 818) 1.0 (Line, 390) 1.0 (Line, 587) 1.0 (Line, 124) 0.0 (Line, 608) 1.0 (Line, 682) 1.0 (Line, 683) 1.0 (Line, 29) 1.0 (Line, 218) 1.0 (Line, 714) 1.0 (Line, 595) 1.0 (Line, 596) 1.0 (Line, 371) 1.0 (Line, 17) 1.0 (Line, 514) 1.0 (Line, 511) 1.0 (Line, 510) 1.0 (Line, 58) 1.0 (Line, 254) 1.0 (Line, 792) 1.0 (Line, 436) 1.0 (Line, 11) 1.0 dtype: float64
In [ ]:
Copied!