Problem 5.4#
Integrated Energy Grids
Problem 5.4
Assume we have three buses (Denmark, Sweden, and Norway) with nominal voltage \(V_{nom}\)= 2000 V connected by three transmission lines. In each of the buses, there is a gas power generator whose variable cost is 50 EUR/MWh and installed capacity is 50 MW. In the Denmark bus, there is a wind generator whose variable cost is zero and whose installed capacity is 200 MW. The transmission lines have a unitary resistance \(r\)=0.01 and reactance \(x\)=0.1, and nominal capacity \(S_{nom}=100\) VA. The demand is 50 MW for Denmark and Sweden and 30 MW for Norway. Using Python for Power System Analysis (PyPSA):
a) Calculate the optimal dispatch that minimizes the total system cost, the energy produced by each generator, and the power flows along the transmission lines using AC power flow representation.
b) Calculate the optimal dispatch that minimizes the total system cost, the energy produced by each generator, and the power flows along the transmission lines using a linearized approximation (also known as DC optimal power flow).
c) Calculate the optimal dispatch that minimizes the total system cost, the energy produced by each generator, and the power flows along the transmission lines using the Net Transfer Capacity (NTC) approach for the transmission lines and discuss the results.
Note
If you have not yet set up Python on your computer, you can execute this tutorial in your browser via Google Colab. Click on the rocket in the top right corner and launch “Colab”. If that doesn’t work download the .ipynb file and import it in Google Colab.
Then install the following packages by executing the following command in a Jupyter cell at the top of the notebook.
!pip install numpy pypsa
import numpy as np
import pypsa
We start by creating the network object and adding the three buses corresponding to Denmark, Sweden, and Norway.
network = pypsa.Network()
for node in ["Denmark", "Sweden", "Norway"]:
network.add("Bus", "bus {}".format(node), v_nom=2000)
network.buses
| v_nom | type | x | y | carrier | unit | location | v_mag_pu_set | v_mag_pu_min | v_mag_pu_max | control | generator | sub_network | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Bus | |||||||||||||
| bus Denmark | 2000.0 | 0.0 | 0.0 | AC | 1.0 | 0.0 | inf | PQ | |||||
| bus Sweden | 2000.0 | 0.0 | 0.0 | AC | 1.0 | 0.0 | inf | PQ | |||||
| bus Norway | 2000.0 | 0.0 | 0.0 | AC | 1.0 | 0.0 | inf | PQ |
We add the three lines connecting the buses
network.add("Line", "line DK-SW", bus0 = "bus Denmark", bus1 = "bus Sweden", s_nom = 100, x=0.1, r=0.01)
network.add("Line", "line DK-NO", bus0 = "bus Denmark", bus1 = "bus Norway", s_nom = 100, x=0.1, r=0.01)
network.add("Line", "line SW-NO", bus0 = "bus Sweden", bus1 = "bus Norway", s_nom = 100, x=0.1, r=0.01)
network.lines
| bus0 | bus1 | type | x | r | g | b | s_nom | s_nom_mod | s_nom_extendable | ... | v_ang_min | v_ang_max | sub_network | x_pu | r_pu | g_pu | b_pu | x_pu_eff | r_pu_eff | s_nom_opt | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Line | |||||||||||||||||||||
| line DK-SW | bus Denmark | bus Sweden | 0.1 | 0.01 | 0.0 | 0.0 | 100.0 | 0.0 | False | ... | -inf | inf | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | ||
| line DK-NO | bus Denmark | bus Norway | 0.1 | 0.01 | 0.0 | 0.0 | 100.0 | 0.0 | False | ... | -inf | inf | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | ||
| line SW-NO | bus Sweden | bus Norway | 0.1 | 0.01 | 0.0 | 0.0 | 100.0 | 0.0 | False | ... | -inf | inf | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 |
3 rows × 31 columns
We add the generators
for node in ["Denmark", "Sweden", "Norway"]:
network.add("Generator", "gas {}".format(node), bus="bus {}".format(node),
p_nom=50,
marginal_cost=50) #EUR/MWh_elec
network.add("Generator",
"wind Denmark",
bus="bus Denmark",
p_nom=200,
marginal_cost=10)
network.generators
| bus | control | type | p_nom | p_nom_mod | p_nom_extendable | p_nom_min | p_nom_max | p_min_pu | p_max_pu | ... | min_up_time | min_down_time | up_time_before | down_time_before | ramp_limit_up | ramp_limit_down | ramp_limit_start_up | ramp_limit_shut_down | weight | p_nom_opt | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Generator | |||||||||||||||||||||
| gas Denmark | bus Denmark | PQ | 50.0 | 0.0 | False | 0.0 | inf | 0.0 | 1.0 | ... | 0 | 0 | 1 | 0 | NaN | NaN | 1.0 | 1.0 | 1.0 | 0.0 | |
| gas Sweden | bus Sweden | PQ | 50.0 | 0.0 | False | 0.0 | inf | 0.0 | 1.0 | ... | 0 | 0 | 1 | 0 | NaN | NaN | 1.0 | 1.0 | 1.0 | 0.0 | |
| gas Norway | bus Norway | PQ | 50.0 | 0.0 | False | 0.0 | inf | 0.0 | 1.0 | ... | 0 | 0 | 1 | 0 | NaN | NaN | 1.0 | 1.0 | 1.0 | 0.0 | |
| wind Denmark | bus Denmark | PQ | 200.0 | 0.0 | False | 0.0 | inf | 0.0 | 1.0 | ... | 0 | 0 | 1 | 0 | NaN | NaN | 1.0 | 1.0 | 1.0 | 0.0 |
4 rows × 37 columns
We add the loads.
for node in ["Denmark", "Sweden"]:
network.add("Load", "load {}".format(node),
bus="bus {}".format(node),
p_set=50)
network.add("Load", "load Norway",
bus="bus Norway",
p_set=30)
network.loads
| bus | carrier | type | p_set | q_set | sign | active | |
|---|---|---|---|---|---|---|---|
| Load | |||||||
| load Denmark | bus Denmark | 50.0 | 0.0 | -1.0 | True | ||
| load Sweden | bus Sweden | 50.0 | 0.0 | -1.0 | True | ||
| load Norway | bus Norway | 30.0 | 0.0 | -1.0 | True |
We optimize searching for the minimum system cost.
network.optimize()
#network.optimize(solver='gurobi', assign_all_duals=True)
WARNING:pypsa.consistency:The following lines have carriers which are not defined:
Index(['line DK-SW', 'line DK-NO', 'line SW-NO'], dtype='object', name='Line')
WARNING:pypsa.consistency:The following buses have carriers which are not defined:
Index(['bus Denmark', 'bus Sweden', 'bus Norway'], dtype='object', name='Bus')
INFO:linopy.model: Solve problem using Highs solver
INFO:linopy.io: Writing time: 0.02s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 7 primals, 18 duals
Objective: 1.30e+03
Solver model: available
Solver message: Optimal
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, Kirchhoff-Voltage-Law were not assigned to the network.
Running HiGHS 1.10.0 (git hash: fd86653): Copyright (c) 2025 HiGHS under MIT licence terms
LP linopy-problem-kgagg25r has 18 rows; 7 cols; 27 nonzeros
Coefficient ranges:
Matrix [3e-03, 1e+00]
Cost [1e+01, 5e+01]
Bound [0e+00, 0e+00]
RHS [3e+01, 2e+02]
Presolving model
4 rows, 7 cols, 13 nonzeros 0s
Dependent equations search running on 4 equations with time limit of 1000.00s
Dependent equations search removed 0 rows and 0 nonzeros in 0.00s (limit = 1000.00s)
4 rows, 7 cols, 13 nonzeros 0s
Presolve : Reductions: rows 4(-14); columns 7(-0); elements 13(-14)
Solving the presolved LP
Using EKK dual simplex solver - serial
Iteration Objective Infeasibilities num(sum)
0 -6.1582334282e-06 Pr: 4(74.75) 0s
4 1.3000000000e+03 Pr: 0(0) 0s
Solving the original LP from the solution after postsolve
Model name : linopy-problem-kgagg25r
Model status : Optimal
Simplex iterations: 4
Objective value : 1.3000000000e+03
Relative P-D gap : 0.0000000000e+00
HiGHS run time : 0.00
Writing the solution to /tmp/linopy-solve-jg1ndexv.sol
('ok', 'optimal')
Now, we can look at what is the optimal dispatch form every genenerator. As expected, the wind generator is producing power to supply the demand in every node.
network.generators_t.p
| Generator | gas Denmark | gas Sweden | gas Norway | wind Denmark |
|---|---|---|---|---|
| snapshot | ||||
| now | -0.0 | -0.0 | -0.0 | 130.0 |
We can see the optimal dispatch in the generators and then solve the non-linear power flow using a Newton-Raphson method.
network.generators_t.p_set = network.generators_t.p
network.pf()
INFO:pypsa.pf:Performing non-linear load-flow on AC sub-network <pypsa.networks.SubNetwork object at 0x7fb85895b8d0> for snapshots Index(['now'], dtype='object', name='snapshot')
{'n_iter': SubNetwork 0
snapshot
now 2,
'error': SubNetwork 0
snapshot
now 7.362638e-09,
'converged': SubNetwork 0
snapshot
now True}
ok, the solution converge, we can check now the active power flow on the lines.
network.lines_t.p0
| line DK-SW | line DK-NO | line SW-NO | |
|---|---|---|---|
| snapshot | |||
| now | 43.333338 | 36.66667 | -6.666667 |
We can also check the voltage angles on the buses
network.buses_t.v_ang * 180 / np.pi
| Bus | bus Denmark | bus Sweden | bus Norway |
|---|---|---|---|
| snapshot | |||
| now | 0.0 | -0.000062 | -0.000053 |
and their per-unit mangitudes
network.buses_t.v_mag_pu
| Bus | bus Denmark | bus Sweden | bus Norway |
|---|---|---|---|
| snapshot | |||
| now | 1.0 | 1.0 | 1.0 |
b) Calculate the optimal dispatch that minimizes the total system cost, the energy produced by each generator, and the power flows along the transmission lines using a linearized approximation (also known as DC optimal power flow)
In this case, since the voltage angles are very small, the linear power flow should be a good approximation. We can calculate the power flows in the line using the linear power flow (lpf) and check that we obtained a very similar result.
network.lpf()
INFO:pypsa.pf:Performing linear load-flow on AC sub-network <pypsa.networks.SubNetwork object at 0x7fb85a5f2ad0> for snapshot(s) Index(['now'], dtype='object', name='snapshot')
network.lines_t.p0
| line DK-SW | line DK-NO | line SW-NO | |
|---|---|---|---|
| snapshot | |||
| now | 43.333333 | 36.666667 | -6.666667 |
c) Calculate the optimal dispatch that minimizes the total system cost, the energy produced by each generator, and the power flows along the transmission lines using Net Transfer Capacity (NTC) approach for the transmission lines and discuss the results.
We can create the problem again and this time use links to represent lines using only their Net Transfer Capacities. By selecting p_min_pu=-1 we make the link reversible.
network = pypsa.Network()
for node in ["Denmark", "Sweden", "Norway"]:
network.add("Bus", "bus {}".format(node), v_nom=2000)
network.add("Link","line DK-SW", bus0 = "bus Denmark", bus1 = "bus Sweden", p_nom = 100, p_min_pu=-1)
network.add("Link","line DK-NO", bus0 = "bus Denmark", bus1 = "bus Norway", p_nom =100, p_min_pu=-1)
network.add("Link","line SW-NO", bus0 = "bus Sweden", bus1 = "bus Norway", p_nom = 100, p_min_pu=-1)
for node in ["Denmark", "Sweden"]:
network.add("Load", "load {}".format(node),
bus="bus {}".format(node),
p_set=50)
network.add("Load", "load Norway",
bus="bus Norway",
p_set=30)
network.loads
for node in ["Denmark", "Sweden", "Norway"]:
network.add("Generator", "gas {}".format(node), bus="bus {}".format(node),
p_nom=50,
marginal_cost=50) #EUR/MWh_elec
network.add("Generator", "wind Denmark", bus="bus Denmark",
p_nom=200,
marginal_cost=10)
network.optimize()
WARNING:pypsa.consistency:The following links have carriers which are not defined:
Index(['line DK-SW', 'line DK-NO', 'line SW-NO'], dtype='object', name='Link')
WARNING:pypsa.consistency:The following buses have carriers which are not defined:
Index(['bus Denmark', 'bus Sweden', 'bus Norway'], dtype='object', name='Bus')
INFO:linopy.model: Solve problem using Highs solver
INFO:linopy.io: Writing time: 0.02s
INFO:linopy.constants: Optimization successful:
Status: ok
Termination condition: optimal
Solution: 7 primals, 17 duals
Objective: 1.30e+03
Solver model: available
Solver message: Optimal
INFO:pypsa.optimization.optimize:The shadow-prices of the constraints Generator-fix-p-lower, Generator-fix-p-upper, Link-fix-p-lower, Link-fix-p-upper were not assigned to the network.
Running HiGHS 1.10.0 (git hash: fd86653): Copyright (c) 2025 HiGHS under MIT licence terms
LP linopy-problem-pkzsdtf7 has 17 rows; 7 cols; 24 nonzeros
Coefficient ranges:
Matrix [1e+00, 1e+00]
Cost [1e+01, 5e+01]
Bound [0e+00, 0e+00]
RHS [3e+01, 2e+02]
Presolving model
3 rows, 7 cols, 10 nonzeros 0s
Dependent equations search running on 3 equations with time limit of 1000.00s
Dependent equations search removed 0 rows and 0 nonzeros in 0.00s (limit = 1000.00s)
3 rows, 7 cols, 10 nonzeros 0s
Presolve : Reductions: rows 3(-14); columns 7(-0); elements 10(-14)
Solving the presolved LP
Using EKK dual simplex solver - serial
Iteration Objective Infeasibilities num(sum)
0 -3.5008301912e-05 Pr: 3(470) 0s
5 1.3000000000e+03 Pr: 0(0) 0s
Solving the original LP from the solution after postsolve
Model name : linopy-problem-pkzsdtf7
Model status : Optimal
Simplex iterations: 5
Objective value : 1.3000000000e+03
Relative P-D gap : 0.0000000000e+00
HiGHS run time : 0.00
Writing the solution to /tmp/linopy-solve-pegnoqny.sol
('ok', 'optimal')
network.generators_t.p
| Generator | gas Denmark | gas Sweden | gas Norway | wind Denmark |
|---|---|---|---|---|
| snapshot | ||||
| now | -0.0 | -0.0 | -0.0 | 130.0 |
network.links_t.p0
| Link | line DK-SW | line DK-NO | line SW-NO |
|---|---|---|---|
| snapshot | |||
| now | 100.0 | -20.0 | 50.0 |
In this case, the power flows are also compatible with the nodal balances, but they are significantly different from those obtained using AC power flow or linearized AC power flow to represent the power flowing through the different lines.