Solving a Power Flow
To follow along, you can download this tutorial as a Julia script (.jl) or Jupyter notebook (.ipynb).
In this tutorial, you'll solve power flows on a 5-bus test system using three different solvers and compare their results.
Building a System
To get started, load the needed packages. We're using a standard test system and want to keep output clean, so we adjust the logging settings to filter out a few precautionary warnings.
using PowerSystemCaseBuilder
using PowerFlows
using PowerSystems
using Logging
disable_logging(Logging.Warn)LogLevel(1001)Create a System from PowerSystemCaseBuilder.jl:
sys = build_system(MatpowerTestSystems, "matpower_case5_sys")| System | |
| Property | Value |
|---|---|
| Name | |
| Description | |
| System Units Base | SYSTEM_BASE |
| Base Power | 100.0 |
| Base Frequency | 60.0 |
| Num Components | 28 |
| Static Components | |
| Type | Count |
|---|---|
| ACBus | 5 |
| Arc | 6 |
| Area | 1 |
| Line | 5 |
| LoadZone | 1 |
| PhaseShiftingTransformer | 2 |
| PowerLoad | 3 |
| ThermalStandard | 5 |
DC Power Flow
DCPowerFlow solves for bus voltage angles using the bus admittance matrix, then computes branch flows from the angle differences. Create a DCPowerFlow solver:
pf_dc = DCPowerFlow()DCPowerFlow(nothing, PowerNetworkMatrices.NetworkReduction[], 1, String[], false, false)Solve the power flow with solve_power_flow. For DC methods, pass a FlowReporting mode so flows are reported on the same arc basis as elsewhere in the package:
dc_results = solve_power_flow(pf_dc, sys, FlowReporting.ARC_FLOWS)Dict{String, Dict{String, DataFrame}} with 1 entry:
"1" => Dict("flow_results"=>6×10 DataFrame…The result is a Dict{String, Dict{String, DataFrame}}. The outer key is the time step name: "1". The inner dictionary stores the power flow results at that time step: "bus_results" for bus data and "flow_results" for AC line data. (There is also a third key, "lcc_results", for HVDC lines, but this system contains no such components, so the matching dataframe will be empty.) Inspect "bus_results":
dc_results["1"]["bus_results"]| bus_number | Vm | θ | P_gen | P_load | P_net | Q_gen | Q_load | Q_net |
|---|---|---|---|---|---|---|---|---|
| Int64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
| 1 | 1.0 | 0.05813712342707429 | 210.0 | 0.0 | 210.0 | 0.0 | 0.0 | 0.0 |
| 2 | 1.0 | -0.010697105096729564 | 0.0 | 300.0 | -300.0 | 0.0 | 0.0 | 0.0 |
| 3 | 1.0 | -0.004752965105551588 | 324.498 | 300.0 | 24.497999999999998 | 0.0 | 0.0 | 0.0 |
| 4 | 1.0 | 0.0 | 0.0 | 400.0 | -400.0 | 0.0 | 0.0 | 0.0 |
| 10 | 1.0 | 0.07261406474007273 | 470.694 | 0.0 | 470.694 | 0.0 | 0.0 | 0.0 |
Notice that Vm (voltage magnitude) is 1.0 for all buses, and Q_gen and Q_load are 0. This is expected for DC power flow, which assumes flat voltage magnitudes and ignores reactive power.
dc_results["1"]["flow_results"]| flow_name | bus_from | bus_to | P_from_to | Q_from_to | P_to_from | Q_to_from | P_losses | Q_losses | angle_difference |
|---|---|---|---|---|---|---|---|---|---|
| String | Int64 | Int64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
| 1-2 | 1 | 2 | 244.9616680920614 | 0.0 | -244.9616680920614 | 0.0 | 1.6861747492479113 | 0.0 | 0.06883422852380386 |
| 1-4 | 1 | 4 | 191.2405399235394 | 0.0 | -191.2405399235394 | 0.0 | 1.1118175009515048 | 0.0 | 0.05813712342707429 |
| 1-10 | 1 | 10 | -226.20220801560063 | 0.0 | 226.20220801560063 | 0.0 | 0.32747160903125155 | 0.0 | -0.01447694131299844 |
| 2-3 | 2 | 3 | -55.03833190793862 | 0.0 | 55.03833190793862 | 0.0 | 0.03271555417545088 | 0.0 | -0.005944139991177976 |
| 3-4 | 3 | 4 | -30.5403319079386 | 0.0 | 30.5403319079386 | 0.0 | 0.03054095074902095 | 0.0 | -0.004752965105551588 |
| 4-10 | 4 | 10 | -244.49179198439933 | 0.0 | 244.49179198439933 | 0.0 | 1.7753542195279608 | 0.0 | -0.07261406474007273 |
Likewise, Q_from_to and Q_to_from (reactive power flow on the line) are zero, for all lines.
PTDF DC Power Flow
PTDFDCPowerFlow computes branch flows directly from bus power injections using the Power Transfer Distribution Factor matrix, without solving for voltage angles as an intermediate step. (This means we can omit the angle computation in contexts where we only care about line flows, though we don't have that option implemented here.) Create a PTDFDCPowerFlow solver:
pf_ptdf = PTDFDCPowerFlow()PTDFDCPowerFlow(nothing, false, PowerNetworkMatrices.NetworkReduction[], 1, String[], false)As before, solve the power flow with solve_power_flow:
ptdf_results = solve_power_flow(pf_ptdf, sys, FlowReporting.ARC_FLOWS)Dict{String, Dict{String, DataFrame}} with 1 entry:
"1" => Dict("flow_results"=>6×10 DataFrame…Look at the bus results:
ptdf_results["1"]["bus_results"]| bus_number | Vm | θ | P_gen | P_load | P_net | Q_gen | Q_load | Q_net |
|---|---|---|---|---|---|---|---|---|
| Int64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
| 1 | 1.0 | 0.05813712342707429 | 210.0 | 0.0 | 210.0 | 0.0 | 0.0 | 0.0 |
| 2 | 1.0 | -0.010697105096729564 | 0.0 | 300.0 | -300.0 | 0.0 | 0.0 | 0.0 |
| 3 | 1.0 | -0.004752965105551588 | 324.498 | 300.0 | 24.497999999999998 | 0.0 | 0.0 | 0.0 |
| 4 | 1.0 | 0.0 | 0.0 | 400.0 | -400.0 | 0.0 | 0.0 | 0.0 |
| 10 | 1.0 | 0.07261406474007273 | 470.694 | 0.0 | 470.694 | 0.0 | 0.0 | 0.0 |
The results match DCPowerFlow, as they should: the two are mathematically equivalent. For very large systems where forming the full PTDF matrix would be too expensive, consider vPTDFDCPowerFlow, which computes the same results without storing the dense matrix.
AC Power Flow
Create an ACPowerFlow solver:
pf_ac = ACPowerFlow()ACPowerFlow{NewtonRaphsonACPowerFlow}(false, nothing, false, false, nothing, true, false, false, false, PowerNetworkMatrices.NetworkReduction[], 1, String[], false, Dict{Symbol, Any}())Solve the power flow:
ac_results = solve_power_flow(pf_ac, sys)Dict{String, DataFrame} with 3 entries:
"flow_results" => 6×10 DataFrame…
"bus_results" => 5×9 DataFrame…
"lcc_results" => 0×13 DataFrame…AC results are returned as a flat Dict{String, DataFrame}, with the same keys as before: "bus_results", "flow_results" (AC lines), and "lcc_results" (HVDC lines). (We don't support multi-period AC power flows yet.) Look at the bus results:
ac_results["bus_results"]| bus_number | Vm | θ | P_gen | P_load | P_net | Q_gen | Q_load | Q_net |
|---|---|---|---|---|---|---|---|---|
| Int64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
| 1 | 1.07762 | 0.05377973212862824 | 210.0 | 0.0 | 210.0 | 155.86199106764073 | 0.0 | 155.86199106764073 |
| 2 | 1.0841946221315348 | 8.152349463538256e-5 | 0.0 | 300.0 | -300.0 | 0.0 | 98.61 | -98.61 |
| 3 | 1.1 | 0.006160589908746271 | 324.498 | 300.0 | 24.497999999999998 | 518.182149466319 | 98.61 | 419.57214946631893 |
| 4 | 1.06414 | 0.0 | 2.2015611506636112 | 400.0 | -397.7984388493364 | -116.00444489674346 | 131.47 | -247.47444489674345 |
| 10 | 1.06907 | 0.06665597012564134 | 470.694 | 0.0 | 470.694 | -164.31698920802899 | 0.0 | -164.31698920802899 |
Notice that Vm now varies across buses (not all 1.0), and Q_gen has non-zero values.
Look at the line flows:
ac_results["flow_results"]| flow_name | bus_from | bus_to | P_from_to | Q_from_to | P_to_from | Q_to_from | P_losses | Q_losses | angle_difference |
|---|---|---|---|---|---|---|---|---|---|
| String | Int64 | Int64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 | Float64 |
| 1-2 | 1 | 2 | 219.04815651765563 | -41.538343627190756 | -217.84617277462783 | 52.726508142521936 | 1.2019837430278013 | 11.188164515331177 | 0.053698208633992855 |
| 1-4 | 1 | 4 | 206.03176235568074 | 32.25263490048205 | -204.89263500836265 | -21.615578151334624 | 1.139127347318114 | 10.637056749147423 | 0.05377973212862824 |
| 1-10 | 1 | 10 | -215.0798358132005 | 165.14902875652393 | 215.48841730273134 | -164.66348802499908 | 0.40858148953084594 | 0.4855407315248472 | -0.012876237997013103 |
| 2-3 | 2 | 3 | -82.153855250951 | -151.33695655178607 | 82.4232748938078 | 151.82170300927572 | 0.2694196428567941 | 0.4847464574896376 | -0.006079066414110889 |
| 3-4 | 3 | 4 | -57.925274893807945 | 267.75044645704486 | 60.60738862294479 | -241.67008678435033 | 2.682113729136848 | 26.08035967269453 | 0.006160589908746271 |
| 4-10 | 4 | 10 | -253.51319246391864 | 15.812083986954987 | 255.20569169372328 | 0.3464988169711386 | 1.6924992298046337 | 16.158582803926123 | -0.06665597012564134 |
Q_from_to and Q_to_from now show reactive power flows, and P_from_to differs from P_to_from due to losses.
When AC Power Flow Fails
Unlike DC power flow, AC power flow is iterative and not guaranteed to converge. Systems with high impedance lines, poor initial voltage profiles, or insufficient reactive power support can cause the solver to fail. When this happens, solve_power_flow returns missing: you'll also see a logged error. If you encounter convergence failures, consider using a more robust solver such as TrustRegionACPowerFlow or RobustHomotopyPowerFlow.