Compensated Load Flow Solutions for Distribution System State Estimation

: State estimation of distribution systems typically relies on measurement sets with very low redundancy levels. In this paper, this fact is exploited by ﬁrst solving a conventional load ﬂow, using exclusively a critical set of measurements, and then compensating the solution to account for the few redundant measurements available. This leads to a suboptimal but sufﬁciently accurate estimate. It is shown how the sparse triangular factorization of the load ﬂow Jacobian matrix can be fully exploited throughout the compensation-based procedure, preventing in this way the ill-conditioning associated with the gain matrix arising in the conventional least-squares formulation. Simulation results are provided for measurement conﬁgurations customarily found in distribution systems, showing the potential advantages of the proposed methodology.


Introduction
Advanced distribution management systems (ADMS) are expected to host very sophisticated operational tools, capable of tackling the numerous challenges posed by upcoming active distribution networks. It is well known that the State Estimator (SE) is the cornerstone on which any other ADMS application is based [1]. Using topological information, network parameters, and the available measurements and pseudo-measurements sets for each loading/generating scenarios, the SE determines the most likely state of the system [2].
Recent studies have carried out a careful revision of the state of the art on Distribution System State Estimation (DSSE) [3][4][5]. Those works analyze the most relevant challenges that utilities have to overcome in order to obtain all the benefits provided by SEs. The most troublesome lies in the need to increase the level of available measurements, most times confined to High Voltage (HV)/ Medium Voltage (MV) substations and a few more measurements scattered along the feeders (usually at MV/MV regulating/switching substations) [6,7] and secondary distribution transformers [8]. This lack of redundancy still continues today, due to the large extension of distribution networks and the huge cost that fully telemetering such large networks would imply. This limitation has been usually circumvented by resorting to so-called load allocation schemes, aimed at inferring critical pseudo-measurements at secondary substations (power flows injections) [9,10]. The load allocation problem has been substantially enriched by the advent of Advanced Measurement Infrastructure (AMI), including smart meters, as this equipment provides an indirect measurement of the complex powers demanded downstream secondary substations, possibly including an estimation of Low Voltage (LV) power losses. References [11][12][13], along with several others reviewed in [4], make use of these new information sources for the sake of improving the observability of distribution systems. Despite this fruitful line of research, mainly focused on improving the quality of pseudo measurements, and, despite the trend to increase the deployment of system monitoring devices, DSSE will continue to commonly face low redundancy measurement systems in the short term, and will remain so in the medium term [5,14].
Another relevant challenge is that associated with the capability of DSSEs to cope with several sources of ill-conditioning, usually linked to distribution systems, such as [2]: a large proportion of injection measurements, very short and long lines simultaneously ending at the same bus, and very large weighting factors used to enforce virtual measurements. This is in addition to the well-known natural ill-conditioning of the gain matrix arising in Weighted Least-Squares (WLS) estimators. The ill-conditioning of the gain matrix may cause both convergence problems and inaccurate estimates, as thoroughly discussed in [15][16][17]. Furthermore, many of the available measurements at these levels are currents rather than powers [18], leading to additional numerical issues typically deteriorating the performance of DSSEs (see chapter 9 of [2] for more details).
Mainly for those reasons, notable research efforts have been directed towards the design of alternative, WLS-based DSSE formulations, specifically tailored to distribution power networks with radial topology, some of them based on branch current state variables instead of the traditional bus voltages [19][20][21][22][23][24][25][26].
The current trend towards fully automated distribution systems, based on robust, yet low redundancy monitoring systems, encourages the search for ad hoc SE solutions, better suited to those boundary constraints. This paper is focused on developing a new procedure for DSSE characterized by its simplicity, while at the same time retaining the capability to cope with all factors encountered at these voltage levels: low redundancy of measurements, handling of current measurements and risk of ill-conditioning. The goal is to reach a good compromise between accuracy and numerical robustness. The basic idea after the proposed solution method, which is also the main contribution of this work, consists of using a simpler and robust power flow tool as a first step, and then compensating the preliminary load flow solution to account for the remaining redundant measurements, not yet used by the load flow solver. This leads to a suboptimal yet accurate enough solution, which can be easily implemented on any type of distribution network, balanced or unbalanced, radially operated or meshed, as long as the measurement redundancy remains relatively low.
The paper is organized as follows: the next section discusses the context in which DSSEs are typically applied, which motivates the proposed methodology. Sections 3 and 4 present the mathematical framework on which the compensation method is based and the resulting algorithm, respectively. Next, Section 5 presents simulation results obtained with several benchmark systems, followed by the concluding remarks.

