Power and Current Analysis

In the following section, we have provided a list of functions that can be utilized for post-processing analysis. Once the voltage values are obtained through power flow analysis, optimal power flow analysis, or state estimation, these functions can be used to calculate power or current values. The specific procedures for computing these values depend on the chosen analysis, which are described in separate manuals for further information.

To load power system model API functionalities into the current scope, utilize the following command:

using JuliaGrid

AC Power Analysis
AC Current Analysis
DC Power Analysis

AC Power Analysis

JuliaGrid.power!Method
power!(system::PowerSystem, analysis::AC)

The function computes the active and reactive powers associated with buses, branches, and generators for AC analysis.

Updates

This function updates the power field of the AC abstract type by computing the following electrical quantities:

  • injection: Active and reactive power bus injections.
  • supply: Active and reactive power bus injections from the generators.
  • shunt: Active and reactive power values associated with shunt element at each bus.
  • from: Active and reactive power flows at the from-bus end of each branch.
  • to: Active and reactive power flows at the to-bus end of each branch.
  • charging: Active and reactive power values linked with branch charging admittances for each branch.
  • series Active and reactive power losses through each branch series impedance.
  • generator: Produced active and reactive power outputs of each generator (not for state estimation).

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
power!(system, analysis)
source
JuliaGrid.injectionPowerMethod
injectionPower(system::PowerSystem, analysis::AC, label)

The function returns the active and reactive power injections associated with a specific bus in the AC framework. The label keyword argument must match an existing bus label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = injectionPower(system, analysis; label = 1)
source
JuliaGrid.supplyPowerMethod
supplyPower(system::PowerSystem, analysis::AC, label)

The function returns the active and reactive power injections from the generators associated with a specific bus in the AC framework. The label keyword argument must match an existing bus label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = supplyPower(system, analysis; label = 1)
source
JuliaGrid.shuntPowerMethod
shuntPower(system::PowerSystem, analysis::AC, label)

The function returns the active and reactive power values of the shunt element associated with a specific bus in the AC framework. The label keyword argument must match an existing bus label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = shuntPower(system, analysis; label = 9)

