AC State Estimation

To perform nonlinear or AC state estimation, the initial requirement is to have the PowerSystem type configured with the AC model, along with the Measurement type storing measurement data. Next, we can develop either the weighted least-squares (WLS) model, utilizing the Gauss-Newton method, or the least absolute value (LAV) model. These models are encapsulated within the ACStateEstimation type:

• gaussNewton,
• acLavStateEstimation.

To obtain bus voltages and solve the state estimation problem, users can either implement an iterative process for the WLS model or simply execute the following function for the LAV model:

• solve!.

After solving the AC state estimation, JuliaGrid provides functions for computing powers and currents:

• power!,
• current!.

Alternatively, instead of designing their own iteration process for the Gauss-Newton method or using the function responsible for solving the LAV model, and computing powers and currents, users can use the wrapper function:

• stateEstimation!.

Users can also access specialized functions for computing specific types of powers or currents for individual buses, branches, or generators within the power system.

After obtaining the bus voltages, when the user employs the WLS model, they can check if the measurement set contains outliers through bad data analysis and remove those measurements using:

• residualTest!.

Bus Type Modification

In AC state estimation, it is necessary to designate a slack bus, where the bus voltage angle is known. Therefore, when establishing the ACStateEstimation type, the initially assigned slack bus is evaluated and may be altered. If the designated slack bus (type = 3) lacks a connected in-service generator, it will be changed to a demand bus (type = 1). Conversely, the first generator bus (type = 2) with an active in-service generator linked to it will be reassigned as the new slack bus (type = 3).

Weighted Least-Squares Estimator

To begin, we will define the PowerSystem and Measurement types:

system = powerSystem()
device = measurement()

addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2", type = 1, active = 0.6, reactive = 0.1)
addBus!(system; label = "Bus 3", type = 1, active = 0.5, reactive = 0.2)

@branch(resistance = 0.14, conductance = 1e-4, susceptance = 0.04)
addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.25)
addBranch!(system; label = "Branch 2", from = "Bus 1", to = "Bus 3", reactance = 0.35)
addBranch!(system; label = "Branch 3", from = "Bus 2", to = "Bus 3", reactance = 0.16)

addGenerator!(system; label = "Generator 1", bus = "Bus 1", active = 1.2, reactive = 0.3)

@voltmeter(label = "Voltmeter ? (!)")
addVoltmeter!(system, device; bus = "Bus 1", magnitude = 1.0, variance = 1e-5)

@ammeter(label = "Ammeter ? (!)")
addAmmeter!(system, device; from = "Branch 3", magnitude = 0.0356, variance = 1e-3)
addAmmeter!(system, device; to = "Branch 2", magnitude = 0.5892, variance = 1e-3)

@wattmeter(label = "Wattmeter ? (!)")
addWattmeter!(system, device; from = "Branch 1", active = 0.7067, variance = 1e-4)
addWattmeter!(system, device; bus = "Bus 2", active = -0.6, variance = 2e-4)

@varmeter(label = "Varmeter ? (!)")
addVarmeter!(system, device; from = "Branch 1", reactive = 0.2125, variance = 1e-4)
addVarmeter!(system, device; bus = "Bus 2", reactive = -0.1, variance = 1e-5)

Next, to establish the AC state estimation model, we will utilize the gaussNewton function:

analysis = gaussNewton(system, device)

analysis = gaussNewton(system, device, LDLt)

Setup Initial Voltages

The initial voltages for the Gauss-Newton method are determined based on the specified initial voltage magnitudes and angles within the buses of the PowerSystem type. These values are then forwarded to the ACStateEstimation during the execution of the gaussNewton function. Therefore, the initial voltages in this example are as follows:

julia> print(system.bus.label, analysis.voltage.magnitude, analysis.voltage.angle)Bus 1: 1.0, 0.0
Bus 2: 1.0, 0.0
Bus 3: 1.0, 0.0

Users have the flexibility to modify these vectors according to their own requirements in order to adjust the initial voltages. For instance, users can conduct an AC power flow analysis and utilize the obtained solution as the initial voltages for AC state estimation:

powerFlow = newtonRaphson(system)
for iteration = 1:10
    mismatch!(system, powerFlow)
    solve!(system, powerFlow)
end

setInitialPoint!(powerFlow, analysis)

State Estimator

To conduct an iterative process using the Gauss-Newton method, it is essential to include the solve! function inside the iteration loop. For example:

for iteration = 1:20
    solve!(system, analysis)
end

Once the state estimator is obtained, users can access the bus voltage magnitudes and angles using:

julia> print(system.bus.label, analysis.voltage.magnitude, analysis.voltage.angle)Bus 1: 0.9999999893867697, 0.0
Bus 2: 0.8551739067243623, -0.16932686045421688
Bus 3: 0.8460326304864767, -0.17123533032549176

Breaking the Iterative Process

The iterative process can be terminated using the solve! function. The following code demonstrates how to utilize this function to break out of the iteration loop:

analysis = gaussNewton(system, device)
for iteration = 1:20
    stopping = solve!(system, analysis)
    if stopping < 1e-8
        println("Solution Found.")
        break
    end
end

Solution Found.

The solve! function returns the maximum absolute values of the state variable increment, which are commonly used as a convergence criterion in the iterative Gauss-Newton algorithm.

Wrapper Function

JuliaGrid provides a wrapper function stateEstimation! for AC state estimation. Hence, it offers a way to solve it using the Gauss-Newton method with reduced implementation effort:

analysis = gaussNewton(system, device)
stateEstimation!(system, analysis; verbose = 3)

Number of wattmeters: 2   Number of varmeters: 2   Number of voltmeters: 1
  In-service:         2     In-service:        2     In-service:         1
  Out-of-service:     0     Out-of-service:    0     Out-of-service:     0

Number of ammeters:   2   Number of PMUs:      0   Number of devices:    7
  In-service:         2     In-service:        0     In-service:         7
  Out-of-service:     0     Out-of-service:    0     Out-of-service:     0

Number of entries in the Jacobian matrix: 29
Number of measurement functions:           7
Number of state variables:                 5
Number of buses:                           3
Number of branches:                        3

Iteration   Maximum Increment   Objective Value
        1          2.0358e-01    8.01893670e+03
        2          4.8274e-02    2.97107558e+02
        3          1.8715e-03    2.16194897e-01
        4          2.5886e-06    4.71061964e-06
        5          1.5800e-10    3.54147541e-07

                     Minimum Value   Maximum Value
Magnitude Increment:   -1.0215e-11      1.1398e-10
Angle Increment:       -1.5800e-10     -7.4174e-11

EXIT: The solution was found using the Gauss-Newton method in 5 iterations.

Inclusion of PMUs in Rectangular Coordinates

In the example above, our focus is solely on solving the AC state estimation using SCADA measurements. However, users have the option to also integrate PMUs into the AC state estimation, either in the rectangular or polar coordinate system.

The default approach is to include PMUs in the rectangular coordinate system:

@pmu(label = "PMU ? (!)")
addPmu!(system, device; to = "Branch 1", magnitude = 0.7466, angle = 2.8011)

In the case of the rectangular system, inclusion resolves ill-conditioned problems arising in polar coordinates due to small values of current magnitudes. However, this approach's main disadvantage is related to measurement errors, as measurement errors correspond to polar coordinates. Therefore, the covariance matrix must be transformed from polar to rectangular coordinates [3]. As a result, measurement errors of a single PMU are correlated, and the covariance matrix does not have a diagonal form. Despite that, the measurement error covariance matrix is usually considered as a diagonal matrix, affecting the accuracy of the state estimation.

In the example above, we specifically include PMUs where measurement error correlations are disregarded. This is evident through the precision matrix, which maintains a diagonal form:

julia> analysis = gaussNewton(system, device);
julia> analysis.method.precision9×9 SparseArrays.SparseMatrixCSC{Float64, Int64} with 9 stored entries:
 100000.0      ⋅       ⋅        ⋅       ⋅   …        ⋅    ⋅          ⋅
       ⋅   1000.0      ⋅        ⋅       ⋅            ⋅    ⋅          ⋅
       ⋅       ⋅   1000.0       ⋅       ⋅            ⋅    ⋅          ⋅
       ⋅       ⋅       ⋅   10000.0      ⋅            ⋅    ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅   5000.0           ⋅    ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅   …        ⋅    ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅      100000.0   ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅            ⋅   1.05192e8   ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅            ⋅    ⋅         1.64807e8

Lastly, we incorporate correlation into our model by adding a new PMU with the desired error correlation:

addPmu!(system, device; bus = "Bus 3", magnitude = 0.846, angle = -0.1712, correlated = true)

Now, we can observe the precision matrix that does not hold a diagonal form:

julia> analysis = gaussNewton(system, device);
julia> analysis.method.precision11×11 SparseArrays.SparseMatrixCSC{Float64, Int64} with 13 stored entries:
 100000.0      ⋅       ⋅        ⋅   …   ⋅          ⋅          ⋅
       ⋅   1000.0      ⋅        ⋅       ⋅          ⋅          ⋅
       ⋅       ⋅   1000.0       ⋅       ⋅          ⋅          ⋅
       ⋅       ⋅       ⋅   10000.0      ⋅          ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅          ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅   …   ⋅          ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅          ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅          ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅      1.64807e8   ⋅          ⋅
       ⋅       ⋅       ⋅        ⋅       ⋅         1.01153e8  6.66801e6
       ⋅       ⋅       ⋅        ⋅   …   ⋅         6.66801e6  1.38567e8

Inclusion of PMUs in Polar Coordinates

The second approach involves incorporating these measurements into the polar coordinate system. For instance:

addPmu!(system, device; from = "Branch 1", magnitude = 0.7379, angle = -0.2921, polar = true)

This inclusion of PMUs provides more accurate state estimates compared to rectangular inclusion, but demands longer computing time. PMUs are handled in the same manner as SCADA measurements. However, this approach is susceptible to ill-conditioned problems arising in polar coordinates due to small values of current magnitudes [3, 4].

Alternative Formulation

The resolution of the WLS state estimation problem using the conventional method typically progresses smoothly. However, it is widely acknowledged that in certain situations common to real-world systems, this method can be vulnerable to numerical instabilities. Such conditions might impede the algorithm from finding a satisfactory solution. In such cases, users may opt for an alternative formulation of the WLS state estimation, namely, employing an approach called orthogonal method [5, Sec. 3.2].

This approach is suitable when measurement errors are uncorrelated, and the precision matrix remains diagonal. Therefore, as a preliminary step, we need to eliminate the correlation, as we did previously:

updatePmu!(system, device; label = "PMU 2 (Bus 3)", correlated = false)

Subsequently, by specifying the Orthogonal argument in the gaussNewton function, JuliaGrid implements a more robust approach to obtain the WLS estimator, which proves particularly beneficial when substantial differences exist among measurement variances:

analysis = gaussNewton(system, device, Orthogonal)
stateEstimation!(system, analysis; verbose = 1)

EXIT: The solution was found using the Gauss-Newton method in 8 iterations.

Print Results in the REPL

Users have the option to print the results in the REPL using any units that have been configured, such as:

@voltage(pu, deg)
printBusData(system, analysis)

|-----------------------------|
| Bus Data                    |
|-----------------------------|
| Label |       Voltage       |
|       |                     |
|   Bus | Magnitude |   Angle |
|       |      [pu] |   [deg] |
|-------|-----------|---------|
| Bus 1 |    0.9999 |  0.0000 |
| Bus 2 |    0.8551 | -9.7028 |
| Bus 3 |    0.8460 | -9.8090 |
|-----------------------------|

Next, users can easily customize the print results for specific buses, for example:

printBusData(system, analysis; label = "Bus 1", header = true)
printBusData(system, analysis; label = "Bus 2")
printBusData(system, analysis; label = "Bus 3", footer = true)

Save Results to a File

Users can also redirect print output to a file. For example, data can be saved in a text file as follows:

open("bus.txt", "w") do file
    printBusData(system, analysis, file)
end

Bad Data Processing

After acquiring the WLS solution using the Gauss-Newton method, users can conduct bad data analysis employing the largest normalized residual test. Continuing with our defined power system and measurement set, let us introduce a new measurement. Upon proceeding to find the solution for this updated state:

addWattmeter!(system, device; from = "Branch 2", active = 31.1)

analysis = gaussNewton(system, device)
stateEstimation!(system, analysis; verbose = 1)

EXIT: The solution was found using the Gauss-Newton method in 11 iterations.

Here, we can observe the impact of the outlier on the solution:

julia> print(system.bus.label, analysis.voltage.magnitude, analysis.voltage.angle)Bus 1: 1.346801341763577, 0.0
Bus 2: 1.1982060394195844, -0.12025815367863564
Bus 3: 0.848589007494443, -0.183265103284885

Following the solution acquisition, we can verify the presence of erroneous data. Detection of such data is determined by the threshold keyword. If the largest normalized residual's value exceeds the threshold, the measurement will be identified as bad data and consequently removed from the AC state estimation model:

outlier = residualTest!(system, device, analysis; threshold = 4.0)

Users can examine the data obtained from the bad data analysis:

julia> outlier.detecttrue
julia> outlier.maxNormalizedResidual3021.9944591812973
julia> outlier.label"Wattmeter 3 (From Branch 2)"

Hence, upon detecting bad data, the detect variable will hold true. The maxNormalizedResidual variable retains the value of the largest normalized residual, while the label contains the label of the measurement identified as bad data. JuliaGrid will mark the respective measurement as out-of-service within the Measurement type.

After removing bad data, a new estimate can be computed without considering this specific measurement. The user has the option to either restart the gaussNewton function or proceed directly to the iteration loop. However, if the latter option is chosen, using voltages obtained with outlier presence as the initial point could significantly impede algorithm convergence. To avoid this undesirable outcome, the user should first establish a new initial point and commence the iteration procedure. For instance:

setInitialPoint!(system, analysis)
stateEstimation!(system, analysis; verbose = 1)

EXIT: The solution was found using the Gauss-Newton method in 8 iterations.

Consequently, we obtain a new solution devoid of the impact of the outlier measurement:

julia> print(system.bus.label, analysis.voltage.magnitude, analysis.voltage.angle)Bus 1: 0.9998752176863828, 0.0
Bus 2: 0.8550575043077884, -0.16934541287631205
Bus 3: 0.8459996766101084, -0.1711998537428169

Least Absolute Value Estimator

The LAV method presents an alternative estimation technique known for its increased robustness compared to WLS. While the WLS method relies on specific assumptions regarding measurement errors, robust estimators like LAV are designed to maintain unbiasedness even in the presence of various types of measurement errors and outliers. This characteristic often eliminates the need for extensive bad data processing procedures [5, Ch. 6]. However, it is important to note that achieving robustness typically involves increased computational complexity.

To obtain an LAV estimator, users need to employ one of the solvers listed in the JuMP documentation. In many common scenarios, the Ipopt solver proves sufficient to obtain a solution:

using Ipopt

analysis = acLavStateEstimation(system, device, Ipopt.Optimizer; verbose = 1)

Setup Initial Primal Values

In JuliaGrid, the assignment of initial primal values for optimization variables takes place when the solve! function is executed. Initial primal values are determined based on the voltage fields within the ACStateEstimation type. By default, these values are established using the initial bus voltage magnitudes and angles from PowerSystem type:

julia> print(system.bus.label, analysis.voltage.magnitude, analysis.voltage.angle)Bus 1: 1.0, 0.0
Bus 2: 1.0, 0.0
Bus 3: 1.0, 0.0

Users have the flexibility to customize these values according to their requirements, and they will be utilized as the initial primal values when executing the solve! function. Additionally, the setInitialPoint! function allows users to configure the initial point as required.

Solution

To solve the formulated LAV state estimation model, simply execute the following function:

solve!(system, analysis)

EXIT: Solved To Acceptable Level.

Upon obtaining the solution, access the bus voltage magnitudes and angles using:

julia> print(system.bus.label, analysis.voltage.magnitude, analysis.voltage.angle)Bus 1: 0.9999955872956292, 0.0
Bus 2: 0.8551679159257191, -0.16932741993230138
Bus 3: 0.8460243689616362, -0.17123746785110638

Measurement Set Update

After establishing the Measurement type using the measurement function, users gain the capability to incorporate new measurement devices or update existing ones.

Once updates are completed, users can seamlessly progress towards generating the ACStateEstimation type using the gaussNewton or acLavStateEstimation function. Ultimately, resolving the AC state estimation is achieved through the utilization of the solve! function:

system = powerSystem()
device = measurement() # <- Initialize the Measurement instance

addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2", type = 1, active = 0.6, reactive = 0.1)
addBus!(system; label = "Bus 3", type = 1, active = 0.5, reactive = 0.2)

@branch(resistance = 0.14, conductance = 1e-4, susceptance = 0.04)
addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.25)
addBranch!(system; label = "Branch 2", from = "Bus 1", to = "Bus 3", reactance = 0.35)
addBranch!(system; label = "Branch 3", from = "Bus 2", to = "Bus 3", reactance = 0.16)

addGenerator!(system; label = "Generator 1", bus = "Bus 1", active = 1.2, reactive = 0.3)

@voltmeter(label = "Voltmeter ? (!)", variance = 1e-3)
addVoltmeter!(system, device; bus = "Bus 1", magnitude = 1.0)

@wattmeter(label = "Wattmeter ? (!)", varianceBus = 1e-2)
addWattmeter!(system, device; from = "Branch 1", active = 0.7067)
addWattmeter!(system, device; bus = "Bus 2", active = -0.6)

@varmeter(label = "Varmeter ? (!)", varianceFrom = 1e-3)
addVarmeter!(system, device; from = "Branch 1", reactive = 0.2125)
addVarmeter!(system, device; bus = "Bus 2", reactive = -0.1)

@pmu(label = "PMU ? (!)")
addPmu!(system, device; bus = "Bus 2", magnitude = 0.8552, angle = -0.1693)

analysis = gaussNewton(system, device) # <- Build ACStateEstimation for the defined model
stateEstimation!(system, analysis)

addWattmeter!(system, device; from = "Branch 3", active = 0.0291)
updateWattmeter!(system, device; label = "Wattmeter 2 (Bus 2)", variance = 1e-4)

addVarmeter!(system, device; to = "Branch 3", reactive = -0.037, variance = 1e-5)
updateVarmeter!(system, device; label = "Varmeter 2 (Bus 2)", reactive = -0.11)

updatePmu!(system, device; label = "PMU 1 (Bus 2)", polar = false)

analysis = gaussNewton(system, device) # <- Build ACStateEstimation for the updated model
stateEstimation!(system, analysis)

State Estimation Update

An advanced methodology involves users establishing the ACStateEstimation type using gaussNewton or acLavStateEstimation just once. After this initial setup, users can seamlessly modify existing measurement devices without the need to recreate the ACStateEstimation type.

This advancement extends beyond the previous scenario where recreating the Measurement type was unnecessary, to now include the scenario where ACStateEstimation also does not need to be recreated.

We can modify the prior example to achieve the same model without establishing ACStateEstimation twice:

system = powerSystem()
device = measurement() # <- Initialize the Measurement instance

addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2", type = 1, active = 0.6, reactive = 0.1)
addBus!(system; label = "Bus 3", type = 1, active = 0.5, reactive = 0.2)

@branch(resistance = 0.14, conductance = 1e-4, susceptance = 0.04)
addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.25)
addBranch!(system; label = "Branch 2", from = "Bus 1", to = "Bus 3", reactance = 0.35)
addBranch!(system; label = "Branch 3", from = "Bus 2", to = "Bus 3", reactance = 0.16)

addGenerator!(system; label = "Generator 1", bus = "Bus 1", active = 1.2, reactive = 0.3)