Motivation
This paper focuses on, but is not restricted to European MV networks, which are typically three-phase balanced systems. Figure 1 illustrates the topology, main components, and most commonly available measurements for these distribution systems.
Three types of substations can be found at this level [27]: the HV/MV primary substation, switching substations, which can include MV/MV regulating transformers in some cases, and MV/LV secondary substations. MV feeders leave HV/MV substations to reach consumers and distributed generators connected to secondary substations through typically shorter LV feeders. Most times, MV systems are planned with a ring or network topology, but they are radially operated for simplicity of protections, even though a meshed operational topology is also possible [28,29].
Regarding the monitored electrical quantities, all readings from smart meters downstream of LV feeders are collected by data concentrators located at secondary substations, constituting a good set of pseudomeasurements for the power delivered on the LV side of transformers. In addition, some utilities are investing in monitoring strategically selected secondary substations, where the power flows at the LV side of transformers are directly measured, rather than being the sum of downstream smart meters [30]. It is most important also to duly consider the exact zero-injection measurements at those nodes where neither consumers nor distributed generators are connected to, which are typically the MV sides of secondary transformers. All these injections (directly or indirectly measured), along with the voltage measurement at the MV bus of the primary substation, constitute a critical set of measurements that are the common input data to run the load flow tool.   [18] can also be available and used by the SE. Anyway, unlike in transmission systems, enjoying redundancy levels easily exceeding 4, current distribution systems barely exceed 1.15 (i.e., only 15% redundant measurements).
Note that all the previously commented pseudo-measurements, exact null injections, and real-time measurements have very different weighting factors associated, which causes ill-conditioning problems in conventional WLS-based SEs. Other relevant sources of problems, as discussed above, are the co-existence of long feeders and very short sections incident at the same bus (for example, some underground/overhead line segments can be a few meters long in urban networks), or the high proportion of power injections in the set of available measurements.

Proposed Method: Theoretical Background
As discussed above, for the low redundancy levels typically arising in distribution systems, it makes sense to divide the entire set of measurements into two disjoint subsets, namely: z c , representing a critical subset of n measurements allowing the n components of the state vector to be computed unambiguously; z r , representing the remaining subset of r measurements, redundant by nature, where typically r < n (unlike in transmission systems). Both sets of measurements are related to the unknown state vector through nonlinear measurement functions, as follows: where c and r are the vectors of measurement errors, assumed Gaussian and independent, characterized by the diagonal covariance matrices R c and R r , respectively. Let us assume initially that the redundant measurements are missing. Then, we could resort to a load flow solver to obtain the state vector, x c , exactly satisfying the set of critical measurements, z c , i.e., By applying the Newton's method, this involves the repeated solution of the following linear system, until ∆x k is sufficiently small: where k represents the iteration counter, H c is the Jacobian of h c , and is the mismatch or residual vector, which becomes null at the solution. When both types of measurements (critical and redundant) are simultaneously considered, the optimal estimatex is obtained by repeatedly solving the well-known Normal Equations (NE) arising from the Gauss-Newton iterative method: is the residual vector of the redundant measurements, and W c = R −1 c , W r = R −1 r represent the customary weighting matrices, accounting for the respective measurement uncertainties. Note that the dependence of the Jacobian matrices on the current state, x k , has been dropped for simplicity of notation.
Let the NE be rearranged as follows: If we replace x k in Equation (7) by the fully converged load flow solution (x k = x c ), then ∆z k c ≈ 0, yielding: where ∆z rc = z r − h r (x c ) is the residual vector of the redundant measurements, computed at the load flow solution, x c , provided by the critical measurements. For convenience, we introduce at this point the following elements: • An auxiliary r × 1 vector: An auxiliary r × n matrix: Based on the above definitions, the system Equation (8) can be more compactly written as follows: Note that matrix K can be very efficiently obtained by rows (i.e., columns of K t ) through repeated forward/backward solutions with the sparse triangular (or LU) factors of H c (the load flow Jacobian available from the last iteration), namely by successively solving: for every row i of H r . Substituting ∆x k from Equation (9) into the definition of y leads to the following reduced system: and, rearranging: The above system involves a dense coefficient matrix, but its size is small (r × r), in accordance with the low redundancy levels available in distribution systems. Given a load flow solution, and the resulting ∆z rc , y can be obtained from Equation (12), which in turn allows computing ∆x k from Equation (9). This constitutes the main contribution of the paper.
Solving Equation (9) can be interpreted as performing an additional load flow iteration, in which the mismatch vector is non null as a consequence of the redundant measurements. Obviously, owing to nonlinearities and measurement noise, the updated (compensated) state vector will not correspond exactly to the optimal solutionx, obtained by successively iterating with Equation (5), but its accuracy may be sufficient for practical purposes, and in any case much better than that obtained when the redundant measurements are fully ignored. Note that, in the ideal or theoretical case in which all measurements were exact, both ∆z rc and y would be null, and there would be no need to compensate the load flow solution.
A major advantage of performing load flow compensated solutions is that they involve factorizing only the Jacobian matrix H c , which is several orders of magnitude better conditioned than the WLS gain matrix (H t W H). See [2] for a description of the main sources affecting the ill-conditioning of the NE.

Practical Implementation: Post-Compensating the Load Flow Solution
The following approximate procedure may provide sufficiently good solutions for practical purposes, particularly when the whole set of measurements can be rearranged so that z c contains the more accurate subset:

1.
Disregarding initially the set of redundant measurements, obtain x c by solving the load flow problem Equation (3). This involves repeatedly applying Equation (4). As a byproduct, the last Jacobian matrix H c and its LU factors are available.

2.
Compute the auxiliary matrix K and the residual vector ∆z rc .

4.
Solve the large sparse system Equation (9). Figure 2 shows a flow diagram to help the reader visualize the major steps of the proposed method, including the iterative process of the load flow problem Equation (4). It is worth noting that, unlike in a regular load flow problem, where practical solution adjustments should be duly considered, such as enforcement of reactive power limits, in this case, such subtleties are not needed.
In order to be able to run the load flow, which is the basic and most important hypothesis on which the proposed method is based, we need at least a pair of independent measurements for each bus. For the feeder head bus, we assume a voltage (V) measurement and take the angle as phase reference (θ = 0). For the remaining buses, we assume that both P and Q measurements are available (possibly in the form of null-injection constraints), but in the presence of a distributed generator we could alternatively have P and V measurements, much like in transmission systems. If any V measurement is missing, we can simply use V = 1 as a pseudo-measurement, or the V measurement of a nearby bus, if available. If a P measurement is rather missing, we can replace it by a forecasted value, as the loads typically follow a repetitive daily pattern (this fact is very well known and widely exploited in power systems, for security analysis, etc.). Finally, if a Q measurement is missing, we can use the typical power factor of loads (e.g., cos φ = 0.97), depending on whether they are residential, industrial, etc., to generate a reasonable pseudo-measurement.

Results
Experiments comparing the proposed method with the conventional WLS SE are presented and discussed in this section. Both algorithms have been programmed in Matlab R2017a (MathWorks, Natick, MA, USA) and run on a 64-bit Intel Core i7 processor (2.8 GHz, 8 GB of RAM) (Intel, Santa Clara, CA, USA).
For the different scenarios, the corresponding sets of measurements have been generated from a load flow solution (considered as the exact network state for comparison purposes), to which Gaussian noise compatible with the standard deviation of each measurement is added.
The scenarios are built by considering the types of measurements usually deployed at the distribution levels: voltage magnitude at the head busbar, plus injections of active and reactive powers at the remaining load nodes. These critical sets have been completed with a variable number of current flow measurements randomly distributed throughout the network branches. For a given number of redundant current measurements, 100 different simulations are carried out, each differing in the branches where those measurements are located.
The standard deviations associated with the measurement errors are 0.0025 pu for the voltage magnitudes and power injections, and 0.01 pu for the current measurements. Virtual measurements (at zero-injection buses), unless otherwise indicated, are modeled as very exact measurements with a standard deviation of 10 −5 .
The convergence threshold is 10 −4 for all iterative processes (load flow and conventional SE). Four networks ranging in size from 85 to over 3000 buses are tested in the sequel. The reader is referred to elsewhere regarding the detailed information fully characterizing those networks [31][32][33]. Table 1 summarizes the relevant parameters which can be helpful to understand the difficulties encountered when running the DSSE, namely: rated voltage, base power, total current flow at the head bus and branch impedances with maximum and minimum absolute magnitudes.

