Note
You can download this example as a Jupyter notebook or start it in interactive mode.
Non-linear power flow after LOPF#
In this example, the dispatch of generators is optimised using the linear OPF, then a non-linear power flow is run on the resulting dispatch.
Data sources#
Grid: based on SciGRID Version 0.2 which is based on OpenStreetMap.
Load size and location: based on Landkreise (NUTS 3) GDP and population.
Load time series: from ENTSO-E hourly data, scaled up uniformly by factor 1.12 (a simplification of the methodology in Schumacher, Hirth (2015)).
Conventional power plant capacities and locations: BNetzA list.
Wind and solar capacities and locations: EEG Stammdaten, based on http://www.energymap.info/download.html, which represents capacities at the end of 2014. Units without PLZ are removed.
Wind and solar time series: REatlas, Andresen et al, “Validation of Danish wind time series from a new global renewable energy atlas for energy system analysis,” Energy 93 (2015) 1074 - 1088.
Where SciGRID nodes have been split into 220kV and 380kV substations, all load and generation is attached to the 220kV substation.
Warnings#
The data behind the notebook is no longer supported. See PyPSA/pypsa-eur for a newer model that covers the whole of Europe.
This dataset is ONLY intended to demonstrate the capabilities of PyPSA and is NOT (yet) accurate enough to be used for research purposes.
Known problems include:
Rough approximations have been made for missing grid data, e.g. 220kV-380kV transformers and connections between close sub-stations missing from OSM.
There appears to be some unexpected congestion in parts of the network, which may mean for example that the load attachment method (by Voronoi cell overlap with Landkreise) isn’t working, particularly in regions with a high density of substations.
Attaching power plants to the nearest high voltage substation may not reflect reality.
There is no proper n-1 security in the calculations - this can either be simulated with a blanket e.g. 70% reduction in thermal limits (as done here) or a proper security constrained OPF (see e.g. http://www.pypsa.org/examples/scigrid-sclopf.ipynb).
The borders and neighbouring countries are not represented.
Hydroelectric power stations are not modelled accurately.
The marginal costs are illustrative, not accurate.
Only the first day of 2011 is in the github dataset, which is not representative. The full year of 2011 can be downloaded at http://www.pypsa.org/examples/scigrid-with-load-gen-trafos-2011.zip.
The ENTSO-E total load for Germany may not be scaled correctly; it is scaled up uniformly by factor 1.12 (a simplification of the methodology in Schumacher, Hirth (2015), which suggests monthly factors).
Biomass from the EEG Stammdaten are not read in at the moment.
Power plant start up costs, ramping limits/costs, minimum loading rates are not considered.
[1]:
import pypsa
import numpy as np
import pandas as pd
import os
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
[2]:
network = pypsa.examples.scigrid_de(from_master=True)
WARNING:pypsa.io:Importing network from PyPSA version v0.17.1 while current version is v0.23.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 scigrid-de.nc has buses, generators, lines, loads, storage_units, transformers
Plot the distribution of the load and of generating tech
[3]:
fig, ax = plt.subplots(
1, 1, subplot_kw={"projection": ccrs.EqualEarth()}, figsize=(8, 8)
)
load_distribution = (
network.loads_t.p_set.loc[network.snapshots[0]].groupby(network.loads.bus).sum()
)
network.plot(bus_sizes=1e-5 * load_distribution, ax=ax, title="Load distribution");
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/conda/latest/lib/python3.11/site-packages/cartopy/mpl/style.py:76: UserWarning: facecolor will have no effect as it has been defined as "never".
warnings.warn('facecolor will have no effect as it has been '

[4]:
network.generators.groupby("carrier")["p_nom"].sum()
[4]:
carrier
Brown Coal 20879.500000
Gas 23913.130000
Geothermal 31.700000
Hard Coal 25312.600000
Multiple 152.700000
Nuclear 12068.000000
Oil 2710.200000
Other 3027.800000
Run of River 3999.100000
Solar 37041.524779
Storage Hydro 1445.000000
Waste 1645.900000
Wind Offshore 2973.500000
Wind Onshore 37339.895329
Name: p_nom, dtype: float64
[5]:
network.storage_units.groupby("carrier")["p_nom"].sum()
[5]:
carrier
Pumped Hydro 9179.5
Name: p_nom, dtype: float64
[6]:
techs = ["Gas", "Brown Coal", "Hard Coal", "Wind Offshore", "Wind Onshore", "Solar"]
n_graphs = len(techs)
n_cols = 3
if n_graphs % n_cols == 0:
n_rows = n_graphs // n_cols
else:
n_rows = n_graphs // n_cols + 1
fig, axes = plt.subplots(
nrows=n_rows, ncols=n_cols, subplot_kw={"projection": ccrs.EqualEarth()}
)
size = 6
fig.set_size_inches(size * n_cols, size * n_rows)
for i, tech in enumerate(techs):
i_row = i // n_cols
i_col = i % n_cols
ax = axes[i_row, i_col]
gens = network.generators[network.generators.carrier == tech]
gen_distribution = (
gens.groupby("bus").sum()["p_nom"].reindex(network.buses.index, fill_value=0.0)
)
network.plot(ax=ax, bus_sizes=2e-5 * gen_distribution)
ax.set_title(tech)
fig.tight_layout()
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/conda/latest/lib/python3.11/site-packages/cartopy/mpl/style.py:76: UserWarning: facecolor will have no effect as it has been defined as "never".
warnings.warn('facecolor will have no effect as it has been '

Run Linear Optimal Power Flow on the first day of 2011.
To approximate n-1 security and allow room for reactive power flows, don’t allow any line to be loaded above 70% of their thermal rating
[7]:
contingency_factor = 0.7
network.lines.s_max_pu = contingency_factor
There are some infeasibilities without small extensions
[8]:
network.lines.loc[["316", "527", "602"], "s_nom"] = 1715
We performing a linear OPF for one day, 4 snapshots at a time.
[9]:
group_size = 4
network.storage_units.state_of_charge_initial = 0.0
for i in range(int(24 / group_size)):
# set the initial state of charge based on previous round
if i:
network.storage_units.state_of_charge_initial = (
network.storage_units_t.state_of_charge.loc[
network.snapshots[group_size * i - 1]
]
)
network.optimize(
network.snapshots[group_size * i : group_size * i + group_size],
solver_name="cbc",
)
WARNING:pypsa.components: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='Transformer')
WARNING:pypsa.components: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='Transformer')
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
INFO:linopy.model: Solve linear problem using Cbc solver
INFO:linopy.io:Writing objective.
Writing constraints.: 100%|██████████| 15/15 [00:00<00:00, 64.21it/s]
Writing variables.: 100%|██████████| 7/7 [00:00<00:00, 231.66it/s]
INFO:linopy.io: Writing time: 0.3s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 9940 primals, 23828 duals
Objective: 1.45e+06
Solver model: not available
Solver message: Optimal - objective value 1449687.25014697
Welcome to the CBC MILP Solver
Version: 2.10.10
Build Date: Apr 19 2023
command line - cbc -printingOptions all -import /tmp/linopy-problem-m02q7lbe.lp -solve -solu /tmp/linopy-solve-6zjwgoop.sol (default strategy 1)
Option for printingOptions changed from normal to all
Presolve 2370 (-21458) rows, 4321 (-5619) columns and 15455 (-28591) elements
Perturbing problem by 0.001% of 2258.0586 - largest nonzero change 0.00099976141 ( 0.017785398%) - largest zero change 0.00099947237
0 Obj 88132.25 Primal inf 5731630.7 (2165) Dual inf 303.14607 (4)
122 Obj 18780.935 Primal inf 3741280.4 (2105)
244 Obj 18782.011 Primal inf 3156892.8 (2057)
366 Obj 18782.868 Primal inf 2974452.7 (2015)
488 Obj 18783.553 Primal inf 2896611.5 (1961)
610 Obj 18785.856 Primal inf 2633354.5 (1910)
732 Obj 18786.768 Primal inf 5489151.4 (1846)
854 Obj 18789.822 Primal inf 2603103.1 (1768)
976 Obj 18793.855 Primal inf 3580883.3 (1724)
1098 Obj 18799.182 Primal inf 2346268.8 (1609)
1220 Obj 18801.571 Primal inf 3652423.9 (1578)
1342 Obj 18804.272 Primal inf 7973281.5 (1484)
1464 Obj 18806.306 Primal inf 8234356.8 (1440)
1586 Obj 22625.93 Primal inf 13245722 (1401)
1708 Obj 22630.526 Primal inf 6509204.9 (1335)
1830 Obj 22635.378 Primal inf 4610184.4 (1164)
1952 Obj 22641.296 Primal inf 9984859.6 (1126)
2074 Obj 24886.471 Primal inf 8314771.3 (1002)
2196 Obj 86626.796 Primal inf 1464418.1 (833)
2318 Obj 170935.97 Primal inf 15776717 (854)
2440 Obj 317367.81 Primal inf 63700418 (723)
2562 Obj 474031.76 Primal inf 945173.56 (440)
2684 Obj 972991.86 Primal inf 992999 (378)
2806 Obj 1276987.6 Primal inf 22645.097 (190)
2928 Obj 1385632.5 Primal inf 13127.236 (119)
3050 Obj 1437639.7 Primal inf 949476.01 (628)
3172 Obj 1449323 Primal inf 116.40819 (21)
3213 Obj 1449696.3
3213 Obj 1449687.3 Dual inf 0.00094475941 (3)
3217 Obj 1449687.3
Optimal - objective value 1449687.3
After Postsolve, objective 1449687.3, infeasibilities - dual 214.20597 (7), primal 1.7945027e-06 (7)
Presolved model was optimal, full model needs cleaning up
0 Obj 1449687.3 Dual inf 2.142059 (7)
End of values pass after 10 iterations
10 Obj 1449687.3
Optimal - objective value 1449687.3
Optimal objective 1449687.25 - 3227 iterations time 0.412, Presolve 0.02
Total time (CPU seconds): 0.61 (Wallclock seconds): 0.57
WARNING:pypsa.components: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='Transformer')
WARNING:pypsa.components: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='Transformer')
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
INFO:linopy.model: Solve linear problem using Cbc solver
INFO:linopy.io:Writing objective.
Writing constraints.: 100%|██████████| 15/15 [00:00<00:00, 65.18it/s]
Writing variables.: 100%|██████████| 7/7 [00:00<00:00, 224.03it/s]
INFO:linopy.io: Writing time: 0.3s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 9940 primals, 23828 duals
Objective: 8.74e+05
Solver model: not available
Solver message: Optimal - objective value 873763.90232907
Welcome to the CBC MILP Solver
Version: 2.10.10
Build Date: Apr 19 2023
command line - cbc -printingOptions all -import /tmp/linopy-problem-cymqu1b4.lp -solve -solu /tmp/linopy-solve-ajf21uvo.sol (default strategy 1)
Option for printingOptions changed from normal to all
Presolve 2391 (-21437) rows, 4359 (-5581) columns and 15576 (-28470) elements
Perturbing problem by 0.001% of 2620.1345 - largest nonzero change 0.00099700892 ( 0.019493681%) - largest zero change 0.00099607519
0 Obj 81156.815 Primal inf 5995243.1 (2192) Dual inf 447.30519 (4)
122 Obj 11805.501 Primal inf 4003831 (2158)
244 Obj 11806.142 Primal inf 3894941.3 (2105)
366 Obj 11808.343 Primal inf 3489517.9 (2054)
488 Obj 11808.697 Primal inf 2879336.3 (1987)
610 Obj 11810.936 Primal inf 2805000.5 (1953)
732 Obj 11812.831 Primal inf 8250220.7 (1906)
854 Obj 11815.308 Primal inf 9179258.7 (1826)
976 Obj 11819.869 Primal inf 9834441 (1778)
1098 Obj 11823.145 Primal inf 10259354 (1723)
1220 Obj 11827.588 Primal inf 9027524.3 (1692)
1342 Obj 11831.177 Primal inf 5661830.6 (1640)
1464 Obj 11836.33 Primal inf 8649644.1 (1515)
1586 Obj 11840.748 Primal inf 25990687 (1446)
1708 Obj 11846.02 Primal inf 23035619 (1364)
1830 Obj 11852.455 Primal inf 7786896.5 (1193)
1952 Obj 11857.707 Primal inf 5105664.6 (1059)
2074 Obj 12212.649 Primal inf 5.0970035e+08 (1027)
2196 Obj 13492.323 Primal inf 1.1483789e+09 (1014)
2318 Obj 54457.313 Primal inf 48601189 (840)
2440 Obj 224493.74 Primal inf 77655088 (845)
2562 Obj 377861.99 Primal inf 13558376 (663)
2684 Obj 499790.01 Primal inf 222833.36 (421)
2806 Obj 663977.22 Primal inf 409323.95 (329)
2928 Obj 813446.17 Primal inf 36340.277 (165)
3050 Obj 857175.94 Primal inf 171789.94 (412)
3172 Obj 873625.51 Primal inf 94.924055 (21)
3209 Obj 873771.66
Optimal - objective value 873763.9
After Postsolve, objective 873763.9, infeasibilities - dual 78.034695 (3), primal 9.1550873e-07 (3)
Presolved model was optimal, full model needs cleaning up
0 Obj 873763.9 Dual inf 0.78034665 (3)
End of values pass after 3 iterations
3 Obj 873763.9
Optimal - objective value 873763.9
Optimal objective 873763.9023 - 3212 iterations time 0.372, Presolve 0.02
Total time (CPU seconds): 0.58 (Wallclock seconds): 0.53
WARNING:pypsa.components: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='Transformer')
WARNING:pypsa.components: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='Transformer')
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
INFO:linopy.model: Solve linear problem using Cbc solver
INFO:linopy.io:Writing objective.
Writing constraints.: 100%|██████████| 15/15 [00:00<00:00, 65.30it/s]
Writing variables.: 100%|██████████| 7/7 [00:00<00:00, 224.02it/s]
INFO:linopy.io: Writing time: 0.3s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 9940 primals, 23828 duals
Objective: 7.91e+05
Solver model: not available
Solver message: Optimal - objective value 790786.59488191
Welcome to the CBC MILP Solver
Version: 2.10.10
Build Date: Apr 19 2023
command line - cbc -printingOptions all -import /tmp/linopy-problem-flw5nwkg.lp -solve -solu /tmp/linopy-solve-dfqqgq80.sol (default strategy 1)
Option for printingOptions changed from normal to all
Presolve 2459 (-21369) rows, 4434 (-5506) columns and 15837 (-28209) elements
Perturbing problem by 0.001% of 2745.1675 - largest nonzero change 0.00099629266 ( 0.010269742%) - largest zero change 0.00099225248
0 Obj -46.556636 Primal inf 5896437.2 (2262)
124 Obj -46.556636 Primal inf 3846928.8 (2202)
248 Obj -44.82498 Primal inf 3891477.9 (2161)
372 Obj -44.156006 Primal inf 3309202.3 (2099)
496 Obj -43.381401 Primal inf 3305884.9 (2059)
620 Obj -41.881269 Primal inf 3984973.3 (2013)
744 Obj -40.263946 Primal inf 3632772.1 (1932)
868 Obj -38.4766 Primal inf 2692217.7 (1845)
992 Obj -36.318847 Primal inf 5247973.4 (1836)
1116 Obj -34.308108 Primal inf 3554862.8 (1797)
1240 Obj -29.841983 Primal inf 6654365.9 (1644)
1364 Obj -27.452571 Primal inf 10423780 (1593)
1488 Obj -23.449055 Primal inf 6977159.9 (1500)
1612 Obj -20.90136 Primal inf 4160718.4 (1331)
1736 Obj 7645.0605 Primal inf 25154037 (1372)
1860 Obj 7648.9785 Primal inf 33121681 (1293)
1984 Obj 7653.6105 Primal inf 36524583 (1193)
2108 Obj 7664.7036 Primal inf 3595244.1 (1058)
2232 Obj 7669.8619 Primal inf 7969548.7 (1105)
2356 Obj 18319.399 Primal inf 12583904 (868)
2480 Obj 156737.75 Primal inf 5.9686087e+08 (621)
2604 Obj 301824.11 Primal inf 7.4529458e+08 (539)
2728 Obj 465315.97 Primal inf 53702.188 (283)
2852 Obj 598987.02 Primal inf 148166.16 (392)
2976 Obj 704133.04 Primal inf 49089.848 (226)
3100 Obj 768065.66 Primal inf 251485.75 (227)
3224 Obj 785981.22 Primal inf 10111.895 (80)
3308 Obj 790794.89
Optimal - objective value 790786.59
After Postsolve, objective 790786.59, infeasibilities - dual 7.9999999 (1), primal 4.8500652e-07 (1)
Presolved model was optimal, full model needs cleaning up
0 Obj 790786.59 Dual inf 0.0799999 (1)
End of values pass after 1 iterations
1 Obj 790786.59
Optimal - objective value 790786.59
Optimal objective 790786.5949 - 3309 iterations time 0.412, Presolve 0.03
Total time (CPU seconds): 0.61 (Wallclock seconds): 0.57
WARNING:pypsa.components: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='Transformer')
WARNING:pypsa.components: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='Transformer')
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
INFO:linopy.model: Solve linear problem using Cbc solver
INFO:linopy.io:Writing objective.
Writing constraints.: 100%|██████████| 15/15 [00:00<00:00, 65.47it/s]
Writing variables.: 100%|██████████| 7/7 [00:00<00:00, 228.83it/s]
INFO:linopy.io: Writing time: 0.3s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 9940 primals, 23828 duals
Objective: 1.46e+06
Solver model: not available
Solver message: Optimal - objective value 1455232.67827293
Welcome to the CBC MILP Solver
Version: 2.10.10
Build Date: Apr 19 2023
command line - cbc -printingOptions all -import /tmp/linopy-problem-g7bybe5q.lp -solve -solu /tmp/linopy-solve-stn_524j.sol (default strategy 1)
Option for printingOptions changed from normal to all
Presolve 2446 (-21382) rows, 4424 (-5516) columns and 15796 (-28250) elements
Perturbing problem by 0.001% of 2193.4858 - largest nonzero change 0.00099891412 ( 0.012052428%) - largest zero change 0.00099818337
0 Obj 25296.904 Primal inf 5950415.8 (2263) Dual inf 118.67629 (1)
123 Obj 7959.0755 Primal inf 3978176.8 (2212)
246 Obj 7960.1705 Primal inf 3508620.9 (2157)
369 Obj 7961.9396 Primal inf 3123948.8 (2093)
492 Obj 7962.9138 Primal inf 2798936.5 (2048)
615 Obj 7964.9705 Primal inf 3357865.7 (1994)
738 Obj 7966.6686 Primal inf 3127552.2 (1948)
861 Obj 7971.3784 Primal inf 2301315.8 (1848)
984 Obj 7974.9175 Primal inf 2372425.7 (1783)
1107 Obj 7977.7474 Primal inf 1934878.7 (1700)
1230 Obj 7981.5501 Primal inf 3937741.9 (1699)
1353 Obj 7985.3938 Primal inf 18263103 (1584)
1476 Obj 7987.4254 Primal inf 3476458.5 (1569)
1599 Obj 12741.997 Primal inf 4981282.2 (1441)
1722 Obj 23727.135 Primal inf 2844103.7 (1342)
1845 Obj 23732.575 Primal inf 14172222 (1288)
1968 Obj 23737.962 Primal inf 5646282.6 (1395)
2091 Obj 23743.25 Primal inf 6416597.8 (1202)
2214 Obj 23869.961 Primal inf 5135444 (943)
2337 Obj 84853.49 Primal inf 724643 (705)
2460 Obj 329083.93 Primal inf 30185874 (846)
2583 Obj 434682.64 Primal inf 297958.23 (508)
2706 Obj 800514.73 Primal inf 405043.6 (494)
2829 Obj 1045136.2 Primal inf 423910.25 (319)
2952 Obj 1336163.6 Primal inf 20098.938 (181)
3075 Obj 1429673 Primal inf 4333.2787 (114)
3198 Obj 1455135.9 Primal inf 14.494187 (17)
3215 Obj 1455244.8
Optimal - objective value 1455232.7
After Postsolve, objective 1455232.7, infeasibilities - dual 25.49824 (1), primal 7.6227303e-08 (1)
Presolved model was optimal, full model needs cleaning up
0 Obj 1455232.7 Dual inf 0.30197075 (2)
End of values pass after 0 iterations
0 Obj 1455232.7 Dual inf 0.30197075 (2) w.o. free dual inf (1)
End of values pass after 12 iterations
12 Obj 1455232.7 Dual inf 0.046988453 (1) w.o. free dual inf (0)
13 Obj 1455232.7
Optimal - objective value 1455232.7
Optimal objective 1455232.678 - 3228 iterations time 0.402, Presolve 0.02
Total time (CPU seconds): 0.60 (Wallclock seconds): 0.57
WARNING:pypsa.components: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='Transformer')
WARNING:pypsa.components: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='Transformer')
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
INFO:linopy.model: Solve linear problem using Cbc solver
INFO:linopy.io:Writing objective.
Writing constraints.: 100%|██████████| 15/15 [00:00<00:00, 63.89it/s]
Writing variables.: 100%|██████████| 7/7 [00:00<00:00, 225.45it/s]
INFO:linopy.io: Writing time: 0.3s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 9940 primals, 23828 duals
Objective: 2.65e+06
Solver model: not available
Solver message: Optimal - objective value 2649026.44753696
Welcome to the CBC MILP Solver
Version: 2.10.10
Build Date: Apr 19 2023
command line - cbc -printingOptions all -import /tmp/linopy-problem-cbwtdjz4.lp -solve -solu /tmp/linopy-solve-ody6fmwr.sol (default strategy 1)
Option for printingOptions changed from normal to all
Presolve 2390 (-21438) rows, 4356 (-5584) columns and 15619 (-28427) elements
Perturbing problem by 0.001% of 2193.4858 - largest nonzero change 0.00099820297 ( 0.010152757%) - largest zero change 0.00099818337
0 Obj 97166.872 Primal inf 5972833.9 (2225) Dual inf 474.70514 (4)
122 Obj 27987.1 Primal inf 4466877.7 (2181)
244 Obj 27987.999 Primal inf 3523912.1 (2116)
366 Obj 27987.999 Primal inf 3307009.2 (2089)
488 Obj 27990.787 Primal inf 3429125.9 (2060)
610 Obj 27993.493 Primal inf 3698479 (1990)
732 Obj 27994.983 Primal inf 6970738.1 (1906)
854 Obj 28002.991 Primal inf 3701372.4 (1812)
976 Obj 28007.381 Primal inf 3743336.9 (1727)
1098 Obj 28009.916 Primal inf 11038730 (1714)
1220 Obj 28012.41 Primal inf 1919136.3 (1604)
1342 Obj 28015.155 Primal inf 2008027.9 (1498)
1464 Obj 44626.704 Primal inf 2004446.4 (1435)
1586 Obj 59481.392 Primal inf 1887845.6 (1288)
1708 Obj 59491.691 Primal inf 7485056.2 (1319)
1830 Obj 66967.925 Primal inf 7564881.8 (1248)
1952 Obj 87265.368 Primal inf 13515530 (1176)
2074 Obj 123053.26 Primal inf 28281720 (1264)
2196 Obj 194531.43 Primal inf 40238455 (1202)
2318 Obj 405537.74 Primal inf 7459420.9 (839)
2440 Obj 511907.66 Primal inf 9970654.1 (860)
2562 Obj 1194775.8 Primal inf 691561.84 (554)
2684 Obj 1884492.2 Primal inf 50732.517 (355)
2806 Obj 2231687.5 Primal inf 251590.65 (594)
2928 Obj 2462901.8 Primal inf 20933.126 (236)
3050 Obj 2592817.7 Primal inf 39739.072 (243)
3172 Obj 2623840.1 Primal inf 11215.152 (173)
3294 Obj 2645834.3 Primal inf 2552.0989 (97)
3382 Obj 2649041.9
Optimal - objective value 2649026.4
After Postsolve, objective 2649026.4, infeasibilities - dual 124.09584 (3), primal 7.3815671e-07 (3)
Presolved model was optimal, full model needs cleaning up
0 Obj 2649026.4 Dual inf 1.2409581 (3)
End of values pass after 12 iterations
12 Obj 2649026.4
Optimal - objective value 2649026.4
Optimal objective 2649026.448 - 3394 iterations time 0.452, Presolve 0.02
Total time (CPU seconds): 0.66 (Wallclock seconds): 0.62
WARNING:pypsa.components: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='Transformer')
WARNING:pypsa.components: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='Transformer')
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
<__array_function__ internals>:200: RuntimeWarning: invalid value encountered in cast
INFO:linopy.model: Solve linear problem using Cbc solver
INFO:linopy.io:Writing objective.
Writing constraints.: 100%|██████████| 15/15 [00:00<00:00, 65.11it/s]
Writing variables.: 100%|██████████| 7/7 [00:00<00:00, 227.35it/s]
INFO:linopy.io: Writing time: 0.3s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 9940 primals, 23828 duals
Objective: 2.14e+06
Solver model: not available
Solver message: Optimal - objective value 2142956.88872811
Welcome to the CBC MILP Solver
Version: 2.10.10
Build Date: Apr 19 2023
command line - cbc -printingOptions all -import /tmp/linopy-problem-uw6x393b.lp -solve -solu /tmp/linopy-solve-5aa4ajzr.sol (default strategy 1)
Option for printingOptions changed from normal to all
Presolve 2391 (-21437) rows, 4374 (-5566) columns and 15653 (-28393) elements
Perturbing problem by 0.001% of 8902.5915 - largest nonzero change 0.00097303945 ( 0.0039378659%) - largest zero change 0.00097207615
0 Obj 92710.129 Primal inf 6096896.4 (2218) Dual inf 1195.3238 (4)
122 Obj 23359.141 Primal inf 3993516.6 (2165)
244 Obj 23359.367 Primal inf 3279783 (2093)
366 Obj 23359.993 Primal inf 3144630.8 (2089)
488 Obj 23360.199 Primal inf 3082814.7 (2047)
610 Obj 23362.366 Primal inf 2788197.4 (1967)
732 Obj 23363.072 Primal inf 4921548.2 (1903)
854 Obj 23364.579 Primal inf 2683838.3 (1804)
976 Obj 23367.964 Primal inf 2691822.5 (1765)
1098 Obj 23368.91 Primal inf 4423665.5 (1640)
1220 Obj 23369.979 Primal inf 2721819.4 (1653)
1342 Obj 23371.108 Primal inf 2933038.8 (1571)
1464 Obj 23373.102 Primal inf 2246191.3 (1420)
1586 Obj 40824.39 Primal inf 1982682.2 (1277)
1708 Obj 40826.919 Primal inf 3724465.8 (1218)
1830 Obj 40829.913 Primal inf 2262983 (1192)
1952 Obj 40833.174 Primal inf 8634958.7 (1156)
2074 Obj 76839.224 Primal inf 5894395.7 (995)
2196 Obj 112207.39 Primal inf 1180179 (832)
2318 Obj 348102.06 Primal inf 841771.56 (686)
2440 Obj 678666.25 Primal inf 2318933.6 (749)
2562 Obj 949060.21 Primal inf 331056.5 (767)
2684 Obj 1436171 Primal inf 65058.383 (287)
2806 Obj 1845544.9 Primal inf 14568.194 (237)
2928 Obj 2045415.3 Primal inf 4098.1153 (157)
3050 Obj 2123896.8 Primal inf 38108.989 (361)
3172 Obj 2142760.5 Primal inf 7.4289097 (19)
3193 Obj 2142962.6
Optimal - objective value 2142956.9
After Postsolve, objective 2142956.9, infeasibilities - dual 138.1249 (4), primal 5.5034699e-07 (4)
Presolved model was optimal, full model needs cleaning up
0 Obj 2142956.9 Dual inf 1.3812486 (4)
End of values pass after 5 iterations
5 Obj 2142956.9
Optimal - objective value 2142956.9
Optimal objective 2142956.889 - 3198 iterations time 0.362, Presolve 0.02
Total time (CPU seconds): 0.54 (Wallclock seconds): 0.52
[10]:
p_by_carrier = network.generators_t.p.groupby(network.generators.carrier, axis=1).sum()
p_by_carrier.drop(
(p_by_carrier.max()[p_by_carrier.max() < 1700.0]).index, axis=1, inplace=True
)
p_by_carrier.columns
[10]:
Index(['Brown Coal', 'Gas', 'Hard Coal', 'Nuclear', 'Run of River', 'Solar',
'Wind Offshore', 'Wind Onshore'],
dtype='object', name='carrier')
[11]:
colors = {
"Brown Coal": "brown",
"Hard Coal": "k",
"Nuclear": "r",
"Run of River": "green",
"Wind Onshore": "blue",
"Solar": "yellow",
"Wind Offshore": "cyan",
"Waste": "orange",
"Gas": "orange",
}
# reorder
cols = [
"Nuclear",
"Run of River",
"Brown Coal",
"Hard Coal",
"Gas",
"Wind Offshore",
"Wind Onshore",
"Solar",
]
p_by_carrier = p_by_carrier[cols]
[12]:
c = [colors[col] for col in p_by_carrier.columns]
fig, ax = plt.subplots(figsize=(12, 6))
(p_by_carrier / 1e3).plot(kind="area", ax=ax, linewidth=4, color=c, alpha=0.7)
ax.legend(ncol=4, loc="upper left")
ax.set_ylabel("GW")
ax.set_xlabel("")
fig.tight_layout()