@voltmeter(label = "Voltmeter ? (!)")
addVoltmeter!(system, device; bus = "Bus 1", magnitude = 1.0)

@wattmeter(label = "Wattmeter ? (!)", varianceBus = 1e-2)
addWattmeter!(system, device; from = "Branch 1", active = 0.7067)
addWattmeter!(system, device; bus = "Bus 2", active = -0.6)
addWattmeter!(system, device; from = "Branch 3", active = 0.0291, status = 0)

@varmeter(label = "Varmeter ? (!)", varianceFrom = 1e-3)
addVarmeter!(system, device; from = "Branch 1", reactive = 0.2125)
addVarmeter!(system, device; bus = "Bus 2", reactive = -0.1)
addVarmeter!(system, device; to = "Branch 3", reactive = -0.037, variance = 1e-5, status = 0)

@pmu(label = "PMU ? (!)")
addPmu!(system, device; bus = "Bus 2", magnitude = 0.8552, angle = -0.1693)

analysis = gaussNewton(system, device) # <- Build ACStateEstimation for the defined model
stateEstimation!(system, analysis)

updateWattmeter!(system, device, analysis; label = "Wattmeter 3 (From Branch 3)", status = 1)
updateWattmeter!(system, device, analysis; label = "Wattmeter 2 (Bus 2)", variance = 1e-4)

updateVarmeter!(system, device, analysis; label = "Varmeter 3 (To Branch 3)", status = 1)
updateVarmeter!(system, device, analysis; label = "Varmeter 2 (Bus 2)", reactive = -0.11)

updatePmu!(system, device, analysis; label = "PMU 1 (Bus 2)", polar = false)

# <- No need for re-build; we have already updated the existing ACStateEstimation instance
stateEstimation!(system, analysis)

Power and Current Analysis

After obtaining the solution from the AC state estimation, we can calculate various electrical quantities related to buses and branches using the power! and current! functions. For instance, let us consider the model for which we obtained the AC state estimation solution:

system = powerSystem()
device = measurement()

addBus!(system; label = "Bus 1", type = 3, susceptance = 0.002)
addBus!(system; label = "Bus 2", type = 1, active = 0.1, reactive = 0.01)
addBus!(system; label = "Bus 3", type = 1, active = 2.5, reactive = 0.2)

@branch(resistance = 0.02, conductance = 1e-4, susceptance = 0.04)
addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.05)
addBranch!(system; label = "Branch 2", from = "Bus 1", to = "Bus 3", reactance = 0.05)
addBranch!(system; label = "Branch 3", from = "Bus 2", to = "Bus 3", reactance = 0.03)

addGenerator!(system; label = "Generator 1", bus = "Bus 1", active = 3.2, reactive = 0.3)

addWattmeter!(system, device; from = "Branch 1", active = 1.046, variance = 1e-2)
addWattmeter!(system, device; bus = "Bus 2", active = -0.1, variance = 1e-3)
addWattmeter!(system, device; from = "Branch 3", active = 0.924, variance = 1e-3)

addVarmeter!(system, device; from = "Branch 1", reactive = 0.059, variance = 1e-3)
addVarmeter!(system, device; bus = "Bus 2", reactive = -0.01, variance = 1e-2)
addVarmeter!(system, device; to = "Branch 3", reactive = -0.044, variance = 1e-3)

analysis = gaussNewton(system, device)
stateEstimation!(system, analysis; verbose = 1)

EXIT: The solution was found using the Gauss-Newton method in 7 iterations.

We can now utilize the provided functions to compute powers and currents:

power!(system, analysis)
current!(system, analysis)

For instance, if we want to show the active power injections and the from-bus current angles, we can employ the following code:

julia> print(system.bus.label, analysis.power.injection.active)Bus 1: 2.6951392445386304
Bus 2: -0.09999730749216829
Bus 3: -2.499287793840309
julia> print(system.branch.label, analysis.current.from.angle)Branch 1: -0.056343895480383585
Branch 2: -0.1541309681169004
Branch 3: -0.08788060713202979