Convergence and Condition Number
The condition number evaluates the numerical ill-conditioning of a linear system of equations. In the poorly conditioned systems (i.e., those with a very large condition number), the numerical errors can impact the solution to the extent that the convergence rate is severely impaired. In extreme cases, the coefficient matrix cannot simply be factored, preventing a solution to be reached.
First, a large 3239-bus network [31], with 200 current flow measurements included (redundancy 1.03), is tested. This network is a structurally poorly-conditioned system, which means that its parameters are the source of ill-conditioning. This happens when very short branches (low impedances) coexist with long branches (large impedances). Table 2 shows the condition numbers arising for the different algorithms. When running a conventional state estimator with this network, ill-conditioning problems appear in the factorization of equation system (5): the condition number reaches a very large value 1.4 × 10 19 , which prevents from an estimated state to be computed in any of the scenarios tested.
However, when running the proposed compensated load flow solution, the condition number is drastically reduced, and the convergence is reached in all tested cases. A total of four iterations are needed by the load flow to converge, with a condition number of 10 10 . In the additional iteration, at the compensation stage Equation (12), the condition number is 6.8 × 10 3 .
Another factor that significantly affects the condition number is the weights assigned to the zero-injection virtual measurements. Since those measurements are 'exact' (null injections), they are incorporated with a very low standard deviation, leading to a much larger weight than those assigned to the conventional measurements. This increases the numerical ill-conditioning of Equation (5).
Next, a real radial distribution network of 690 buses is tested. This network contains 44.7% of zero-injection buses [33]. The critical set of load flow measurements has been supplemented with 30 current flow measurements (redundancy 1.02). Three different standard deviations for the virtual (exact) measurements are simulated, namely σ zv = 10 −4 , σ zv = 10 −5 and σ zv = 10 −6 .
When running the conventional state estimator with the lowest weight (σ zv = 10 −4 ), all simulated scenarios converge, requiring an average of 3.9 iterations, with a condition number of 10 16 . However, when the weights of virtual measurements are increased (σ zv = 10 −5 ), convergence and ill-conditioning problems appear when solving Equation (5). In this case, when the algorithm converges, the average number of iterations increases to 8.7, and the condition number rises up to 10 18 . In 9% of the simulations, the WLS problem cannot be solved. If the standard deviation of the virtual measurements is further reduced to σ zv = 10 −6 , the ill-conditioning generalizes and none of the scenarios can be solved.
As opposed to the conventional WLS state estimator, the proposed compensated load flow converges for all the simulated scenarios, regardless of the weights associated with the virtual measurements. The load flow stage converges in four iterations, plus the final extra iteration of the compensation stage. Once again, the condition number of the equation systems gets significantly reduced. For the load flow stage, the condition number is 6 × 10 6 ; for the small equation system Equation (12) of the final stage, the condition number is just 81.
In addition to the networks tested in this section, Table 2 gathers the condition numbers of two smaller systems (85 and 117 buses) that will be tested in the following sections. For the latter networks, 10 current flow measurements have been considered when computing the condition number. As can be seen, in all cases, the proposed methodology gives rise to condition numbers significantly lower than those obtained with a conventional state estimator.

