# 🔥 Power flow
VeraGrid has the most power flow features in any open-source software.
This is the power flow settings page:
🎥 [Watch the IEEE14 demo on YouTube](https://www.youtube.com/watch?v=d9dTyO4wels)
Below, a list of the parameters and solvers available.
Settings:
- **Solver**: The power flow solver to use.
- Newton-Raphson in power:
- Newton-Raphson in current:
- Newton-Raphson-Iwamoto:
- Levenberg-Marquardt:
- Fast-Decoupled:
- Holomorphic-Embedding:
- Linear "AC" approximation:
- "DC" approximation
The following table relates which power flow controls are available in each solver:
| | Newton
Raphson |Powell
Dog-leg|Levenberg-
Marquardt|Iwamoto|Fast-
decoupled|Gauss-
seidel|Holomorphic
embedding|Linear
without
voltage modules|Linear
with voltage
modules|
|-----------------------------------------------------------------------------------------|---|---|---|---|---|---|---|---|---|
| Local voltage
control using
a Generator. | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ |
| Remote voltage
control using
a Generator. | ✅ | ✅ | ✅ | ✅ | ✅ | | | | |
| Generator
reactive power
limits. | ✅ | ✅ | ✅ | ✅ | ✅ | ✅ | | | |
| Local and
remote voltage
control using
a transformer's
tap changer. | ✅ | ✅ | ✅ | | | | | | |
| Local active
power control
using a
transformer's
tap changer. | ✅ | ✅ | ✅ | | | | | | |
| Local reactive
power control
using a
transformer's
tap changer. | ✅ | ✅ | ✅ | | | | | | |
| Local and
remote
AC and DC
voltage control
using a
converter. | ✅ | ✅ | ✅ | | | | | | |
| Local
AC and DC
active power
control using
a converter. |✅ | ✅ | ✅ | | | | | | |
| Local AC
reactive power
control using
a converter. | ✅ | ✅ | ✅ | | | | | | |
| 3-phase
unbalanced. | ✅ | ✅ | ✅ | | | | | | |
- **Tolerance**: per-unit error tolerance to use in the solver. Exponent of the numerical precision. i.e. `4` corresponds to `1e-4` MW in p.u. of precision
- **Automatic precision (find)**:
The precision to use for the numerical solvers depends on the magnitude of the power injections.
If we are dealing with hundreds of MW, the precision may be `1e-3`, but if we are dealing with Watts, the precision has
to be greater. The automatic precision checks the loading for a suitable precision such that the results are fine.
- **Trust radius**: For newton-like methods this is the solution trust radius to use
- **Max. Iterations**: Maximum number of iterations in the solver.
- **Verbosity**: Level of verbosity (0: none, 1: some, >=2: all)
Flags:
- **Retry with other methods is failed**:
This option tries other numerical solvers to try to find a power flow solution.
This option is relevant because different numerical algorithms may be more suited to certain grid configurations.
In general the Newton-Raphson implementation in VeraGrid includes back-tracing and other innovations that make it
a very competitive method to consider by default.
- **Use initial guess**: In the buses there are two properties `Vm0`and `Va0` both serve to
specify the initial voltage guess. If selected those values are used to initialize the power flow solution.
- **Ignore single island nodes**: If selected, the islands of 1 one are considered in blackout and not calculated.
- **3-phase**: If selected a three phase power flow is run. For now only supporting simple local voltage controls.
- **Distributed slack**: If selected, the slack power is distributed simply among the existing generators.
- **Control Q limits**: This is the mode of reactive power control for the generators that are set in PV mode.
This is a node-level control, so for this to work, there must not be more than one generator per bus. Otherwise a
per-node repartition will be made, but the limits respecting cannot be ensured. For that use the AC-OPF.
- **Control tap module**: If selected the tap module controls are possible, if active at device level.
- **Control tap phase**: If selected the tap phase controls are possible, if active at device level.
- **Control remote voltages**: If selected, the remote voltage controls are possible if set at device level.
- **Orthogonalize controls**: If selected, the controls are "orthogonalized" to their device specified steps.
This applies to transformer tap changers and non-linear shunts.
- **Initialize angles**: If selected the power flow angles are initialized with a linear power flow.
For Holomorphic embedding this is not necessary.
- **Apply temperature correction**: When selected the branches apply the correction of the resistance due to the temperature.
- **Apply impedance tolerances**:
- **Add report**: Inspect the results looking for violations to report in the power flow logger.
## Registered Result Properties
### `PowerFlowTimeSeriesResults` registered properties
The power flow time-series result stores the solved network quantities for every time step.
| Property | Type | Description |
|----------|------|-------------|
| `bus_names` | `StrVec` | Names aligned with bus-indexed result arrays. |
| `branch_names` | `StrVec` | Names aligned with branch-indexed result arrays. |
| `hvdc_names` | `StrVec` | Names aligned with HVDC line-indexed result arrays. |
| `bus_types` | `IntVec` | Bus type code used by the solved numerical model. |
| `time_array` | `DateVec` | Time stamps represented by the result object. |
| `F` | `IntVec` | Branch from-bus index for each branch. |
| `T` | `IntVec` | Branch to-bus index for each branch. |
| `hvdc_F` | `IntVec` | HVDC from-bus index for each HVDC line. |
| `hvdc_T` | `IntVec` | HVDC to-bus index for each HVDC line. |
| `bus_area_indices` | `IntVec` | Area index assigned to each bus. |
| `area_names` | `IntVec` | Area names or area identifiers used for inter-area aggregation. |
| `S` | `CxMat` | Complex bus power result matrix. |
| `voltage` | `CxMat` | Complex bus voltage solution. |
| `Sf` | `CxMat` | Complex branch power flow at the from side. |
| `St` | `CxMat` | Complex branch power flow at the to side. |
| `If` | `CxMat` | Complex branch current at the from side. |
| `It` | `CxMat` | Complex branch current at the to side. |
| `tap_module` | `Mat` | Transformer tap module used in the solved state. |
| `tap_angle` | `Mat` | Transformer tap angle used in the solved state. |
| `Vbranch` | `CxMat` | Complex branch voltage result used by branch reports. |
| `loading` | `CxMat` | Branch loading result. |
| `losses` | `CxMat` | Complex branch losses. |
| `hvdc_losses` | `Mat` | HVDC result field `hvdc_losses`. |
| `hvdc_Pf` | `Mat` | HVDC result field `hvdc_Pf`. |
| `hvdc_Pt` | `Mat` | HVDC result field `hvdc_Pt`. |
| `hvdc_loading` | `Mat` | HVDC result field `hvdc_loading`. |
### `PowerFlowResults3Ph` registered properties
The three-phase power flow result stores phase-resolved bus, branch, HVDC, VSC, generator, battery, and shunt quantities.
| Property | Type | Description |
|----------|------|-------------|
| `bus_names` | `StrVec` | Names aligned with bus-indexed result arrays. |
| `branch_names` | `StrVec` | Names aligned with branch-indexed result arrays. |
| `hvdc_names` | `StrVec` | Names aligned with HVDC line-indexed result arrays. |
| `gen_names` | `StrVec` | Names aligned with generator-indexed result arrays. |
| `batt_names` | `StrVec` | Names aligned with battery-indexed result arrays. |
| `sh_names` | `StrVec` | Names aligned with shunt-indexed result arrays. |
| `load_names` | `StrVec` | Names aligned with load-indexed result arrays. |
| `bus_types` | `IntVec` | Bus type code used by the solved numerical model. |
| `F` | `IntVec` | Branch from-bus index for each branch. |
| `T` | `IntVec` | Branch to-bus index for each branch. |
| `hvdc_F` | `IntVec` | HVDC from-bus index for each HVDC line. |
| `hvdc_T` | `IntVec` | HVDC to-bus index for each HVDC line. |
| `bus_area_indices` | `IntVec` | Area index assigned to each bus. |
| `area_names` | `IntVec` | Area names or area identifiers used for inter-area aggregation. |
| `Sbus_N` | `CxVec` | Complex bus power injection for neutral. |
| `Sbus_A` | `CxVec` | Complex bus power injection for phase A. |
| `Sbus_B` | `CxVec` | Complex bus power injection for phase B. |
| `Sbus_C` | `CxVec` | Complex bus power injection for phase C. |
| `voltage_N` | `CxVec` | Complex bus voltage solution for neutral. |
| `voltage_A` | `CxVec` | Complex bus voltage solution for phase A. |
| `voltage_B` | `CxVec` | Complex bus voltage solution for phase B. |
| `voltage_C` | `CxVec` | Complex bus voltage solution for phase C. |
| `Sf_A` | `CxVec` | Complex branch power flow at the from side for phase A. |
| `Sf_B` | `CxVec` | Complex branch power flow at the from side for phase B. |
| `Sf_C` | `CxVec` | Complex branch power flow at the from side for phase C. |
| `St_A` | `CxVec` | Complex branch power flow at the to side for phase A. |
| `St_B` | `CxVec` | Complex branch power flow at the to side for phase B. |
| `St_C` | `CxVec` | Complex branch power flow at the to side for phase C. |
| `If_N` | `CxVec` | Complex branch current at the from side for neutral. |
| `If_A` | `CxVec` | Complex branch current at the from side for phase A. |
| `If_B` | `CxVec` | Complex branch current at the from side for phase B. |
| `If_C` | `CxVec` | Complex branch current at the from side for phase C. |
| `It_N` | `CxVec` | Complex branch current at the to side for neutral. |
| `It_A` | `CxVec` | Complex branch current at the to side for phase A. |
| `It_B` | `CxVec` | Complex branch current at the to side for phase B. |
| `It_C` | `CxVec` | Complex branch current at the to side for phase C. |
| `tap_module` | `Vec` | Transformer tap module used in the solved state. |
| `tap_angle` | `Vec` | Transformer tap angle used in the solved state. |
| `Vbranch_A` | `CxVec` | Complex branch voltage result for phase A. |
| `Vbranch_B` | `CxVec` | Complex branch voltage result for phase B. |
| `Vbranch_C` | `CxVec` | Complex branch voltage result for phase C. |
| `loading_A` | `CxVec` | Loading result for phase A. |
| `loading_B` | `CxVec` | Loading result for phase B. |
| `loading_C` | `CxVec` | Loading result for phase C. |
| `losses_A` | `CxVec` | Losses result for phase A. |
| `losses_B` | `CxVec` | Losses result for phase B. |
| `losses_C` | `CxVec` | Losses result for phase C. |
| `losses_hvdc` | `Vec` | Losses result. |
| `Pf_hvdc_A` | `Vec` | Registered result field `Pf_hvdc_A`. |
| `Pf_hvdc_B` | `Vec` | Registered result field `Pf_hvdc_B`. |
| `Pf_hvdc_C` | `Vec` | Registered result field `Pf_hvdc_C`. |
| `Pt_hvdc_A` | `Vec` | Registered result field `Pt_hvdc_A`. |
| `Pt_hvdc_B` | `Vec` | Registered result field `Pt_hvdc_B`. |
| `Pt_hvdc_C` | `Vec` | Registered result field `Pt_hvdc_C`. |
| `loading_hvdc` | `Vec` | Loading result. |
| `losses_vsc` | `Vec` | Losses result. |
| `Pfp_vsc` | `Vec` | VSC result field `Pfp_vsc`. |
| `St_vsc_A` | `CxVec` | Complex branch power flow at the to side for phase A. |
| `St_vsc_B` | `CxVec` | Complex branch power flow at the to side for phase B. |
| `St_vsc_C` | `CxVec` | Complex branch power flow at the to side for phase C. |
| `If_vsc` | `Vec` | Complex branch current at the from side. |
| `It_vsc_A` | `CxVec` | Complex branch current at the to side for phase A. |
| `It_vsc_B` | `CxVec` | Complex branch current at the to side for phase B. |
| `It_vsc_C` | `CxVec` | Complex branch current at the to side for phase C. |
| `loading_vsc` | `Vec` | Loading result. |
| `gen_q_A` | `Vec` | Generator reactive power output for phase A. |
| `gen_q_B` | `Vec` | Generator reactive power output for phase B. |
| `gen_q_C` | `Vec` | Generator reactive power output for phase C. |
| `battery_q_A` | `Vec` | Battery reactive power output for phase A. |
| `battery_q_B` | `Vec` | Battery reactive power output for phase B. |
| `battery_q_C` | `Vec` | Battery reactive power output for phase C. |
| `shunt_q_A` | `Vec` | Shunt reactive power output for phase A. |
| `shunt_q_B` | `Vec` | Shunt reactive power output for phase B. |
| `shunt_q_C` | `Vec` | Shunt reactive power output for phase C. |
## API
Using the simplified API:
```python
import os
import VeraGridEngine as vg
folder = os.path.join('..', 'Grids_and_profiles', 'grids')
fname = os.path.join(folder, 'IEEE39_1W.veragrid')
main_circuit = vg.open_file(fname)
results = vg.power_flow(main_circuit)
print(main_circuit.name)
print('Converged:', results.converged, 'error:', results.error)
print(results.get_bus_df())
print(results.get_branch_df())
```
Using the more complex library objects:
```python
import os
import VeraGridEngine as vg
folder = os.path.join('..', 'Grids_and_profiles', 'grids')
fname = os.path.join(folder, 'IEEE14_from_raw.veragrid')
main_circuit = vg.open_file(fname)
options = vg.PowerFlowOptions(vg.SolverType.NR, verbose=False)
power_flow = vg.PowerFlowDriver(main_circuit, options)
power_flow.run()
print(main_circuit.name)
print('Converged:', power_flow.results.converged, 'error:', power_flow.results.error)
print(power_flow.results.get_bus_df())
print(power_flow.results.get_branch_df())
```
### `PowerFlowResults` registered properties
The `PowerFlowResults` object stores the following registered properties for persistence and programmatic access.
The arrays use the grid ordering used by the numerical circuit: buses follow `bus_names`, branches follow
`branch_names`, HVDC links follow `hvdc_names`, generators follow `gen_names`, batteries follow `batt_names`,
and shunts follow `sh_names`.
| Property | Type | Description |
|----------|------|-------------|
| `bus_names` | `StrVec` | Bus names aligned with bus-indexed result arrays. |
| `branch_names` | `StrVec` | Branch names aligned with branch-indexed result arrays. |
| `hvdc_names` | `StrVec` | HVDC line names aligned with HVDC-indexed result arrays. |
| `gen_names` | `StrVec` | Generator names aligned with generator-indexed result arrays. |
| `batt_names` | `StrVec` | Battery names aligned with battery-indexed result arrays. |
| `sh_names` | `StrVec` | Shunt names aligned with shunt-indexed result arrays. |
| `bus_types` | `IntVec` | Bus type code used by the solved numerical model. |
| `F` | `IntVec` | Branch "from" bus index for each branch. |
| `T` | `IntVec` | Branch "to" bus index for each branch. |
| `hvdc_F` | `IntVec` | HVDC "from" bus index for each HVDC line. |
| `hvdc_T` | `IntVec` | HVDC "to" bus index for each HVDC line. |
| `bus_area_indices` | `IntVec` | Area index assigned to each bus. |
| `area_names` | `IntVec` | Area names or area identifiers used for inter-area aggregation. |
| `Sbus` | `CxVec` | Complex bus power injection. |
| `voltage` | `CxVec` | Complex bus voltage solution. |
| `Sf` | `CxVec` | Complex branch power flow at the "from" side. |
| `St` | `CxVec` | Complex branch power flow at the "to" side. |
| `If` | `CxVec` | Complex branch current at the "from" side. |
| `It` | `CxVec` | Complex branch current at the "to" side. |
| `tap_module` | `Vec` | Transformer tap module used in the solved state. |
| `tap_angle` | `Vec` | Transformer tap angle used in the solved state. |
| `Vbranch` | `CxVec` | Complex branch voltage magnitude reference used by branch reports. |
| `loading` | `CxVec` | Branch loading result. |
| `losses` | `CxVec` | Complex branch losses. |
| `losses_hvdc` | `Vec` | HVDC active power losses. |
| `Pf_hvdc` | `Vec` | HVDC active power flow at the "from" side. |
| `Pt_hvdc` | `Vec` | HVDC active power flow at the "to" side. |
| `loading_hvdc` | `Vec` | HVDC loading result. |
| `losses_vsc` | `Vec` | VSC active power losses. |
| `Pfp_vsc` | `Vec` | VSC positive-pole active power flow at the "from" side. |
| `Pfn_vsc` | `Vec` | VSC negative-pole active power flow at the "from" side. |
| `St_vsc` | `CxVec` | Complex VSC power flow at the AC terminal. |
| `If_vsc` | `Vec` | VSC current at the "from" side. |
| `It_vsc` | `CxVec` | Complex VSC current at the AC terminal. |
| `loading_vsc` | `Vec` | VSC loading result. |
| `gen_q` | `Vec` | Generator reactive power output. |
| `battery_q` | `Vec` | Battery reactive power output. |
| `shunt_q` | `Vec` | Shunt reactive power output. |
| `three_phase` | `bool` | Flag indicating whether these results belong to a three-phase power flow. |
Output:
```text
IEEE14_from_raw
Converged: True error: 5.98e-08
Bus resuts:
Vm Va P Q
BUS 1 1.06 0.00 232.39 -16.55
BUS 2 1.04 -4.98 18.30 30.86
BUS 3 1.01 -12.73 -94.20 6.08
BUS 4 1.02 -10.31 -47.80 3.90
BUS 5 1.02 -8.77 -7.60 -1.60
BUS 6 1.07 -14.22 -11.20 5.23
BUS 7 1.06 -13.36 0.00 0.00
BUS 8 1.09 -13.36 0.00 17.62
BUS 9 1.06 -14.94 -29.50 -16.60
BUS 10 1.05 -15.10 -9.00 -5.80
BUS 11 1.06 -14.79 -3.50 -1.80
BUS 12 1.06 -15.08 -6.10 -1.60
BUS 13 1.05 -15.16 -13.50 -5.80
BUS 14 1.04 -16.03 -14.90 -5.00
Branch results:
Pf Qf Pt Qt loading Ploss Qloss
1_2_1 156.882887 -20.404291 -152.585286 27.676248 15688288652036.908203 4.297600 7.271957
1_5_1 75.510380 3.854989 -72.747507 2.229360 7551037982438.064453 2.762872 6.084349
2_3_1 73.237578 3.560203 -70.914309 1.602232 7323757808601.912109 2.323269 5.162436
2_4_1 56.131495 -1.550352 -54.454837 3.020689 5613149456668.273438 1.676658 1.470337
2_5_1 41.516214 1.170996 -40.612460 -2.099032 4151621353697.657715 0.903753 -0.928036
3_4_1 -23.285691 4.473114 23.659136 -4.835650 -2328569062725.765625 0.373445 -0.362537
4_5_1 -61.158231 15.823642 61.672651 -14.201004 -6115823108351.800781 0.514420 1.622637
6_11_1 7.353277 3.560471 -7.297904 -3.444512 735327693069.753418 0.055373 0.115959
6_12_1 7.786067 2.503414 -7.714258 -2.353959 778606687855.751465 0.071809 0.149455
6_13_1 17.747977 7.216574 -17.535891 -6.798912 1774797671583.112793 0.212085 0.417662
7_8_1 -0.000000 -17.162967 0.000000 17.623448 -0.001718 0.000000 0.460481
7_9_1 28.074179 5.778690 -28.074179 -4.976621 2807417855964.891602 0.000000 0.802069
9_10_1 5.227551 4.219139 -5.214676 -4.184938 522755058212.680359 0.012875 0.034201
9_14_1 9.426380 3.610007 -9.310226 -3.362932 942638030136.208130 0.116154 0.247075
10_11_1 -3.785324 -1.615061 3.797906 1.644513 -378532426869.186707 0.012581 0.029451
12_13_1 1.614258 0.753959 -1.607959 -0.748260 161425771970.211853 0.006298 0.005698
13_14_1 5.643852 1.747172 -5.589774 -1.637068 564385175482.526855 0.054078 0.110105
4_7_1 28.074176 -9.681066 -28.074176 11.384281 2807417645485.176270 0.000000 1.703214
4_9_1 16.079758 -0.427611 -16.079758 1.732322 1607975830176.256104 0.000000 1.304711
5_6_1 44.087319 12.470682 -44.087319 -8.049520 4408731875605.579102 0.000000 4.421161
```
Getting the Jacobian:
```python
import os
import VeraGridEngine as vg
from VeraGridEngine.Simulations.Derivatives.ac_jacobian import AC_jacobian
folder = os.path.join('..', 'Grids_and_profiles', 'grids')
fname = os.path.join(folder, 'IEEE39_1W.veragrid')
main_circuit = vg.open_file(fname)
nc = vg.compile_numerical_circuit_at(circuit=main_circuit)
adm = nc.get_admittance_matrices()
sim = nc.get_simulation_indices()
J = AC_jacobian(
Ybus=adm.Ybus,
V=nc.bus_data.Vbus, # The voltage initla guess, substitute by any other compatible voltage
pvpq=sim.no_slack,
pq=sim.pq
)
print(J.toarray())
```
## Benchmark
### Matpower grids
Matpower's excellent formulations and consistency has allowed this and other
projects to develop, relying on its sound math. That is why VeraGrid reads Matpower
cases out of the box, without you having to do anything special.
And of course, VeraGrid solves all Matpower 8 provided grids,
solving the continental USA case in about 1 second:
| name | n_buses | n_branches | P imbalance (%) | Flat start | converged | error (p.u.) | iterations | time (ms) |
|-------------------------|---------|------------|-----------------|------------|-----------|--------------|------------|-----------|
| case_SyntheticUSA.m | 82000 | 104121 | -0.12 | FALSE | TRUE | 2.03E-08 | 13 | 3181.41 |
| case_ACTIVSg70k.m | 70000 | 88207 | 0.64 | FALSE | TRUE | 8.00E-07 | 4 | 1170.59 |
| case_ACTIVSg25k.m | 25000 | 32230 | -2.72 | TRUE | TRUE | 7.77E-10 | 13 | 993.16 |
| case13659pegase.m | 13659 | 20467 | 2411.01 | TRUE | TRUE | 1.66E-07 | 8 | 284.45 |
| case_ACTIVSg10k.m | 10000 | 12706 | -7.61 | FALSE | TRUE | 3.81E-11 | 5 | 152.66 |
| case9241pegase.m | 9241 | 16049 | 683.53 | TRUE | TRUE | 9.97E-11 | 11 | 319.74 |
| case8387pegase.m | 8387 | 14561 | -44.18 | TRUE | TRUE | 1.27E-11 | 15 | 362.14 |
| case6515rte.m | 6515 | 9037 | -47.51 | FALSE | TRUE | 9.36E-08 | 7 | 132.49 |
| case6495rte.m | 6495 | 9019 | -48.91 | FALSE | TRUE | 1.08E-07 | 6 | 111.79 |
| case6470rte.m | 6470 | 9005 | -47.97 | FALSE | TRUE | 7.89E-09 | 7 | 130.71 |
| case6468rte.m | 6468 | 9000 | -46.27 | FALSE | TRUE | 1.20E-09 | 7 | 145.08 |
| case3375wp.m | 3374 | 4161 | -73.29 | FALSE | TRUE | 2.26E-09 | 5 | 52.88 |
| case3120sp.m | 3120 | 3693 | -99.96 | TRUE | TRUE | 1.56E-08 | 12 | 142.90 |
| case3012wp.m | 3012 | 3572 | -98.93 | FALSE | TRUE | 2.85E-10 | 6 | 67.21 |
| case2869pegase.m | 2869 | 4582 | 561.42 | TRUE | TRUE | 1.22E-11 | 10 | 78.79 |
| case2868rte.m | 2868 | 3808 | -46.03 | TRUE | TRUE | 5.05E-07 | 15 | 119.37 |
| case2848rte.m | 2848 | 3776 | -41.34 | TRUE | TRUE | 1.30E-11 | 18 | 153.64 |
| case2746wp.m | 2746 | 3514 | -95.83 | TRUE | TRUE | 9.13E-09 | 9 | 72.19 |
| case2746wop.m | 2746 | 3514 | -96.55 | TRUE | TRUE | 2.78E-07 | 9 | 128.23 |
| case2737sop.m | 2737 | 3506 | -94.20 | TRUE | TRUE | 3.01E-10 | 11 | 94.70 |
| case2736sp.m | 2736 | 3504 | -95.17 | TRUE | TRUE | 1.74E-09 | 10 | 114.32 |
| case2383wp.m | 2383 | 2896 | -97.43 | TRUE | TRUE | 5.17E-12 | 10 | 87.65 |
| case_ACTIVSg2000.m | 2000 | 3206 | 10.84 | TRUE | TRUE | 6.16E-10 | 12 | 136.57 |
| case1951rte.m | 1951 | 2596 | -45.97 | FALSE | TRUE | 6.67E-07 | 4 | 31.74 |
| case1888rte.m | 1888 | 2531 | -47.08 | FALSE | TRUE | 1.56E-07 | 5 | 39.21 |
| case1354pegase.m | 1354 | 1991 | 862.85 | TRUE | TRUE | 7.43E-09 | 8 | 56.87 |
| case_ACTIVSg500.m | 500 | 597 | 2.78 | TRUE | TRUE | 4.17E-09 | 7 | 19.51 |
| case300.m | 300 | 411 | -38.64 | TRUE | TRUE | 1.54E-09 | 8 | 7.69 |
| case_ACTIVSg200.m | 200 | 245 | 6.48 | TRUE | TRUE | 3.07E-10 | 5 | 11.00 |
| case145.m | 145 | 453 | -100.00 | TRUE | TRUE | 1.04E-09 | 8 | 14.71 |
| case141.m | 141 | 140 | -100.00 | TRUE | TRUE | 4.69E-09 | 2 | 1.15 |
| case136ma.m | 136 | 156 | -100.00 | TRUE | TRUE | 1.14E-08 | 2 | 1.11 |
| case118zh.m | 118 | 132 | -100.00 | TRUE | TRUE | 1.46E-08 | 2 | 1.35 |
| case118.m | 118 | 186 | -28.33 | TRUE | TRUE | 1.94E-07 | 7 | 3.05 |
| case94pi.m | 94 | 93 | -100.00 | TRUE | TRUE | 2.08E-11 | 2 | 1.32 |
| case89pegase.m | 89 | 210 | 4.19 | TRUE | TRUE | 2.81E-09 | 4 | 2.54 |
| case85.m | 85 | 84 | -100.00 | TRUE | TRUE | 7.90E-12 | 2 | 1.13 |
| case74ds.m | 74 | 73 | -100.00 | TRUE | TRUE | 8.74E-07 | 1 | 0.58 |
| case_RTS_GMLC.m | 73 | 120 | -80.55 | TRUE | TRUE | 1.62E-07 | 7 | 11.74 |
| case70da.m | 70 | 76 | -100.00 | TRUE | TRUE | 2.17E-12 | 2 | 0.89 |
| case69.m | 69 | 68 | -100.00 | TRUE | TRUE | 7.20E-09 | 2 | 0.78 |
| case60nordic.m | 60 | 88 | 96.26 | TRUE | TRUE | 5.15E-08 | 4 | 1.35 |
| case57.m | 57 | 80 | -100.00 | TRUE | TRUE | 2.82E-10 | 7 | 11.47 |
| case51ga.m | 51 | 50 | -100.00 | TRUE | TRUE | 1.85E-12 | 2 | 0.86 |
| case51he.m | 51 | 50 | -100.00 | TRUE | TRUE | 6.16E-07 | 1 | 0.40 |
| case39.m | 39 | 46 | -26.13 | TRUE | TRUE | 1.93E-11 | 7 | 10.79 |
| case38si.m | 38 | 37 | -100.00 | TRUE | TRUE | 7.26E-12 | 2 | 1.06 |
| case34sa.m | 34 | 33 | -100.00 | TRUE | TRUE | 8.24E-13 | 2 | 0.74 |
| case33bw.m | 33 | 37 | -100.00 | TRUE | TRUE | 7.38E-09 | 2 | 1.09 |
| case33mg.m | 33 | 37 | -100.00 | TRUE | TRUE | 7.46E-12 | 2 | 0.70 |
| case30.m | 30 | 41 | -39.59 | TRUE | TRUE | 9.57E-10 | 3 | 1.13 |
| case_ieee30.m | 30 | 41 | -3.24 | TRUE | TRUE | 5.18E-08 | 3 | 0.92 |
| case30pwl.m | 30 | 41 | -39.59 | TRUE | TRUE | 9.57E-10 | 3 | 0.89 |
| case30Q.m | 30 | 41 | -39.59 | TRUE | TRUE | 9.57E-10 | 3 | 0.92 |
| case28da.m | 28 | 27 | -100.00 | TRUE | TRUE | 6.85E-07 | 1 | 0.53 |
| case24_ieee_rts.m | 24 | 38 | -70.52 | TRUE | TRUE | 1.63E-08 | 5 | 9.04 |
| case22.m | 22 | 21 | -100.00 | TRUE | TRUE | 2.13E-07 | 1 | 0.47 |
| case18.m | 18 | 17 | -100.00 | TRUE | TRUE | 1.27E-08 | 3 | 0.74 |
| case18nbr.m | 18 | 17 | -100.00 | TRUE | TRUE | 1.35E-07 | 2 | 0.84 |
| case17me.m | 17 | 16 | -100.00 | TRUE | TRUE | 3.19E-08 | 3 | 0.91 |
| case16ci.m | 16 | 16 | -100.00 | TRUE | TRUE | 1.38E-09 | 2 | 0.47 |
| case15da.m | 15 | 14 | -100.00 | TRUE | TRUE | 8.57E-07 | 1 | 0.68 |
| case15nbr.m | 15 | 14 | -100.00 | TRUE | TRUE | 5.43E-08 | 2 | 0.63 |
| case16am.m | 15 | 14 | -100.00 | TRUE | TRUE | 1.22E-07 | 2 | 0.89 |
| case14.m | 14 | 20 | 1.26 | TRUE | TRUE | 5.98E-08 | 3 | 0.61 |
| case12da.m | 12 | 11 | -100.00 | TRUE | TRUE | 2.71E-07 | 1 | 1.36 |
| case10ba.m | 10 | 9 | -100.00 | TRUE | TRUE | 6.22E-08 | 2 | 0.52 |
| case9target.m | 9 | 9 | -41.08 | TRUE | TRUE | 1.40E-07 | 5 | 2.12 |
| case9_gurobi_test.m | 9 | 9 | 1.68 | TRUE | TRUE | 3.42E-07 | 3 | 1.61 |
| case9.m | 9 | 9 | 1.68 | TRUE | TRUE | 3.42E-07 | 3 | 0.67 |
| case9Q.m | 9 | 9 | -21.27 | TRUE | TRUE | 5.71E-07 | 3 | 0.64 |
| case6ww.m | 6 | 11 | -47.62 | TRUE | TRUE | 2.09E-10 | 3 | 1.43 |
| case5.m | 5 | 6 | -36.35 | TRUE | TRUE | 6.42E-11 | 3 | 1.87 |
| case4_dist.m | 4 | 3 | -100.00 | TRUE | TRUE | 4.63E-11 | 6 | 8.03 |
| case4gs.m | 4 | 4 | -100.00 | TRUE | TRUE | 6.59E-14 | 6 | 1.14 |
| contab_ACTIVSg200.m | 0 | 0 | 0.00 | TRUE | TRUE | 0.00E+00 | 0 | 0.00 |
| contab_ACTIVSg500.m | 0 | 0 | 0.00 | TRUE | TRUE | 0.00E+00 | 0 | 0.00 |
| scenarios_ACTIVSg200.m | 0 | 0 | 0.00 | TRUE | TRUE | 0.00E+00 | 0 | 0.00 |
| contab_ACTIVSg10k.m | 0 | 0 | 0.00 | TRUE | TRUE | 0.00E+00 | 0 | 0.00 |
| contab_ACTIVSg2000.m | 0 | 0 | 0.00 | TRUE | TRUE | 0.00E+00 | 0 | 0.00 |
| scenarios_ACTIVSg2000.m | 0 | 0 | 0.00 | TRUE | TRUE | 0.00E+00 | 0 | 0.00 |
_Results simulated with AMD 9750x and 64 GB of RAM under Ubuntu 24.04.
All solved using Newton-Raphson, and only using the provided solution
that comes with the files when the flat start fails.
Cool right? This is the code to reproduce the results:
```python
import os
import pandas as pd
import multiprocessing as mp
import VeraGridEngine as gce
folder = "[some path...]/matpower8.0b1/data"
def run_grid(fname):
grid = gce.open_file(fname)
name = os.path.basename(fname)
if grid.get_bus_number() > 0:
res = gce.power_flow(
grid=grid,
options=gce.PowerFlowOptions(solver_type=gce.SolverType.NR,
retry_with_other_methods=False,
use_stored_guess=False)
)
flat_start = True
if not res.converged:
# if it does not converge, retry with the provided solution
res = gce.power_flow(
grid=grid,
options=gce.PowerFlowOptions(solver_type=gce.SolverType.NR,
retry_with_other_methods=False,
use_stored_guess=True)
)
flat_start = False
info = {
"name": name,
"n_buses": grid.get_bus_number(),
"n_branches": grid.get_branch_number(),
"P imbalance (%)": grid.get_imbalance() * 100.0,
"Flat start": flat_start,
"converged": res.converged,
"error (p.u.)": res.error,
"iterations": res.iterations,
"time (ms)": res.elapsed * 1000.0,
}
else:
info = {
"name": name,
"n_buses": grid.get_bus_number(),
"n_branches": grid.get_branch_number(),
"P imbalance (%)": 0.0,
"Flat start": True,
"converged": True,
"error (p.u.)": 0,
"iterations": 0,
"time (ms)": 0,
}
return info
# run this one to compile all JIT routines and have fair measurements
gce.power_flow(gce.open_file(os.path.join(folder, "case_ieee30.m")))
data = list()
files_list = list()
for root, dirs, files in os.walk(folder):
for file in files:
if file.endswith(".m"):
path = os.path.join(root, file)
files_list.append(path)
with mp.Pool(mp.cpu_count()) as p:
data = p.map(run_grid, files_list)
df = pd.DataFrame(data).sort_values(by='n_buses', ascending=False)
df.to_excel("All matpower grids.xlsx", index=False)
```
## Three-phase Unbalanced Power Flow
The power flow problem, or load flow problem, consists of computing the flow of electrical power in a steady-state
system. In practice, this amounts to calculating the voltage phasor in each bus of the power network.
Once the bus voltages are known, the other electrical variables are easy to compute. Mathematical equations for the
power flow problem can be obtained by combining the complex power equation with Ohm's law.
The injected complex power $\vec{S_i}$ can be written as the product of the voltage phasor $\vec{U_i}$ and the
conjugated current phasor $\vec{I_i}^*$ in the bus $i$.
$$
\vec{S_i} = \vec{U_i} \vec{I_i}^*
$$
Ohm's law relates the current $\vec{I_i}$ flowing through the bus $i$ with the bus voltage $\vec{U_i}$ and its
impedance $\vec{Z_i}$. We can express it in terms of admittance $\vec{Y_i}$, which is the impedance inverse.
$$
\vec{I_i} = \frac{\vec{U_i}}{\vec{Z_i}} = \vec{Y_i}\vec{U_i}
$$
Therefore, by substituting the current from the first equation into the second equation, the fundamental expression to
be solved is obtained as:
$$
\vec{S_i} = \vec{U_i}(\vec{Y_i}\vec{U_i})^*
$$
The following magnitudes are defined for each bus $i$:
- $U$: Voltage module.
- $\delta$: Voltage angle.
- $P$: Active power injection or consumption.
- $Q$: Reactive power injection or consumption.
In some buses we artificially define the voltage module and angle (slack buses), in some buses we just know the
consumption power (load buses or PQ), and in some others we know the active power injection and the voltage module
(generation buses or PV). Thus, the different bus types are summarised in the following table.
At least one slack bus is required in order to solve the power flow problem.
| | $U$ | $\delta$ | $P$ | $Q$ |
|-------|----------|----------|----------|----------|
| Slack | Set | Set | Computed | Computed |
| PQ | Computed | Computed | Set | Set |
| PV | Set | Computed | Set | Computed |
### The Newton-Raphson Method
The Newton-Raphson method is used to solve the power flow problem, which consists of a zero-finding algorithm that
produces successively better approximations to the zeros of a function. A zero of a function $f$ is the $x$ value such
that $f(x)=0$. The idea is to start with an initial guess $x_0$, then to approximate the function by its tangent line,
and finally to compute the x-axis interception of this tangent line.
As shown in the following figure, the $x_1$ interception will typically be a better approximation to the
function's real zero $x_n$ than the initial guess $x_0$, and the method can be iterated.
If a function $y=f(x)$ is expanded about an operating point $x_0$ using a Taylor series expansion and ignoring all
terms of higher order than 1, the following equation can be obtained:
$$
y = f(x_0) + \frac{df(x_0)}{dx} (x - x_0)
$$
Returning to the power flow problem and considering again the fundamental equation $\vec{S_i} = \vec{U_i}(\vec{Y_i}\vec{U_i})^*$
to be solved, it can be noted that it may be expressed as $\vec{S}_i = f(\vec{U}_i)$, in the same way as a generic function $y = f(x)$.
Therefore, it can be solved using the Newton-Raphson method.
### Solving the Three-phase Unbalanced Power Flow Problem
From the previous consideration, the solution of the function $\vec{S_i} = f(\vec{U_i})$ can be written in the form of
the expression above, the fundamental Newton-Raphson expression for the solution of the power flow problem,
where $J$ is the Jacobian matrix:
$$
\varDelta \vec{S} = -J \cdot \vec{U}
$$
If the previous equation is developed such that the complex power is expressed in its respective active $P$ and
reactive $Q$ components, $J$ is the Jacobian matrix composed by the partial derivatives, and the voltage is isolated
and expressed in polar coordinates with angle $\delta$ and module $U$. Then, its matrix full representation is as follows:
$$
\begin{bmatrix}
\delta \\
U
\end{bmatrix}
=-
\begin{bmatrix}
\dfrac{\partial P}{\partial \delta} & \dfrac{\partial P}{\partial U} \\ \\
\dfrac{\partial Q}{\partial \delta} & \dfrac{\partial Q}{\partial U}
\end{bmatrix}^{-1}
\begin{bmatrix}
\varDelta P \\
\varDelta Q
\end{bmatrix}
$$
Therefore, the system will be solved iteratively until the difference between the specified power $S_\text{spec}$ and
the calculated power $S_\text{calc}$ is negligible, meaning that the zero of the function has been found:
$$
\varDelta \vec{S} = \vec{S}_\text{calc} - \vec{S}_\text{spec}
$$
The calculated power can be obtained using the voltage and the admittance matrices of the lines and
transformers, as previously developed:
$$
\vec{S}_\text{calc} = \vec{U} \, (\vec{Y} \, \vec{U})^*
$$
In contrast, the specified power incorporates the injection values of the different elements defined by
the user. As previously discussed in the modelling section, these elements may be defined in terms of specified power
$\vec{S}_0$, current $\vec{I}_0$, or admittance $\vec{Y}_0$. Then, the equivalent specified power is calculated using
the following expression:
$$
\vec{S}_\text{spec} = \vec{S}_0 + (\vec{I}_0 + \vec{Y}_0 \, \vec{U})^* \, \vec{U}
$$
The voltage term $\vec{U}$ is a vector composed by the voltage values at each of the system's buses, and therefore has
a size equal to $n$, that is, the total number of buses. The admittance term $\vec{Y}$ is a matrix comprising the
admittance of the lines and transformers, therefore the elements that interconnect the buses with one another.
As a result, it has a size of $n \times n$.
On the other hand, the specified power and current injection vectors, $\vec{S}_0$ and $\vec{I}_0$, have the same length
as the voltage vector $\vec{U}$, while the specified admittance $\vec{Y}_0$ has the same dimensions as the branch
admittance matrix $\vec{Y}$.
One might assume that, in order to formulate a three-phase power flow, it would be sufficient to simply multiply the
number of buses by three, to account for the $abc$ phases. While this assumption could be valid in transmission systems,
it does not hold true when we enter the domain of distribution networks. In such systems, all kinds of elements can be
present, with three-phase, two-phase, and single-phase components.
Since VeraGrid aims to address both transmission and distribution systems, the assumption of simply tripling the number
of buses is not sufficient. A more advanced method must be implemented to accurately detect the real number of buses,
taking into account the presence or absence of the three phases.
This system for detecting the existing phases at each bus must begin with the branches (lines and transformers).
These are the elements where the presence of phases transporting electrical energy is specified. Based on this
information, a binary vector, referred to as the bus **mask**, is constructed. This vector contains ones for the phases
that are present and zeros for those positions where one of the three standard phases is absent.
This mask $M$ must be applied to several variables. The first is the admittance matrix, which will be processed by
multiplying both its rows and columns by the bus mask. In this way, the entries corresponding to non-existent phases
are removed. As a result, the matrix remains of size equal to $n \times n$, that this time corresponds to the total
number of **active buses**.
$$
\vec{Y}_{\text{masked}} = M \cdot \vec{Y} \cdot M^{\top}
$$
Considering a simple example of two buses connected through a line that only has phases $a$ and $c$ active, phase $b$
cannot exist at the second bus. Therefore, the bus mask will contain a zero in that position, effectively removing the
corresponding rows and columns from the admittance matrix:
Naturally, the bus mask must also be applied to the voltage vector at each bus in order to eliminate the entries
corresponding to non-existent buses:
$$
\vec{U}_{\text{masked}} = M^{\top} \cdot \vec{U}
$$
And the same operation must be carried out equivalently for the current injection vector and the power vector
at each bus:
$$
\vec{I}_{\text{masked}} = M^{\top} \cdot \vec{I}
$$
$$
\vec{S}_{\text{masked}} = M^{\top} \cdot \vec{S}
$$
In this way, all the elements required for the power flow are now correctly dimensioned, and the algorithm can be
applied in a manner equivalent to the traditional single-phase power flow.
### Benchmark - IEEE 13 Node Test Feeder
To validate the implementation in VeraGrid of the three-phase power flow algorithm, the IEEE 13 Node Test Feeder has
been selected as a reference benchmark. This distribution system model, published by the IEEE Distribution System
Analysis Subcommittee, represents a typical North American urban radial distribution network.
The network, depicted in the figure above, includes a variety of distribution elements such as:
- Unbalanced overhead power lines with different phasing configurations (three-phase, two-phase, and single-phase).
- A distribution transformer between buses 633 and 634.
- Multiple balanced and unbalanced loads modelled as constant impedance (Z), constant power (P), and constant current (I), with three-phase connections in star (Y) and delta (D), connections between two phases, and also between one phase and the ground.
- Shunt capacitor banks for reactive power compensation, also with different phasing configurations (three-phase and single-phase).
The base voltage levels of the system are 4,16 kV on the distribution side, with a
downstream low-voltage level of 0,48 kV at bus 634. The total nominal load is around 3,8 MVA, distributed across
residential, commercial, and industrial consumers. Due to the network's unbalanced and mixed-phase nature, this test feeder is widely recognised as a challenging yet
representative case for validating unbalanced power flow solvers.
The following code can be used to model the IEEE 13 Node Test Feeder in VeraGrid and to perform the three-phase
unbalanced power flow, obtaining the voltage results depicted bellow.
```python
import VeraGridEngine.api as gce
from VeraGridEngine import WindingType, ShuntConnectionType
import numpy as np
logger = gce.Logger()
grid = gce.MultiCircuit()
grid.fBase = 60
# ----------------------------------------------------------------------------------------------------------------------
# Buses
# ----------------------------------------------------------------------------------------------------------------------
bus_632 = gce.Bus(name='632', Vnom=4.16, xpos=0, ypos=0)
bus_632.is_slack = True
grid.add_bus(obj=bus_632)
gen = gce.Generator(vset=1.0)
grid.add_generator(bus=bus_632, api_obj=gen)
bus_633 = gce.Bus(name='633', Vnom=4.16, xpos=100 * 5, ypos=0)
grid.add_bus(obj=bus_633)
bus_634 = gce.Bus(name='634', Vnom=0.48, xpos=200 * 5, ypos=0)
grid.add_bus(obj=bus_634)
bus_645 = gce.Bus(name='645', Vnom=4.16, xpos=-100 * 5, ypos=0)
grid.add_bus(obj=bus_645)
bus_646 = gce.Bus(name='646', Vnom=4.16, xpos=-200 * 5, ypos=0)
grid.add_bus(obj=bus_646)
bus_652 = gce.Bus(name='652', Vnom=4.16, xpos=-100 * 5, ypos=200 * 5)
grid.add_bus(obj=bus_652)
bus_671 = gce.Bus(name='671', Vnom=4.16, xpos=0, ypos=100 * 5)
grid.add_bus(obj=bus_671)
bus_675 = gce.Bus(name='675', Vnom=4.16, xpos=200 * 5, ypos=100 * 5)
grid.add_bus(obj=bus_675)
bus_611 = gce.Bus(name='611', Vnom=4.16, xpos=-200 * 5, ypos=100 * 5)
grid.add_bus(obj=bus_611)
bus_680 = gce.Bus(name='680', Vnom=4.16, xpos=0, ypos=200 * 5)
grid.add_bus(obj=bus_680)
bus_684 = gce.Bus(name='684', Vnom=4.16, xpos=-100 * 5, ypos=100 * 5)
grid.add_bus(obj=bus_684)
# ----------------------------------------------------------------------------------------------------------------------
# Impedances [Ohm/km]
# ----------------------------------------------------------------------------------------------------------------------
z_601 = np.array([
[0.3465 + 1j * 1.0179, 0.1560 + 1j * 0.5017, 0.1580 + 1j * 0.4236],
[0.1560 + 1j * 0.5017, 0.3375 + 1j * 1.0478, 0.1535 + 1j * 0.3849],
[0.1580 + 1j * 0.4236, 0.1535 + 1j * 0.3849, 0.3414 + 1j * 1.0348]
], dtype=complex) / 1.60934
z_602 = np.array([
[0.7526 + 1j * 1.1814, 0.1580 + 1j * 0.4236, 0.1560 + 1j * 0.5017],
[0.1580 + 1j * 0.4236, 0.7475 + 1j * 1.1983, 0.1535 + 1j * 0.3849],
[0.1560 + 1j * 0.5017, 0.1535 + 1j * 0.3849, 0.7436 + 1j * 1.2112]
], dtype=complex) / 1.60934
z_603 = np.array([
[1.3294 + 1j * 1.3471, 0.2066 + 1j * 0.4591],
[0.2066 + 1j * 0.4591, 1.3238 + 1j * 1.3569]
], dtype=complex) / 1.60934
z_604 = np.array([
[1.3238 + 1j * 1.3569, 0.2066 + 1j * 0.4591],
[0.2066 + 1j * 0.4591, 1.3294 + 1j * 1.3471]
], dtype=complex) / 1.60934
z_605 = np.array([
[1.3292 + 1j * 1.3475]
], dtype=complex) / 1.60934
z_606 = np.array([
[0.7982 + 1j * 0.4463, 0.3192 + 1j * 0.0328, 0.2849 + 1j * -0.0143],
[0.3192 + 1j * 0.0328, 0.7891 + 1j * 0.4041, 0.3192 + 1j * 0.0328],
[0.2849 + 1j * -0.0143, 0.3192 + 1j * 0.0328, 0.7982 + 1j * 0.4463]
], dtype=complex) / 1.60934
z_607 = np.array([
[1.3425 + 1j * 0.5124]
], dtype=complex) / 1.60934
# ----------------------------------------------------------------------------------------------------------------------
# Admittances [S/km]
# ----------------------------------------------------------------------------------------------------------------------
y_601 = np.array([
[1j * 6.2998, 1j * -1.9958, 1j * -1.2595],
[1j * -1.9958, 1j * 5.9597, 1j * -0.7417],
[1j * -1.2595, 1j * -0.7417, 1j * 5.6386]
], dtype=complex) / 10 ** 6 / 1.60934
y_602 = np.array([
[1j * 5.6990, 1j * -1.0817, 1j * -1.6905],
[1j * -1.0817, 1j * 5.1795, 1j * -0.6588],
[1j * -1.6905, 1j * -0.6588, 1j * 5.4246]
], dtype=complex) / 10 ** 6 / 1.60934
y_603 = np.array([
[1j * 4.7097, 1j * -0.8999],
[1j * -0.8999, 1j * 4.6658]
], dtype=complex) / 10 ** 6 / 1.60934
y_604 = np.array([
[1j * 4.6658, 1j * -0.8999],
[1j * -0.8999, 1j * 4.7097]
], dtype=complex) / 10 ** 6 / 1.60934
y_605 = np.array([
[1j * 4.5193]
], dtype=complex) / 10 ** 6 / 1.60934
y_606 = np.array([
[1j * 96.8897, 1j * 0.0000, 1j * 0.0000],
[1j * 0.0000, 1j * 96.8897, 1j * 0.0000],
[1j * 0.0000, 1j * 0.0000, 1j * 96.8897]
], dtype=complex) / 10 ** 6 / 1.60934
y_607 = np.array([
[1j * 88.9912]
], dtype=complex) / 10 ** 6 / 1.60934
# ----------------------------------------------------------------------------------------------------------------------
# Loads
# ----------------------------------------------------------------------------------------------------------------------
load_634 = gce.Load(P1=0.160,
Q1=0.110,
P2=0.120,
Q2=0.090,
P3=0.120,
Q3=0.090)
load_634.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_634, api_obj=load_634)
load_645 = gce.Load(P1=0.0,
Q1=0.0,
P2=0.170,
Q2=0.125,
P3=0.0,
Q3=0.0)
load_645.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_645, api_obj=load_645)
load_646 = gce.Load(G1=0.0,
B1=0.0,
G2=0.230,
B2=-0.132,
G3=0.0,
B3=0.0)
load_646.conn = ShuntConnectionType.Delta
grid.add_load(bus=bus_646, api_obj=load_646)
load_652 = gce.Load(G1=0.128,
B1=-0.086,
G2=0.0,
B2=0.0,
G3=0.0,
B3=0.0)
load_652.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_652, api_obj=load_652)
load_671 = gce.Load(P1=0.385,
Q1=0.220,
P2=0.385,
Q2=0.220,
P3=0.385,
Q3=0.220)
load_671.conn = ShuntConnectionType.Delta
grid.add_load(bus=bus_671, api_obj=load_671)
load_671_692 = gce.Load(Ir1=0.0,
Ii1=0.0,
Ir2=0.0,
Ii2=0.0,
Ir3=0.170,
Ii3=0.151)
load_671_692.conn = ShuntConnectionType.Delta
grid.add_load(bus=bus_671, api_obj=load_671_692)
load_611 = gce.Load(Ir1=0.0,
Ii1=0.0,
Ir2=0.0,
Ii2=0.0,
Ir3=0.170,
Ii3=0.080)
load_611.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_611, api_obj=load_611)
load_632_distrib = gce.Load(P1=0.017 / 2,
Q1=0.010 / 2,
P2=0.066 / 2,
Q2=0.038 / 2,
P3=0.117 / 2,
Q3=0.068 / 2)
load_632_distrib.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_632, api_obj=load_632_distrib)
load_671_distrib = gce.Load(P1=0.017 / 2,
Q1=0.010 / 2,
P2=0.066 / 2,
Q2=0.038 / 2,
P3=0.117 / 2,
Q3=0.068 / 2)
load_671_distrib.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_671, api_obj=load_671_distrib)
load_675 = gce.Load(P1=0.485,
Q1=0.190,
P2=0.068,
Q2=0.060,
P3=0.290,
Q3=0.212)
load_675.conn = ShuntConnectionType.GroundedStar
grid.add_load(bus=bus_675, api_obj=load_675)
# ----------------------------------------------------------------------------------------------------------------------
# Capacitors
# ----------------------------------------------------------------------------------------------------------------------
cap_675 = gce.Shunt(B1=0.2,
B2=0.2,
B3=0.2)
cap_675.conn = ShuntConnectionType.GroundedStar
grid.add_shunt(bus=bus_675, api_obj=cap_675)
cap_611 = gce.Shunt(B1=0.0,
B2=0.0,
B3=0.1)
cap_611.conn = ShuntConnectionType.GroundedStar
grid.add_shunt(bus=bus_611, api_obj=cap_611)
# ----------------------------------------------------------------------------------------------------------------------
# Line Configurations
# ----------------------------------------------------------------------------------------------------------------------
config_601 = gce.create_known_abc_overhead_template(name='Config. 601',
z_abc=z_601,
ysh_abc=y_601,
phases=np.array([1, 2, 3]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_601)
config_602 = gce.create_known_abc_overhead_template(name='Config. 602',
z_abc=z_602,
ysh_abc=y_602,
phases=np.array([1, 2, 3]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_602)
config_603 = gce.create_known_abc_overhead_template(name='Config. 603',
z_abc=z_603,
ysh_abc=y_603,
phases=np.array([2, 3]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_603)
config_604 = gce.create_known_abc_overhead_template(name='Config. 604',
z_abc=z_604,
ysh_abc=y_604,
phases=np.array([1, 3]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_604)
config_605 = gce.create_known_abc_overhead_template(name='Config. 605',
z_abc=z_605,
ysh_abc=y_605,
phases=np.array([3]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_605)
config_606 = gce.create_known_abc_overhead_template(name='Config. 606',
z_abc=z_606,
ysh_abc=y_606,
phases=np.array([1, 2, 3]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_606)
config_607 = gce.create_known_abc_overhead_template(name='Config. 607',
z_abc=z_607,
ysh_abc=y_607,
phases=np.array([1]),
Vnom=4.16,
frequency=60)
grid.add_overhead_line(config_607)
# ----------------------------------------------------------------------------------------------------------------------
# Lines and Transformers
# ----------------------------------------------------------------------------------------------------------------------
line_632_645 = gce.Line(bus_from=bus_632,
bus_to=bus_645,
length=500 * 0.0003048)
line_632_645.apply_template(config_603, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_632_645)
line_645_646 = gce.Line(bus_from=bus_645,
bus_to=bus_646,
length=300 * 0.0003048)
line_645_646.apply_template(config_603, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_645_646)
line_632_633 = gce.Line(bus_from=bus_632,
bus_to=bus_633,
length=500 * 0.0003048)
line_632_633.apply_template(config_602, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_632_633)
XFM_1 = gce.Transformer2W(name='XFM-1',
bus_from=bus_633,
bus_to=bus_634,
HV=4.16,
LV=0.48,
nominal_power=0.5,
rate=0.5,
r=1.1 * 2,
x=2 * 2)
XFM_1.conn_f = WindingType.GroundedStar
XFM_1.conn_t = WindingType.GroundedStar
grid.add_transformer2w(XFM_1)
line_632_671 = gce.Line(bus_from=bus_632,
bus_to=bus_671,
length=2000 * 0.0003048)
line_632_671.apply_template(config_601, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_632_671)
line_671_684 = gce.Line(bus_from=bus_671,
bus_to=bus_684,
length=300 * 0.0003048)
line_671_684.apply_template(config_604, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_671_684)
line_684_611 = gce.Line(bus_from=bus_684,
bus_to=bus_611,
length=300 * 0.0003048)
line_684_611.apply_template(config_605, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_684_611)
line_684_652 = gce.Line(bus_from=bus_684,
bus_to=bus_652,
length=800 * 0.0003048)
line_684_652.apply_template(config_607, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_684_652)
line_671_680 = gce.Line(bus_from=bus_671,
bus_to=bus_680,
length=1000 * 0.0003048)
line_671_680.apply_template(config_601, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_671_680)
line_671_675 = gce.Line(bus_from=bus_671,
bus_to=bus_675,
length=500 * 0.0003048)
line_671_675.apply_template(config_606, grid.Sbase, grid.fBase, logger)
grid.add_line(obj=line_671_675)
# ----------------------------------------------------------------------------------------------------------------------
# Run power flow
# ----------------------------------------------------------------------------------------------------------------------
res = gce.power_flow(grid=grid, options=gce.PowerFlowOptions(three_phase_unbalanced=True))
# ----------------------------------------------------------------------------------------------------------------------
# Show the results
# ----------------------------------------------------------------------------------------------------------------------
print(res.get_bus_df())
# ----------------------------------------------------------------------------------------------------------------------
# Save the grid
# ----------------------------------------------------------------------------------------------------------------------
gce.save_file(grid, "IEEE 13 Node Test Feeder.veragrid")
```
#### Phase A
Voltage magnitude:
Voltage angle:
#### Phase B:
Voltage magnitude:
Voltage angle:
#### Phase C:
Voltage magnitude:
Voltage angle:
