DC State Estimation

To perform the DC state estimation, we first need to have the PowerSystem type that has been created with the DC model, alongside the Measurement type that retains measurement data. Subsequently, we can formulate either the weighted least-squares (WLS) or the least absolute value (LAV) DC state estimation model encapsulated within the type DcStateEstimation using:

• dcStateEstimation,
• dcLavStateEstimation.

To obtain bus voltage angles and solve the DC state estimation problem, users can use the wrapper function:

• solve!.

After solving the DC state estimation, JuliaGrid provides function for computing powers:

• power!.

Alternatively, instead of using functions responsible for solving state estimation and computing powers, users can use the wrapper function:

• stateEstimation!.

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

Weighted Least-Squares Estimator

To solve the DC state estimation and derive WLS estimates using JuliaGrid, the process initiates by defining PowerSystem and Measurement types. Here is an illustrative example:

system, monitoring = ems()

addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2")
addBus!(system; label = "Bus 3")

addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.5)
addBranch!(system; label = "Branch 2", from = "Bus 1", to = "Bus 3", reactance = 0.2)
addBranch!(system; label = "Branch 3", from = "Bus 2", to = "Bus 3", reactance = 0.3)

@wattmeter(label = "Wattmeter ?")
addWattmeter!(monitoring; bus = "Bus 1", active = 0.6, variance = 1e-3)
addWattmeter!(monitoring; bus = "Bus 3", active = -0.4, variance = 1e-2)
addWattmeter!(monitoring; from = "Branch 1", active = 0.18, variance = 1e-4)
addWattmeter!(monitoring; to = "Branch 2", active = -0.42, variance = 1e-4)

The dcStateEstimation function serves to establish the DC state estimation problem:

analysis = dcStateEstimation(monitoring)

analysis = dcStateEstimation(monitoring, LDLt)

To obtain the bus voltage angles, the solve! function can be invoked as shown:

solve!(analysis)

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

julia> print(system.bus.label, analysis.voltage.angle)Bus 1: 0.0
Bus 2: -0.09000000000000001
Bus 3: -0.08399999999999999

Wrapper Function

JuliaGrid provides a wrapper function for DC state estimation analysis and also supports the computation of powers using the stateEstimation! function:

analysis = dcStateEstimation(monitoring)
stateEstimation!(analysis; verbose = 2)

Number of entries in the coefficient matrix: 10
Number of measurement functions:              4
Number of state variables:                    2

EXIT: The solution of the DC state estimation was found.

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 factorization [5, Sec. 3.2].

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

analysis = dcStateEstimation(monitoring, Orthogonal)
stateEstimation!(analysis)

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(analysis)

|-----------------|
| Bus Data        |
|-----------------|
| Label | Voltage |
|       |         |
|   Bus |   Angle |
|       |   [deg] |
|-------|---------|
| Bus 1 |  0.0000 |
| Bus 2 | -5.1566 |
| Bus 3 | -4.8128 |
|-----------------|

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

printBusData(analysis; label = "Bus 1", header = true)
printBusData(analysis; label = "Bus 2")
printBusData(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(analysis, file)
end

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 analysis 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 = dcLavStateEstimation(monitoring, Ipopt.Optimizer)

Setup Initial Primal Values

In JuliaGrid, the assignment of initial primal values for optimization variables takes place when the solve! function is executed. These values are derived from the voltage angles stored in the PowerSystem type and are assigned to the corresponding voltage field within the DcStateEstimation type:

julia> print(system.bus.label, analysis.voltage.angle)Bus 1: 0.0
Bus 2: 0.0
Bus 3: 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. One practical approach is to obtaine WLS estimator and then apply the resulting solution as the starting point for state estimation:

wls = dcStateEstimation(monitoring)
stateEstimation!(wls)

setInitialPoint!(analysis, wls)

As a result, the initial primal values will now reflect the outcome of the WLS state estimation solution:

julia> print(system.bus.label, analysis.voltage.angle)Bus 1: 0.0
Bus 2: -0.09000000000000001
Bus 3: -0.08399999999999999

Solution

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

stateEstimation!(analysis)

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

julia> print(system.bus.label, analysis.voltage.angle)Bus 1: 0.0
Bus 2: -0.08999999999999997
Bus 3: -0.08399999999999999

Measurement Set Update

We begin by creating the PowerSystem and Measurement types with the ems function. The DC model is then configured using dcModel! function. After that, we initialize the DcStateEstimation type through the dcStateEstimation function and solve the resulting state estimation problem:

system, monitoring = ems()

addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2")

addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.5)

dcModel!(system)

@wattmeter(label = "Wattmeter ?")
addWattmeter!(monitoring; bus = "Bus 2", active = -0.11, variance = 1e-3)
addWattmeter!(monitoring; from = "Branch 1", active = 0.09, variance = 1e-4)

analysis = dcStateEstimation(monitoring)
stateEstimation!(analysis)