[13]:
fig, ax = plt.subplots(figsize=(12, 6))
p_storage = network.storage_units_t.p.sum(axis=1)
state_of_charge = network.storage_units_t.state_of_charge.sum(axis=1)
p_storage.plot(label="Pumped hydro dispatch", ax=ax, linewidth=3)
state_of_charge.plot(label="State of charge", ax=ax, linewidth=3)
ax.legend()
ax.grid()
ax.set_ylabel("MWh")
ax.set_xlabel("")
fig.tight_layout()

[14]:
now = network.snapshots[4]
With the linear load flow, there is the following per unit loading:
[15]:
loading = network.lines_t.p0.loc[now] / network.lines.s_nom
loading.describe()
[15]:
count 852.000000
mean -0.002877
std 0.258526
min -0.700000
25% -0.129689
50% 0.003084
75% 0.122643
max 0.700000
dtype: float64
[16]:
fig, ax = plt.subplots(subplot_kw={"projection": ccrs.EqualEarth()}, figsize=(9, 9))
network.plot(
ax=ax,
line_colors=abs(loading),
line_cmap=plt.cm.jet,
title="Line loading",
bus_sizes=1e-3,
bus_alpha=0.7,
)
fig.tight_layout();
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/conda/latest/lib/python3.11/site-packages/cartopy/mpl/style.py:76: UserWarning: facecolor will have no effect as it has been defined as "never".
warnings.warn('facecolor will have no effect as it has been '

Let’s have a look at the marginal prices
[17]:
network.buses_t.marginal_price.loc[now].describe()
[17]:
count 585.000000
mean 15.737598
std 10.941995
min -10.397824
25% 6.992120
50% 15.841190
75% 25.048186
max 52.150120
Name: 2011-01-01 04:00:00, dtype: float64
[18]:
fig, ax = plt.subplots(subplot_kw={"projection": ccrs.PlateCarree()}, figsize=(8, 8))
plt.hexbin(
network.buses.x,
network.buses.y,
gridsize=20,
C=network.buses_t.marginal_price.loc[now],
cmap=plt.cm.jet,
zorder=3,
)
network.plot(ax=ax, line_widths=pd.Series(0.5, network.lines.index), bus_sizes=0)
cb = plt.colorbar(location="bottom")
cb.set_label("Locational Marginal Price (EUR/MWh)")
fig.tight_layout()
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/conda/latest/lib/python3.11/site-packages/cartopy/mpl/style.py:76: UserWarning: facecolor will have no effect as it has been defined as "never".
warnings.warn('facecolor will have no effect as it has been '

Curtailment variable#
By considering how much power is available and how much is generated, you can see what share is curtailed:
[19]:
carrier = "Wind Onshore"
capacity = network.generators.groupby("carrier").sum().at[carrier, "p_nom"]
p_available = network.generators_t.p_max_pu.multiply(network.generators["p_nom"])
p_available_by_carrier = p_available.groupby(network.generators.carrier, axis=1).sum()
p_curtailed_by_carrier = p_available_by_carrier - p_by_carrier
[20]:
p_df = pd.DataFrame(
{
carrier + " available": p_available_by_carrier[carrier],
carrier + " dispatched": p_by_carrier[carrier],
carrier + " curtailed": p_curtailed_by_carrier[carrier],
}
)
p_df[carrier + " capacity"] = capacity
[21]:
p_df["Wind Onshore curtailed"][p_df["Wind Onshore curtailed"] < 0.0] = 0.0
[22]:
fig, ax = plt.subplots(figsize=(10, 4))
p_df[[carrier + " dispatched", carrier + " curtailed"]].plot(
kind="area", ax=ax, linewidth=3
)
p_df[[carrier + " available", carrier + " capacity"]].plot(ax=ax, linewidth=3)
ax.set_xlabel("")
ax.set_ylabel("Power [MW]")
ax.set_ylim([0, 40000])
ax.legend()
fig.tight_layout()

Non-Linear Power Flow#
Now perform a full Newton-Raphson power flow on the first hour. For the PF, set the P to the optimised P.
[23]:
network.generators_t.p_set = network.generators_t.p
network.storage_units_t.p_set = network.storage_units_t.p
Set all buses to PV, since we don’t know what Q set points are
[24]:
network.generators.control = "PV"
# Need some PQ buses so that Jacobian doesn't break
f = network.generators[network.generators.bus == "492"]
network.generators.loc[f.index, "control"] = "PQ"
Now, perform the non-linear PF.
[25]:
info = network.pf();
INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network SubNetwork 0 for snapshots DatetimeIndex(['2011-01-01 00:00:00', '2011-01-01 01:00:00',
'2011-01-01 02:00:00', '2011-01-01 03:00:00',
'2011-01-01 04:00:00', '2011-01-01 05:00:00',
'2011-01-01 06:00:00', '2011-01-01 07:00:00',
'2011-01-01 08:00:00', '2011-01-01 09:00:00',
'2011-01-01 10:00:00', '2011-01-01 11:00:00',
'2011-01-01 12:00:00', '2011-01-01 13:00:00',
'2011-01-01 14:00:00', '2011-01-01 15:00:00',
'2011-01-01 16:00:00', '2011-01-01 17:00:00',
'2011-01-01 18:00:00', '2011-01-01 19:00:00',
'2011-01-01 20:00:00', '2011-01-01 21:00:00',
'2011-01-01 22:00:00', '2011-01-01 23:00:00'],
dtype='datetime64[ns]', name='snapshot', freq=None)
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045870 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.047222 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.046564 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045386 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045401 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045926 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045623 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045584 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045581 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.046002 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044911 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044681 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045854 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044866 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045277 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045268 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044901 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045786 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045293 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.045301 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044935 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044666 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044946 seconds
INFO:pypsa.pf:Newton-Raphson solved in 4 iterations with error of 0.000000 in 0.044325 seconds
Any failed to converge?
[26]:
(~info.converged).any().any()
[26]:
False
With the non-linear load flow, there is the following per unit loading of the full thermal rating.
[27]:
(network.lines_t.p0.loc[now] / network.lines.s_nom).describe()
[27]:
count 852.000000
mean 0.000449
std 0.260873
min -0.813322
25% -0.125479
50% 0.003032
75% 0.126691
max 0.827140
dtype: float64
Let’s inspect the voltage angle differences across the lines have (in degrees)
[28]:
df = network.lines.copy()
for b in ["bus0", "bus1"]:
df = pd.merge(
df, network.buses_t.v_ang.loc[[now]].T, how="left", left_on=b, right_index=True
)
s = df[str(now) + "_x"] - df[str(now) + "_y"]
(s * 180 / np.pi).describe()
[28]:
count 852.000000
mean -0.022763
std 2.373578
min -12.158659
25% -0.463000
50% 0.001597
75% 0.537443
max 17.959258
dtype: float64
Plot the reactive power
[29]:
fig, ax = plt.subplots(subplot_kw={"projection": ccrs.EqualEarth()}, figsize=(9, 9))
q = network.buses_t.q.loc[now]
bus_colors = pd.Series("r", network.buses.index)
bus_colors[q < 0.0] = "b"
network.plot(
bus_sizes=1e-4 * abs(q),
ax=ax,
bus_colors=bus_colors,
title="Reactive power feed-in (red=+ve, blue=-ve)",
);
/home/docs/checkouts/readthedocs.org/user_builds/pypsa/conda/latest/lib/python3.11/site-packages/cartopy/mpl/style.py:76: UserWarning: facecolor will have no effect as it has been defined as "never".
warnings.warn('facecolor will have no effect as it has been '