```

source
JuliaGrid.fromPowerMethod
fromPower(system::PowerSystem, analysis::AC; label)

The function returns the active and reactive power flows at the from-bus end associated with a specific branch in the AC framework. The label keyword argument must match an existing branch label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = fromPower(system, analysis; label = 2)
source
JuliaGrid.toPowerMethod
toPower(system::PowerSystem, analysis::AC; label)

The function returns the active and reactive power flows at the to-bus end associated with a specific branch in the AC framework. The label keyword argument must match an existing branch label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = toPower(system, analysis; label = 2)
source
JuliaGrid.seriesPowerMethod
seriesPower(system::PowerSystem, analysis::AC; label)

The function returns the active and reactive power losses across the series impedance of a specific branch within the AC framework. The label keyword argument should correspond to an existing branch label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = seriesPower(system, analysis; label = 2)
source
JuliaGrid.chargingPowerMethod
chargingPower(system::PowerSystem, analysis::AC; label)

The function returns the active and reactive power values associated with the charging admittances of a specific branch in the AC framework. The label keyword argument must correspond to an existing branch label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = chargingPower(system, analysis; label = 2)
source
JuliaGrid.generatorPowerMethod
generatorPower(system::PowerSystem, analysis::AC)

The function returns the active and reactive powers associated with a specific generator in the AC framework. The label keyword argument must match an existing generator label.

Example

system = powerSystem("case14.h5")
acModel!(system)

analysis = newtonRaphson(system)
for i = 1:10
    stopping = mismatch!(system, analysis)
    if all(stopping .< 1e-8)
        break
    end
    solve!(system, analysis)
end
active, reactive = generatorPower(system, analysis; label = 1)
source

AC Current Analysis

JuliaGrid.current!Method
current!(system::PowerSystem, analysis::AC)

The function computes the currents in the polar coordinate system associated with buses and branches in the AC framework.

Updates

This function calculates various electrical quantities in the polar coordinate system:

  • injection: Current injections at each bus.
  • from: Current flows at each from-bus end of the branch.
  • to: Current flows at each to-bus end of the branch.
  • series: Current flows through the series impedance of the branch in the direction from the from-bus end to the to-bus end of the branch.

Example

using Ipopt

system = powerSystem("case14.h5")
acModel!(system)

analysis = acOptimalPowerFlow(system, Ipopt.Optimizer)
solve!(system, analysis)
current!(system, analysis)
source
JuliaGrid.injectionCurrentMethod
injectionCurrent(system::PowerSystem, analysis::AC; label)

The function returns the current injection in the polar coordinate system associated with a specific bus in the AC framework. The label keyword argument must match an existing bus label.

Example

using Ipopt

system = powerSystem("case14.h5")
acModel!(system)

analysis = acOptimalPowerFlow(system, Ipopt.Optimizer)
solve!(system, analysis)
magnitude, angle = injectionCurrent(system, analysis; label = 1)
source
JuliaGrid.fromCurrentMethod
fromCurrent(system::PowerSystem, analysis::AC; label)

The function returns the current in the polar coordinate system at the from-bus end associated with a specific branch in the AC framework. The label keyword argument must match an existing branch label.

Example

using Ipopt

system = powerSystem("case14.h5")
acModel!(system)

analysis = acOptimalPowerFlow(system, Ipopt.Optimizer)
solve!(system, analysis)
magnitude, angle = fromCurrent(system, analysis; label = 2)
source
JuliaGrid.toCurrentMethod
toCurrent(system::PowerSystem, analysis::AC; label)

The function returns the current in the polar coordinate system at the to-bus end associated with a specific branch in the AC framework. The label keyword argument must match an existing branch label.

Example

using Ipopt

system = powerSystem("case14.h5")
acModel!(system)

analysis = acOptimalPowerFlow(system, Ipopt.Optimizer)
solve!(system, analysis)
magnitude, angle = toCurrent(system, analysis; label = 2)
source
JuliaGrid.seriesCurrentMethod
seriesCurrent(system::PowerSystem, analysis::AC; label)

The function returns the current in the polar coordinate system through series impedance associated with a specific branch in the direction from the from-bus end to the to-bus end of the branch within the AC framework. The label keyword argument must match an existing branch label.

Example

using Ipopt

system = powerSystem("case14.h5")
acModel!(system)

analysis = acOptimalPowerFlow(system, Ipopt.Optimizer)
solve!(system, analysis)
magnitude, angle = seriesCurrent(system, analysis; label = 2)
source

DC Power Analysis

JuliaGrid.power!Method
power!(system::PowerSystem, analysis::DC)

The function calculates the active power values related to buses, branches, and generators within the DC analysis framework.

Updates

This function updates the power field of the DC abstract type by computing the following electrical quantities:

  • injection: Active power injections at each bus.
  • supply: Active power injections from the generators at each bus.
  • from: Active power flows at each from-bus end of the branch.
  • to: Active power flows at each to-bus end of the branch.
  • generator: Output active powers of each generator (excluding for state estimation).

Example

system = powerSystem("case14.h5")
dcModel!(system)

analysis = dcPowerFlow(system)
solve!(system, analysis)
power!(system, analysis)
source
JuliaGrid.injectionPowerMethod
injectionPower(system::PowerSystem, analysis::DC; label)

The function returns the active power injection associated with a specific bus in the DC framework. The label keyword argument must match an existing bus label.

Example

system = powerSystem("case14.h5")
dcModel!(system)

analysis = dcPowerFlow(system)
solve!(system, analysis)
injection = injectionPower(system, analysis; label = 2)
source
JuliaGrid.supplyPowerMethod
supplyPower(system::PowerSystem, analysis::DC; label)

The function returns the active power injection from the generators associated with a specific bus in the DC framework. The label keyword argument must match an existing bus label.

Example

system = powerSystem("case14.h5")
dcModel!(system)

analysis = dcPowerFlow(system)
solve!(system, analysis)
supply = supplyPower(system, analysis; label = 2)
source
JuliaGrid.fromPowerMethod
fromPower(system::PowerSystem, analysis::DC; label)

The function returns the active power flow at the from-bus end associated with a specific branch in the DC framework. The label keyword argument must match an existing branch label.

Example

system = powerSystem("case14.h5")
dcModel!(system)

analysis = dcPowerFlow(system)
solve!(system, analysis)
from = fromPower(system, analysis; label = 2)
source
JuliaGrid.toPowerMethod
toPower(system::PowerSystem, analysis::DC; label)

The function returns the active power flow at the to-bus end associated with a specific branch in the DC framework. The label keyword argument must match an existing branch label.

Example

system = powerSystem("case14.h5")
dcModel!(system)

analysis = dcPowerFlow(system)
solve!(system, analysis)
to = toPower(system, analysis; label = 2)
source
JuliaGrid.generatorPowerMethod
generatorPower(system::PowerSystem, analysis::DC; label)

This function returns the output active power associated with a specific generator in the DC framework. The label keyword argument must match an existing generator label.

Example

system = powerSystem("case14.h5")
dcModel!(system)

analysis = dcPowerFlow(system)
solve!(system, analysis)
generator = generatorPower(system, analysis; label = 1)
source