Next, we modify the existing Measurement type using add and update functions. Then, we create the new DcStateEstimation type based on the modified system and solve the state estimation problem:

addWattmeter!(monitoring; to = "Branch 1", active = -0.12, variance = 1e-4)
updateWattmeter!(monitoring; label = "Wattmeter 1", status = 0)
updateWattmeter!(monitoring; label = "Wattmeter 2", active = 0.1, noise = false)

analysis = dcStateEstimation(monitoring)
stateEstimation!(analysis)

State Estimation Update

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

This advancement extends beyond the previous scenario where recreating the Measurement type was unnecessary, to now include the scenario where DcStateEstimation also does not need to be recreated. Such efficiency can be particularly advantageous in cases where JuliaGrid can reuse gain matrix factorization.

Let us now revisit our defined PowerSystem, Measurement and DcStateEstimation types:

system, monitoring = ems()

addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2")

addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.5)

dcModel!(system)

@wattmeter(label = "Wattmeter ?")
addWattmeter!(monitoring; bus = "Bus 2", active = -0.11, variance = 1e-3)
addWattmeter!(monitoring; from = "Branch 1", active = 0.09, variance = 1e-4)
addWattmeter!(monitoring; to = "Branch 1", active = -0.12, variance = 1e-4, status = 0)

analysis = dcStateEstimation(monitoring)
stateEstimation!(analysis)

Next, we modify the existing Measurement type as well as the DcStateEstimation type using add and update functions. We then immediately proceed to solve the state estimation problem:

updateWattmeter!(analysis; label = "Wattmeter 1", status = 0)
updateWattmeter!(analysis; label = "Wattmeter 2", active = 0.1)
updateWattmeter!(analysis; label = "Wattmeter 3", status = 1, noise = false)

stateEstimation!(analysis)

Reusing Weighted Least-Squares Matrix Factorization

Drawing from the preceding example, our focus now shifts to finding a solution involving modifications that entail adjusting the measurement value of the Wattmeter 2. It is important to note that these adjustments do not impact the variance or status of the measurement device, which can affect the gain matrix. To resolve this updated system, users can simply execute the solve! function:

updateWattmeter!(analysis; label = "Wattmeter 2", active = 0.091)

stateEstimation!(analysis)

Power Analysis

After obtaining the solution from the DC state estimation, calculating powers related to buses and branches is facilitated by using the power! function. For instance, let us consider the model for which we obtained the DC state estimation solution:

system, monitoring = ems()

addBus!(system; label = "Bus 1", type = 3, conductance = 1e-3)
addBus!(system; label = "Bus 2")
addBus!(system; label = "Bus 3")

addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.5)
addBranch!(system; label = "Branch 2", from = "Bus 1", to = "Bus 3", reactance = 0.2)
addBranch!(system; label = "Branch 3", from = "Bus 2", to = "Bus 3", reactance = 0.3)

addWattmeter!(monitoring; bus = "Bus 1", active = 0.6, variance = 1e-3)
addWattmeter!(monitoring; bus = "Bus 3", active = -0.4, variance = 1e-2)
addWattmeter!(monitoring; from = "Branch 1", active = 0.18, variance = 1e-4)
addWattmeter!(monitoring; to = "Branch 2", active = -0.42, variance = 1e-4)

analysis = dcStateEstimation(monitoring)
stateEstimation!(analysis)

We can compute active powers using the following function:

power!(analysis)

For example, active power injections corresponding to buses are:

julia> print(system.bus.label, analysis.power.injection.active)Bus 1: 0.6008333333333332
Bus 2: -0.19983333333333342
Bus 3: -0.4

Print Results in the REPL

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

@power(MW, pu)
printWattmeterData(analysis)

|---------------------------------------------------------------|
| Wattmeter Data                                                |
|---------------------------------------------------------------|
| Label |                     Active Power                      |
|       |                                                       |
|       | Measurement | Variance | Estimate | Residual | Status |
|       |        [MW] |     [MW] |     [MW] |     [MW] |        |
|-------|-------------|----------|----------|----------|--------|
| 1     |     60.0000 | 1.00e-01 |  60.0833 |  -0.0833 |      1 |
| 2     |    -40.0000 | 1.00e+00 | -40.0000 |   0.0000 |      1 |
| 3     |     18.0000 | 1.00e-02 |  17.9917 |   0.0083 |      1 |
| 4     |    -42.0000 | 1.00e-02 | -41.9917 |  -0.0083 |      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(analysis, io; style = false)
CSV.write("bus.csv", CSV.File(take!(io); delim = "|"))

Active Power Injection

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

julia> active = injectionPower(analysis; label = "Bus 1")0.6008333333333332

Active Power Injection from Generators

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

julia> active = supplyPower(analysis; label = "Bus 1")0.6008333333333332

Active Power Flow

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

julia> active = fromPower(analysis; label = "Branch 1")0.17991666666666667
julia> active = toPower(analysis; label = "Branch 1")-0.17991666666666667