Print Results in the REPL

Users can utilize any of the print functions outlined in the Print API. For example, to print state estimation data related to wattmeters, we can use:

@power(MW, pu)
printWattmeterData(system, device, analysis)

|---------------------------------------------------------------|
| Wattmeter Data                                                |
|---------------------------------------------------------------|
| Label |                     Active Power                      |
|       |                                                       |
|       | Measurement | Variance | Estimate | Residual | Status |
|       |        [MW] |     [MW] |     [MW] |     [MW] |        |
|-------|-------------|----------|----------|----------|--------|
| 1     |    104.6000 | 1.00e+00 | 104.6026 |  -0.0026 |      1 |
| 2     |    -10.0000 | 1.00e-01 |  -9.9997 |  -0.0003 |      1 |
| 3     |     92.4000 | 1.00e-01 |  92.3997 |   0.0003 |      1 |
|---------------------------------------------------------------|

Save Results to a CSV File

For CSV output, users should first generate a simple table with style = false, and then save it to a CSV file:

using CSV

io = IOBuffer()
printWattmeterData(system, device, analysis, io; style = false)
CSV.write("bus.csv", CSV.File(take!(io); delim = "|"))

Active and Reactive Power Injection

To calculate the active and reactive power injection associated with a specific bus, the function can be used:

julia> active, reactive = injectionPower(system, analysis; label = "Bus 1")(2.6951392445386304, 0.3132044957940553)

Active and Reactive Power Injection from Generators

To calculate the active and reactive power injection from the generators at a specific bus, the function can be used:

julia> active, reactive = supplyPower(system, analysis; label = "Bus 1")(2.6951392445386304, 0.3132044957940553)

Active and Reactive Power at Shunt Element

To calculate the active and reactive power associated with shunt element at a specific bus, the function can be used:

julia> active, reactive = shuntPower(system, analysis; label = "Bus 1")(0.0, -0.0020066379525219783)

Active and Reactive Power Flow

Similarly, we can compute the active and reactive power flow at both the from-bus and to-bus ends of the specific branch by utilizing the provided functions below:

julia> active, reactive = fromPower(system, analysis; label = "Branch 2")(1.6491130561317449, 0.2562114962034585)
julia> active, reactive = toPower(system, analysis; label = "Branch 2")(-1.5932872079562321, -0.15530958620622967)

Active and Reactive Power at Charging Admittances

To calculate the active and reactive power linked with branch charging admittances of the particular branch, the function can be used:

julia> active, reactive = chargingPower(system, analysis; label = "Branch 1")(9.800363372293906e-5, -0.03920145348917562)

Active powers indicate active losses within the branch's charging admittances. Moreover, charging admittances injected reactive powers into the power system due to their capacitive nature, as denoted by a negative sign.

Active and Reactive Power at Series Impedance

To calculate the active and reactive power across the series impedance of the branch, the function can be used:

julia> active, reactive = seriesPower(system, analysis; label = "Branch 2")(0.05572979175205785, 0.13932447938014464)

The active power also considers active losses originating from the series resistance of the branch, while the reactive power represents reactive losses resulting from the impedance's inductive characteristics.

Current Injection

To calculate the current injection associated with a specific bus, the function can be used:

julia> magnitude, angle = injectionCurrent(system, analysis; label = "Bus 1")(2.708785622827216, -0.11569193435017798)

Current Flow

We can compute the current flow at both the from-bus and to-bus ends of the specific branch by utilizing the provided functions below:

julia> magnitude, angle = fromCurrent(system, analysis; label = "Branch 2")(1.6661346611128287, -0.1541309681169004)
julia> magnitude, angle = toCurrent(system, analysis; label = "Branch 2")(1.6709804057015962, 2.9641707920912417)

Current Through Series Impedance

To calculate the current passing through the series impedance of the branch in the direction from the from-bus end to the to-bus end, we can use the following function:

julia> magnitude, angle = seriesCurrent(system, analysis; label = "Branch 2")(1.6692781636392697, -0.16599467621735334)