AC Power Flow
In this example, we perform several AC power flow analyses using the power system shown in Figure 1. These analyses simulate quasi-steady-state conditions where the system undergoes parameter and topology changes, demonstrating JuliaGrid's efficiency in handling such scenarios.
Figure 1: The 4-bus power system.
Users can download a Julia script containing the scenarios from this section using the following link.
We begin by defining the units for active and reactive powers, as well as voltage magnitudes and angles, which will be used throughout this example:
@power(MW, MVAr)
@voltage(pu, deg)Next, we define the bus parameters for AC power flow analysis. This includes specifying the type of each bus, the connected active and reactive power loads, and shunt capacitor banks with conductance and susceptance values. The bus voltage magnitude and angle serve as initial values for the iterative power flow algorithm. Note that for the slack bus (type = 3), the angle is fixed to the specified value. With these definitions, we can start to build the power system model:
system = powerSystem()
@bus(magnitude = 1.1, angle = -5.7)
addBus!(system; label = "Bus 1", type = 3, angle = 0.0)
addBus!(system; label = "Bus 2", type = 2, active = 20.2, reactive = 10.5)
addBus!(system; label = "Bus 3", type = 1, conductance = 0.1, susceptance = 8.2)
addBus!(system; label = "Bus 4", type = 1, active = 50.8, reactive = 23.1)Next, we refine the transmission line parameters by adding resistance, reactance, and susceptance values. Additionally, for transformers, we specify the off-nominal turns ratio using the turnsRatio keyword:
@branch(label = "Branch ?", reactance = 0.22)
addBranch!(system; from = "Bus 1", to = "Bus 3", resistance = 0.02, susceptance = 0.05)
addBranch!(system; from = "Bus 1", to = "Bus 2", resistance = 0.05, susceptance = 0.04)
addBranch!(system; from = "Bus 2", to = "Bus 3", resistance = 0.04, susceptance = 0.04)
addBranch!(system; from = "Bus 3", to = "Bus 4", turnsRatio = 0.98)Finally, we define the active and reactive power outputs of the generators and the voltage magnitude setpoints. These setpoints fix the voltage magnitudes for the slack bus (type = 3) and generator buses (type = 2):
@generator(label = "Generator ?")
addGenerator!(system; bus = "Bus 1", active = 60.1, reactive = 40.2, magnitude = 0.98)
addGenerator!(system; bus = "Bus 2", active = 18.2, magnitude = 1.01)After defining the power system data, we generate an AC model with the vectors and matrices required for analysis, including the nodal admittance matrix. This model is automatically updated when the data changes and can be reused across analyses:
acModel!(system)Display Data Settings
Before running simulations, we set verbose from silent mode (0) to basic output (1):
@config(verbose = 1)For more detailed solver output, verbose can also be set when calling the power flow solver. Next, we configure the data display settings, including the selection of displayed data elements and the numeric format for relevant power flow values.
For bus-related data, we set:
show1 = Dict("Power Injection" => false)
fmt1 = Dict("Power Generation" => "%.2f", "Power Demand" => "%.2f", "Shunt Power" => "%.2f")Similarly, for branch-related data, we choose:
show2 = Dict("Shunt Power" => false, "Status" => false)
fmt2 = Dict("From-Bus Power" => "%.2f", "To-Bus Power" => "%.2f", "Series Power" => "%.2f")Base Case Analysis
At the start, we create a fast Newton-Raphson XB model:
fnr = fastNewtonRaphsonXB(system)Then, we solve the AC power flow and compute bus voltages together with active and reactive powers:
powerFlow!(fnr; power = true, verbose = 2)Number of entries in the Jacobians: 11
Active Power: 7
Reactive Power: 4
Number of state variables: 5
---------------------------------------------------------------
Iteration Maximum Active Mismatch Maximum Reactive Mismatch
---------------------------------------------------------------
0 4.61818182e-01 9.63142262e-01
1 1.77748104e-01 3.36304103e-02
2 1.91631619e-02 3.40094647e-03
3 2.44565654e-03 4.00326340e-04
4 3.06356277e-04 4.86320320e-05
5 3.76577214e-05 5.92151978e-06
6 4.59272704e-06 7.20154616e-07
7 5.58645482e-07 8.75181686e-08
8 6.78952866e-08 1.06328286e-08
9 8.24944624e-09 1.29169581e-09
Minimum Value Maximum Value
Magnitude Increment: 5.6585e-10 5.7490e-09
Angle Increment: 1.5264e-09 4.9430e-09
EXIT: The solution was found using the fast Newton-Raphson method in 9 iterations.Once the AC power flow is solved, we can inspect the bus results. For instance:
printBusData(fnr; show = show1, fmt = fmt1)|------------------------------------------------------------------------------------------|
| Bus Data |
|------------------------------------------------------------------------------------------|
| Label | Voltage | Power Generation | Power Demand | Shunt Power |
| | | | | |
| Bus | Magnitude | Angle | Active | Reactive | Active | Reactive | Active | Reactive |
| | [pu] | [deg] | [MW] | [MVAr] | [MW] | [MVAr] | [MW] | [MVAr] |
|-------|-----------|----------|--------|----------|--------|----------|--------|----------|
| Bus 1 | 0.9800 | 0.0000 | 53.73 | -16.12 | 0.00 | 0.00 | 0.00 | -0.00 |
| Bus 2 | 1.0100 | -2.8924 | 18.20 | 43.22 | 20.20 | 10.50 | 0.00 | -0.00 |
| Bus 3 | 0.9641 | -4.5972 | 0.00 | 0.00 | 0.00 | 0.00 | 0.09 | -7.62 |
| Bus 4 | 0.9211 | -11.6813 | 0.00 | 0.00 | 50.80 | 23.10 | 0.00 | -0.00 |
|------------------------------------------------------------------------------------------|Similarly, we can inspect the branch results:
printBranchData(fnr; show = show2, fmt = fmt2)|------------------------------------------------------------------------------------------|
| Branch Data |
|------------------------------------------------------------------------------------------|
| Label | From-Bus Power | To-Bus Power | Series Power |
| | | | |
| Branch | From-Bus | To-Bus | Active | Reactive | Active | Reactive | Active | Reactive |
| | | | [MW] | [MVAr] | [MW] | [MVAr] | [MW] | [MVAr] |
|----------|----------|--------|--------|----------|--------|----------|--------|----------|
| Branch 1 | Bus 1 | Bus 3 | 34.90 | 2.87 | -34.64 | -4.74 | 0.26 | 2.85 |
| Branch 2 | Bus 1 | Bus 2 | 18.82 | -18.99 | -18.49 | 16.51 | 0.34 | 1.48 |
| Branch 3 | Bus 2 | Bus 3 | 16.49 | 16.21 | -16.25 | -18.81 | 0.24 | 1.30 |
| Branch 4 | Bus 3 | Bus 4 | 50.80 | 31.17 | -50.80 | -23.10 | 0.00 | 8.07 |
|------------------------------------------------------------------------------------------|The resulting active and reactive power flows are shown in Figure 2.
(a) Active powers.
(b) Reactive powers.
Figure 2: Power flows in the 4-bus power system for the base case scenario.
The branch data also shows the series power losses caused by the series resistance and reactance of each branch. Note that the active power at the from-bus and to-bus ends of a branch differs by the active power loss. Reactive power does not follow the same pattern because branch susceptances provide partial compensation.
Modifying Supplies and Demands
We now modify the active and reactive power outputs of the generators and the active and reactive demands. Rather than creating a new model, we update the power system and fast Newton-Raphson models simultaneously:
updateBus!(fnr; label = "Bus 2", active = 25.5, reactive = 15.0)
updateBus!(fnr; label = "Bus 4", active = 42.0, reactive = 20.0)
updateGenerator!(fnr; label = "Generator 1", active = 58.0, reactive = 20.0)
updateGenerator!(fnr; label = "Generator 2", active = 23.1, reactive = 20.0)Next, we run the AC power flow again without recreating the fast Newton-Raphson model. This call also warm-starts the fast Newton-Raphson method because the initial voltage magnitudes and angles come from the base-case solution:
powerFlow!(fnr; power = true)EXIT: The solution was found using the fast Newton-Raphson method in 8 iterations.Since these changes do not affect the Jacobian matrices, JuliaGrid reuses the Jacobian factorizations from the base case and reduces the computational cost.
Finally, we display the updated branch data:
printBranchData(fnr; show = show2, fmt = fmt2)|------------------------------------------------------------------------------------------|
| Branch Data |
|------------------------------------------------------------------------------------------|
| Label | From-Bus Power | To-Bus Power | Series Power |
| | | | |
| Branch | From-Bus | To-Bus | Active | Reactive | Active | Reactive | Active | Reactive |
| | | | [MW] | [MVAr] | [MW] | [MVAr] | [MW] | [MVAr] |
|----------|----------|--------|--------|----------|--------|----------|--------|----------|
| Branch 1 | Bus 1 | Bus 3 | 29.23 | -0.96 | -29.05 | -1.85 | 0.18 | 1.96 |
| Branch 2 | Bus 1 | Bus 2 | 15.87 | -18.46 | -15.59 | 15.71 | 0.27 | 1.20 |
| Branch 3 | Bus 2 | Bus 3 | 13.19 | 12.67 | -13.04 | -15.76 | 0.15 | 0.84 |
| Branch 4 | Bus 3 | Bus 4 | 42.00 | 25.37 | -42.00 | -20.00 | 0.00 | 5.37 |
|------------------------------------------------------------------------------------------|Compared with the base case, the active power flow directions remain unchanged, but their magnitudes differ. The reactive power values also change, and the flow at the Bus 1 side of Branch 1 reverses, as shown in Figure 3.
(a) Active powers.
(b) Reactive powers.
Figure 3: Power flows in the 4-bus power system with modified supplies and demands.
Modifying Network Topology
Next, we take Branch 3 out-of-service. Although we could update the power system and fast Newton-Raphson models simultaneously, we use the Newton-Raphson method here to demonstrate flexibility. Therefore, we update only the power system model:
updateBranch!(system; label = "Branch 3", status = 0)Next, we create the Newton-Raphson model:
nr = newtonRaphson(system)When the model is created, the Newton-Raphson method is initialized with the voltage magnitudes and angles defined in the original power system model. To warm-start it from the fast Newton-Raphson solution, we transfer the voltage magnitudes and angles:
setInitialPoint!(nr, fnr)We can now solve the AC power flow for this scenario:
powerFlow!(nr; power = true)EXIT: The solution was found using the Newton-Raphson method in 3 iterations.To inspect the power flow distribution after the outage, we print the branch data:
printBranchData(nr; show = show2, fmt = fmt2)|------------------------------------------------------------------------------------------|
| Branch Data |
|------------------------------------------------------------------------------------------|
| Label | From-Bus Power | To-Bus Power | Series Power |
| | | | |
| Branch | From-Bus | To-Bus | Active | Reactive | Active | Reactive | Active | Reactive |
| | | | [MW] | [MVAr] | [MW] | [MVAr] | [MW] | [MVAr] |
|----------|----------|--------|--------|----------|--------|----------|--------|----------|
| Branch 1 | Bus 1 | Bus 3 | 42.56 | 19.70 | -42.09 | -18.98 | 0.48 | 5.27 |
| Branch 2 | Bus 1 | Bus 2 | 2.50 | -15.82 | -2.40 | 12.31 | 0.10 | 0.46 |
| Branch 3 | Bus 2 | Bus 3 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 |
| Branch 4 | Bus 3 | Bus 4 | 42.00 | 26.01 | -42.00 | -20.00 | 0.00 | 6.01 |
|------------------------------------------------------------------------------------------|Compared with the previous cases, the reactive power flow at the Bus 1 side of Branch 1 reverses because of the Branch 3 outage, as shown in Figure 4.
(a) Active powers.
(b) Reactive powers.
Figure 4: Power flows in the 4-bus power system with modified network topology.