Accuracy Analysis
This section compares the accuracy of the compensated load flow solution with that of the conventional state estimator.
For this purpose, two MV distribution feeders are analyzed, namely an 85-bus network [32], and a 117-bus system [31]. In the 85-bus network, the number of current flow (redundant) measurements included in the different scenarios is increased from 1 to 30, whereas, in the 117-bus system, the redundant measurements are increased in 10 steps, each comprising 10 current measurements, up to 100 measurements. For each redundancy level, 1000 different tests (with random current measurement location) are carried out. Figure 3 shows the average absolute differences of the estimated voltage magnitudes and phase angles (in radians), with respect to the exact simulated values, for the 85-bus network scenarios. Figure 4 presents similar results for the 117-bus system. In order to show the improvement achieved by the correcting stage of the proposed algorithm, the average absolute errors arising after the first stage (just load flow computations based on critical sets) are shown with the red line.
The results clearly show the increased accuracy attained when more and more redundant current measurements are added, which somehow saturates after a threshold is exceeded. They also confirm that the suboptimal estimate obtained with the compensated solution is almost indistinguishable from the one provided by the conventional state estimator.

Computational Cost
In this section, the computational effort of both the conventional state estimator and the compensated load flow is compared. Table 3 presents the execution time per iteration required by each algorithm, as well as the number of iterations (in brackets). As mentioned before, the conventional state estimator could not reach a solution for the 3239-bus network due to ill conditioning problems. Further experiments have been performed for the 85-bus and the 117-bus systems. Regarding the convergence rate, Figure 5 shows the average number of iterations employed by a conventional state estimator, as the number of redundant current measurements increases. It can be observed that the higher the redundancy, the larger the number of iterations, reaching four iterations for the 85-bus network, and up to six iterations for the 117-bus system.
On the other hand, the extended load flow involves only three iterations for the load flow stage, plus the additional iteration compensating the load flow solution, in all simulated cases. It is worth stressing that, in the proposed method, the current measurements are incorporated at this last stage and, therefore, the number of iterations of the preliminary load flow solution is not affected by the redundancy level. Figure 6 shows the solution times required by both methodologies. As can be observed, there are no significant differences between both algorithms, regardless of the amount of current measurements considered. The compensated load flow method linearly increases the solution time with the number of redundant measurements. In any case, the total solution times are sufficiently small for this factor not to become the most relevant one in practice.

Incorporation of Additional Measurements
In previous sections, the scenarios considered have been generated by adding a set of current flow measurements to the critical set of measurements used by the load flow. The proposed methodology is, however, not limited to the case in which the redundant set is composed of current measurements only. In this section, scenarios incorporating also extra voltage measurements are analyzed.
For the 85-bus network, in addition to the critical measurements associated with the load flow problem, 10 voltage magnitudes and 10 current flow measurements are incorporated, randomly distributed throughout the network. One hundred different scenarios, with different random locations, are performed.
The conventional estimator shows a condition number of 7.2 × 10 14 , while the condition number of the compensated load flow keeps at 3.2 × 10 5 at the first stage, and 9.3 × 10 4 at the second. Both methodologies employ a similar number of iterations: an average of 3.95 iterations for the conventional state estimator, and four iterations for the proposed method. The solution time is also similar: 9.5 ms for the conventional state estimator versus 9 ms for the compensated load flow. The accuracy of both methodologies is also comparable, being the average absolute error 0.46 × 10 −4 for the voltage magnitudes, and 2.2 × 10 −5 radians for the phase angles.

Conclusions
This paper has examined the possibility of enhancing the results provided by a load flow solver, at minimum complexity and computational cost, with the information provided by a relatively small set of redundant measurements, a situation commonly arising in distribution systems. The goal is to attain estimates of similar accuracy as those provided by a conventional WLS estimator, while at the same time circumventing the severe ill-conditioning problems frequently found at those voltage levels. Simulation results have shown that the compensated load flow proposed in this work significantly reduces the condition number of the equation systems involved in the state estimation process. This way, the proposed methodology can successfully manage extremely ill-conditioned scenarios, where a conventional state estimator fails badly to obtain a solution. Moreover, the accuracy of the solution provided by the proposed algorithm has proved to be comparable to that of the conventional state estimator, while the computational cost is not an issue as long as the measurement redundancy remains low.