DC State Estimation
In this example, we monitor a 6-bus power system, shown in Figure 1, and estimate bus voltage angles using DC state estimation.
Figure 1: The 6-bus power system.
Users can download a Julia script containing the scenarios from this section using the following link.
We start by defining the units for voltage angles, which will be used throughout this example:
@voltage(pu, deg)Next, the power system is defined by specifying buses, branches, and generators with cost functions:
system = powerSystem()
addBus!(system; label = "Bus 1", type = 3)
addBus!(system; label = "Bus 2", type = 1, active = 0.217)
addBus!(system; label = "Bus 3", type = 1, active = 0.478)
addBus!(system; label = "Bus 4", type = 2, active = 0.076)
addBus!(system; label = "Bus 5", type = 1, active = 0.112)
addBus!(system; label = "Bus 6", type = 2, active = 0.295)
addBranch!(system; label = "Branch 1", from = "Bus 1", to = "Bus 2", reactance = 0.05)
addBranch!(system; label = "Branch 2", from = "Bus 2", to = "Bus 3", reactance = 0.23)
addBranch!(system; label = "Branch 3", from = "Bus 3", to = "Bus 4", reactance = 0.19)
addBranch!(system; label = "Branch 4", from = "Bus 4", to = "Bus 5", reactance = 0.17)
addBranch!(system; label = "Branch 5", from = "Bus 5", to = "Bus 6", reactance = 0.04)
addBranch!(system; label = "Branch 6", from = "Bus 1", to = "Bus 6", reactance = 0.21)
addBranch!(system; label = "Branch 7", from = "Bus 2", to = "Bus 6", reactance = 0.13)
addBranch!(system; label = "Branch 8", from = "Bus 5", to = "Bus 2", reactance = 0.34)
addGenerator!(system; label = "Generator 1", bus = "Bus 1", active = 0.8, maxActive = 2.3)
addGenerator!(system; label = "Generator 2", bus = "Bus 4", active = 0.4, maxActive = 2.3)
cost!(system; generator = "Generator 1", active = 2, polynomial = [1100.0; 500.0; 150.0])
cost!(system; generator = "Generator 2", active = 2, polynomial = [1500.0; 700.0; 140.0])After defining the power system data, a DC model is generated, including key matrices and vectors for analysis:
dcModel!(system)Measurement Model
Next, we define the measurements. The question is how to obtain measurement values. In this example, synthetic measurements are generated using DC optimal power flow results.
To start, we initialize the measurement variable:
device = measurement()DC Optimal Power Flow
To obtain bus voltage angless, we solve the DC optimal power flow. Using these values, we compute the active power associated with buses and branches:
powerFlow = dcOptimalPowerFlow(system, Ipopt.Optimizer)
solve!(system, powerFlow)
power!(system, powerFlow)Active Power Injection Measurements
Active power injection measurements will be obtained from the DC optimal power flow analysis:
printBusData(system, powerFlow)|---------------------------------------------------------------------|
| Bus Data |
|---------------------------------------------------------------------|
| Label | Voltage | Power Generation | Power Demand | Power Injection |
| | | | | |
| Bus | Angle | Active | Active | Active |
| | [deg] | [pu] | [pu] | [pu] |
|-------|---------|------------------|--------------|-----------------|
| Bus 1 | 0.0000 | 0.7181 | 0.0000 | 0.7181 |
| Bus 2 | -1.5346 | 0.0000 | 0.2170 | -0.2170 |
| Bus 3 | -4.1563 | 0.0000 | 0.4780 | -0.4780 |
| Bus 4 | -1.1185 | 0.4599 | 0.0760 | 0.3839 |
| Bus 5 | -2.1399 | 0.0000 | 0.1120 | -0.1120 |
| Bus 6 | -2.1948 | 0.0000 | 0.2950 | -0.2950 |
|---------------------------------------------------------------------|Next, these measurements are defined:
@wattmeter(label = "Wattmeter ?")
for (label, idx) in system.bus.label
Pᵢ = powerFlow.power.injection.active[idx]
addWattmeter!(system, device; bus = label, active = Pᵢ, variance = 1e-4, noise = true)
endEnabling noise = true adds white Gaussian noise with a variance of 1e-4 to the exact values, generating the final measurement values.
Active Power Flow Measurements
Next, we will include a certain number of active power flow measurements using the results from the DC optimal power flow analysis:
printBranchData(system, powerFlow)|-----------------------------------------------------------------------|
| Branch Data |
|-----------------------------------------------------------------------|
| Label | From-Bus Power | To-Bus Power | Status |
| | | | |
| Branch | From-Bus | To-Bus | Active | Active | |
| | | | [pu] | [pu] | |
|----------|----------|--------|----------------|--------------|--------|
| Branch 1 | Bus 1 | Bus 2 | 0.5357 | -0.5357 | 1 |
| Branch 2 | Bus 2 | Bus 3 | 0.1989 | -0.1989 | 1 |
| Branch 3 | Bus 3 | Bus 4 | -0.2791 | 0.2791 | 1 |
| Branch 4 | Bus 4 | Bus 5 | 0.1049 | -0.1049 | 1 |
| Branch 5 | Bus 5 | Bus 6 | 0.0239 | -0.0239 | 1 |
| Branch 6 | Bus 1 | Bus 6 | 0.1824 | -0.1824 | 1 |
| Branch 7 | Bus 2 | Bus 6 | 0.0886 | -0.0886 | 1 |
| Branch 8 | Bus 5 | Bus 2 | -0.0311 | 0.0311 | 1 |
|-----------------------------------------------------------------------|Thus, two active power flow measurements are added:
addWattmeter!(system, device; from = "Branch 1", active = powerFlow.power.from.active[1])
addWattmeter!(system, device; from = "Branch 4", active = powerFlow.power.from.active[4])Here, noise is not set, keeping the measurement values exact.
Active Power Measurements
Finally, the complete set of measurements is displayed:
printWattmeterData(system, device)|-----------------------------------------------|
| Wattmeter Data |
|-----------------------------------------------|
| Label | Active Power |
| | |
| | Measurement | Variance | Status |
| | [pu] | [pu] | |
|-------------|-------------|----------|--------|
| Wattmeter 1 | 0.7117 | 1.00e-04 | 1 |
| Wattmeter 2 | -0.2129 | 1.00e-04 | 1 |
| Wattmeter 3 | -0.4806 | 1.00e-04 | 1 |
| Wattmeter 4 | 0.3862 | 1.00e-04 | 1 |
| Wattmeter 5 | -0.1281 | 1.00e-04 | 1 |
| Wattmeter 6 | -0.2893 | 1.00e-04 | 1 |
| Wattmeter 7 | 0.5357 | 1.00e-04 | 1 |
| Wattmeter 8 | 0.1049 | 1.00e-04 | 1 |
|-----------------------------------------------|Figure 2 illustrates this measurement configuration, which includes active power injection measurements at all buses and two active power flow measurements.
Figure 2: The 6-bus power system with active power measurement configuration.
Base Case Analysis
After obtaining the measurements, the DC state estimation model is created:
analysis = dcStateEstimation(system, device)Next, the model is solved to determine the WLS estimator for bus voltage angles, and the results are used to compute power values:
solve!(system, analysis)
power!(system, analysis)This allows users to observe the estimated bus voltages along with the corresponding power values:
printBusData(system, analysis)|---------------------------------------------------------------------|
| Bus Data |
|---------------------------------------------------------------------|
| Label | Voltage | Power Generation | Power Demand | Power Injection |
| | | | | |
| Bus | Angle | Active | Active | Active |
| | [deg] | [pu] | [pu] | [pu] |
|-------|---------|------------------|--------------|-----------------|
| Bus 1 | 0.0000 | 0.7157 | 0.0000 | 0.7157 |
| Bus 2 | -1.5275 | 0.0059 | 0.2170 | -0.2111 |
| Bus 3 | -4.1606 | -0.0014 | 0.4780 | -0.4794 |
| Bus 4 | -1.1171 | 0.4629 | 0.0760 | 0.3869 |
| Bus 5 | -2.1624 | -0.0134 | 0.1120 | -0.1254 |
| Bus 6 | -2.1956 | 0.0084 | 0.2950 | -0.2866 |
|---------------------------------------------------------------------|Additionally, data related to measurement devices can be examined:
printWattmeterData(system, device, analysis)|---------------------------------------------------------------------|
| Wattmeter Data |
|---------------------------------------------------------------------|
| Label | Active Power |
| | |
| | Measurement | Variance | Estimate | Residual | Status |
| | [pu] | [pu] | [pu] | [pu] | |
|-------------|-------------|----------|----------|----------|--------|
| Wattmeter 1 | 0.7117 | 1.00e-04 | 0.7157 | -0.0040 | 1 |
| Wattmeter 2 | -0.2129 | 1.00e-04 | -0.2111 | -0.0018 | 1 |
| Wattmeter 3 | -0.4806 | 1.00e-04 | -0.4794 | -0.0012 | 1 |
| Wattmeter 4 | 0.3862 | 1.00e-04 | 0.3869 | -0.0007 | 1 |
| Wattmeter 5 | -0.1281 | 1.00e-04 | -0.1254 | -0.0027 | 1 |
| Wattmeter 6 | -0.2893 | 1.00e-04 | -0.2866 | -0.0026 | 1 |
| Wattmeter 7 | 0.5357 | 1.00e-04 | 0.5332 | 0.0025 | 1 |
| Wattmeter 8 | 0.1049 | 1.00e-04 | 0.1073 | -0.0024 | 1 |
|---------------------------------------------------------------------|Modifying Measurement Data
Let us now modify the measurement values. Instead of recreating the measurement set and the DC state estimation model from scratch, both are updated simultaneously:
updateWattmeter!(system, device, analysis; label = "Wattmeter 7", active = 1.1)
updateWattmeter!(system, device, analysis; label = "Wattmeter 8", active = 1.6)By changing these measurement values, two outliers are introduced into the dataset, which affects the estimated values.
Next, the DC state estimation is solved again to compute the updated estimate:
solve!(system, analysis)
power!(system, analysis)Bus-related data can now be examined:
printBusData(system, analysis)|---------------------------------------------------------------------|
| Bus Data |
|---------------------------------------------------------------------|
| Label | Voltage | Power Generation | Power Demand | Power Injection |
| | | | | |
| Bus | Angle | Active | Active | Active |
| | [deg] | [pu] | [pu] | [pu] |
|-------|---------|------------------|--------------|-----------------|
| Bus 1 | 0.0000 | 0.8554 | 0.0000 | 0.8554 |
| Bus 2 | -1.7703 | -0.2454 | 0.2170 | -0.4624 |
| Bus 3 | -1.3733 | 0.1520 | 0.4780 | -0.3260 |
| Bus 4 | 2.5037 | 0.9506 | 0.0760 | 0.8746 |
| Bus 5 | -2.5460 | -0.3107 | 0.1120 | -0.4227 |
| Bus 6 | -2.8567 | -0.2239 | 0.2950 | -0.5189 |
|---------------------------------------------------------------------|With the modified measurement values for Wattmeter 7 and Wattmeter 8, the estimated results deviate more significantly from the exact values obtained through DC optimal power flow, as the altered measurements no longer align with their corresponding values.
Now, instead of using the WLS estimator, we compute the LAV estimator:
analysis = dcLavStateEstimation(system, device, Ipopt.Optimizer)
solve!(system, analysis)
power!(system, analysis)Bus-related data can be examined:
printBusData(system, analysis)|---------------------------------------------------------------------|
| Bus Data |
|---------------------------------------------------------------------|
| Label | Voltage | Power Generation | Power Demand | Power Injection |
| | | | | |
| Bus | Angle | Active | Active | Active |
| | [deg] | [pu] | [pu] | [pu] |
|-------|---------|------------------|--------------|-----------------|
| Bus 1 | 0.0000 | 0.7246 | 0.0000 | 0.7246 |
| Bus 2 | -1.5454 | 0.0041 | 0.2170 | -0.2129 |
| Bus 3 | -4.1981 | -0.0026 | 0.4780 | -0.4806 |
| Bus 4 | -1.1579 | 0.4622 | 0.0760 | 0.3862 |
| Bus 5 | -2.1996 | -0.0161 | 0.1120 | -0.1281 |
| Bus 6 | -2.2281 | 0.0057 | 0.2950 | -0.2893 |
|---------------------------------------------------------------------|As observed, the estimates obtained using the LAV method are closer to the exact values from the DC optimal power flow, as LAV is more robust to outliers compared to WLS.
Modifying Measurement Set
Let us proceed with the LAV state estimation model and set two measurements to out-of-service:
updateWattmeter!(system, device, analysis; label = "Wattmeter 1", status = 0)
updateWattmeter!(system, device, analysis; label = "Wattmeter 5", status = 0)Recompute the LAV estimator and active power values:
solve!(system, analysis)
power!(system, analysis)Bus-related data can now be examined:
printBusData(system, analysis)|---------------------------------------------------------------------|
| Bus Data |
|---------------------------------------------------------------------|
| Label | Voltage | Power Generation | Power Demand | Power Injection |
| | | | | |
| Bus | Angle | Active | Active | Active |
| | [deg] | [pu] | [pu] | [pu] |
|-------|---------|------------------|--------------|-----------------|
| Bus 1 | 0.0000 | 1.5930 | 0.0000 | 1.5930 |
| Bus 2 | -3.1513 | 0.0041 | 0.2170 | -0.2129 |
| Bus 3 | -7.1481 | -0.0026 | 0.4780 | -0.4806 |
| Bus 4 | -5.2183 | 0.4622 | 0.0760 | 0.3862 |
| Bus 5 | -7.2534 | -0.8844 | 0.1120 | -0.9964 |
| Bus 6 | -5.9313 | 0.0057 | 0.2950 | -0.2893 |
|---------------------------------------------------------------------|As observed, while the LAV approach is more robust than WLS in handling outliers, the accuracy of the estimated values still depends on factors such as the magnitude of outliers, their number, and the positioning of meters within the power system. Removing two accurate measurements while keeping outliers in the system shows that even the LAV method cannot fully compensate for the loss of reliable data, leading to less accurate estimates.