A rigorous matrix procedure for calculating the line constants and wave parameters of uniform MTLs using SMT/PMU

Summary Real ‐ time monitoring of the operation state of wide ‐ area transmission line links has become possible with the help of phasor meter units. Synchronized information acquired by phasor meter units needs to be adequately processed to permit the accurate estimation of the line constants of the transmission link. In this paper, a novel general rigorous compact procedure for correctly processing the measured voltage and current phasors of uniform multiconductor transmission line systems is proposed. The procedure based on the ABCD matrix and on modal analysis techniques applies to transposed or untransposed multiconductor transmission lines, with arbi-trary geometry, number of conductors, and length. The proposed algorithm, adequate for multiport structures, avoids the approximations usually required by ordinary methods mostly focused on lumped parameters and on 2 ‐ port approaches. The proposed matrix procedure is illustrated and validated using simulation results.


| INTRODUCTION
Smart grids digital communication technology is able to provide accurate real-time measurement of voltages and currents at the ends of multiconductor transmission line (MTL) systems, even when the ends of the line are many hundreds of miles apart. This became possible with the advent of synchrophasor measurement technology (SMT) based on the utilization of phasor measurement units (PMUs), firstly developed in the 1980s by Arun Phadke and collaborators. 1 A historical perspective of SMT/PMU can be found elsewhere. 2,3 Very recent review papers on the subject, covering state-of-the-art applications and offering exhaustive bibliographic references, are available. 4,5 Among the various applications of SMT/PMU, 6,7 we are interested here in the important topic of remotely monitoring the operation state of high-or extrahigh-voltage long overhead power lines (MTL systems). This can be done by resorting to synchronized measurements of MTL voltage and current phasors, enabling the estimation of the power line constants (per unit length impedance Z and per unit length admittance Y). Once the computed estimates are compared with reference values stored in a database the state of operation of the MTL system can be checked, allowing planners to make judicious decisions.
Although not new, the topic of the estimation of line constants using SMT/PMU data has been attracting more and more attention. [8][9][10][11][12][13][14][15][16][17][18][19] In particular, the work in Lowe 19 provides a very good review and discussion of the main methods used to estimate the line constants of perfectly transposed and untransposed 3-phase lines, with emphasis on the untransposed case. However, the analyses or approaches used in the past are, in our view, not entirely satisfactory. On the one hand, they use questionable assumptions and approximations, and on the other, they do not make full use of the properties of the physical system under observation. Some flaws encountered in the literature are 1. Zero, positive and negative sequence modes are assumed to decouple the MTL into equivalent 2-port networks on which the analysis is based. The problem is that none of those "modes" exist; they do not coincide with actual traveling wave modes. 2. The referred equivalent 2-port networks are often given the shape of pi-circuits where the horizontal branch is related to Z and the vertical branches related to Y/2. The problem is that such circuits can only be used in transmission lines of very small length. 3. In some cases, where long MTLs are tackled, the analysis becomes too complicated, involving the solution of unnecessary nonlinear equations; 4. Most of the solutions proposed use lumped parameter approaches, avoiding the more physically adequate distributed parameter approach and associated modal analysis techniques; 5. For untransposed lines, the set of math equations required to describe the solution of the problem is unnecessarily heavy and quite cumbersome; 6. The particular physical properties of the MTL system (reciprocity and symmetry), which can be used to simplify the solution of parameters estimation, are not taken in account.
In this work, we establish a novel general rigorous compact matrix method that, using PMU information, permits the evaluation not only of the transmission line constants (p.u.l. impedance and admittance matrices) but also of the line wave parameters (characteristic wave impedance matrix, propagation constants, and mode eigenvectors). The method applies to uniform MTLs, with any number of conductors, any geometry, single or double circuit, of any length, for any frequency.
We will show that a minimum of n independent samples of the voltage and current phasors at the sending and receiving ends of an MTL with n+1 conductors suffice to obtain all its line constants. The only assumption required is that the MTL system can be considered uniform. Otherwise, if the system were considered nonuniform, the definition/concept of "line constants" itself would make no sense.
One can argue that MTLs are indeed nonuniform systems. However, the chain connection of many line segments (span lengths) acts in a way to homogenize the disturbing effect of the local nonuniformities occurring along the entire line length. This is so because, for 50 or 60 Hz, the size of the perturbations is negligibly small at the scale of the wavelength (6 or 5 Mm). As an elucidative example, take the case of overhead power line analysis considering the nonuniform problem arising from the effect of conductors sagging between towers. 20 When high frequency regimes are studied, the sag effect is important and manifests itself through resonance phenomena; however, when 50 or 60 Hz frequencies are considered, the sag effect produces no noticeable consequences; the analysis can be conducted successfully considering that the system is uniform, assigning an average height to each conductor. 20 This paper is organized into 6 sections, the first of which is introductory. Section 2 is dedicated to background material on the chain matrix description of MTL systems. The proposed general procedure for the evaluation of the chain matrix, based on PMU information, is presented in Section 3. The determination of the MTL line constants, using modal analysis, based on the computed chain matrix, is addressed in Section 4. The effectiveness of the developed method is validated in Section 5 using numerical simulation results. At last, conclusions are outlined in Section 6.

| THE CHAIN MATRIX
Consider a uniform MTL section with n + 1 conductors, one of them being the ground return (reference conductor #0).
As shown in Figure 1, the MTL connects two n-port networks, the generator side and the load side, both enforcing the MTL boundary conditions at x = 0 (sending end) and x = l (receiving end). Due, mainly, to minute changes in the load network, the MTL boundary conditions can fluctuate over time. Voltage and current phasors are acquired, from synchronized measurements, at the sending (S) and receiving (R) ends of the system.
The purpose of this research work is the presentation of a general compact procedure for the determination of the following line constants: • per unit length series impedance matrix, Z = R + jX and • per unit length shunt admittance matrix, Y = G + jS. • characteristic wave impedance matrix Z w , • modal propagation constants γ = α + jβ, and • modal voltage eigenvectors t (transformation matrix).
Since the MTL linear system is a 2n-port network it can be adequately described in the frequency domain by its chain matrix or ABCD matrix, from which all of the above mentioned line constants and wave parameters can be computed. Therefore, the key problem to be solved is the calculation of the entries of the A, B, C, and D matrices from PMU measurements at the MTL ports. Note, however, that the chain matrix does not depend on the particular boundary conditions of the MTL; the chain matrix is an intrinsic property of the MTL system.
Most of the research on the subject of the estimation of line constants using PMU does not take into account the physical properties of the MTL and the corresponding mathematical properties of the ABCD matrix. If that is done, several simplifications can be used, and a lot of work can be saved.
The chain matrix is a 2n × 2n matrix that relates voltages and currents phasors at the sending port with those at the receiving port: The vector phasor quantities U S , I S , U R , and I R are determined by the boundary conditions of the MTL system, which vary over time.
The MTLs are 2n-port networks containing linear passive nonanisotropic components; therefore, they are reciprocal. Consequently, the following properties must be obeyed: 21 where superscript T stands for transposition and 1 is the identity matrix of order n.
Moreover, for uniform MTLs, the network is symmetrical (the S and R ports can be interchanged) and, in this case, the following properties add on: 21 Substituting Equation 6 into Equations 3, 4, and 5 leads to According to the above properties, the unimodular chain matrix can be fully described by the pair A and B, from where the remaining matrices C and D may be obtained through

| EVALUATION OF THE A AND B MATRICES
The evaluation of A = [a jk ] and B = [b jk ] requires only 2 matrix equations. From Equation 1, and from Equations 1 and 8, 3.1 | Case n = 1 For n = 1, all the matrices and vectors in Equations 10 and 11 are simple scalars. The solution of Equations 10 and 11 is trivial: From Equations 12 and 9, all the entries of the chain matrix can be found through which complies with the general properties in Equation 6.

| General case
Things are more intricate when A and B are not scalars but matrices. For generic n, the A and B matrices have n 2 entries; therefore, 2n 2 scalar complex equations are required for their evaluation. B being symmetric would permit the elimination of n(n − 1)/2 equations, but we will not take advantage of that possibility. Reasons are twofold. On the one hand, that would complicate the computational algorithm, that is, more steps would be required to get lesser final equations to be solved (lesser unknowns) but that would break the formal simplicity and compactness of the proposed matrix method, without gained benefit for small n, as is the case of 3-phase power lines. On the other hand, the computed results of matrix B will allow, at an intermediate stage, to control/ check its symmetry-a strong asymmetry of B will be a symptom that the MTL under concern clearly deviates from a uniform system, precluding, in that case, reliable conclusions about the estimated line constants.
For a 2n-port MTL, the results in Equations 10 and 11 provide 2n simultaneous scalar equations, but we need 2n 2 scalar equations to find the matrix entries of A and B. This means that Equations 10 and 11 must be used not once but n times to allow a complete formulation of the problem. A set of n independent phasor measurements are required to be performed at n time instants.
The time gap between measurements should be long enough to permit significant changes of the MTL boundary conditions at one or both ports (for example, considering different load conditions at the receiving port). The meaning of "significant changes of the boundary conditions" will be discussed in Section 5.
To understand the procedure for the general case, we start by introducing the column vectors with all the voltages and currents phasors acquired in measurement (k), with k running from 1 to n: ; ; which, according to Equation 10, must satisfy the matrix equation: Likewise, according to Equation 11, we have Next, 4 square matrices n × n, gathering all the measured quantities, are defined: and by using them, the results in Equations 15 and 16 can be recast into a single equation: which, finally, allows the determination of the A and B matrices:

| EVALUATION OF MTL PARAMETERS
The key matrices A and B belonging to the chain matrix are extracted according to the general compact procedure established in Section 3.2, with the final result given in Equation 19. However, although not required, individual solutions for A and for B may be found from Equation 19. In fact, after some matrix algebra, one may write: At this point, 2 remarks are convenient. Firstly, one may notice that in the limit case n = 1, all of the intervening matrices in Equation 20 turn into scalars, yielding which confirms the conclusion found in Equation 13.
As a second remark, for the general case n > 1, one may notice that the direct application of the result in Equation 20, obtained using partitioned matrix techniques, requires that I SS , U RR , N 1A , and N 1B are nonsingular matrices. However, it must be underlined that these conditions, although sufficient, are not necessary for calculating [A B] from Equation 19. To make clear that I SS and U RR can be singular and the problem can still be solved, we can go back to Equation 18 and, using Pauli matrices, rewrite the equation as from where we get Using again partitioned matrix techniques, individual solutions for A and B, similar to Equation 20, but with different shape, can be found: In this case, the direct application of Equation 23 would not require the nonsingularity of I SS and U RR but would require the nonsingularity of I RR and U SS .
This shows that the calculation of A and B, from In addition, it can be noted from the results in Section 4.1 (Equations 24 and 28) that, for actual MTLs with losses, one can be sure that nonsingular A and B matrices do always exist, since none of their eigenvalues can be 0.
Is the process fool proof? No, of course. If, for example, the set of enforced sending end voltages is kept constant during measurements, and the load impedance at the receiving end does not vary when 2 or more of the measurements are executed, then all the submatrices I RR , I SS , U SS , and U RR turn out to be singular, and the problem cannot be solved. If the load impedance varies only very slightly then the problem can be solved but with large errors. In short, the success of the method relies on the degree of linear independence of the set of n measurements-see discussion in Section 5.2.

| Wave propagation parameters
From the theory of modal analysis of uniform MTL, 22,23 we know that where γ is a diagonal matrix gathering the modal propagation constants and T is the modal transformation matrix, whose columns provide the voltages distribution among conductors for each and every traveling mode, Therefore, the diagonalization of the A matrix, gives, at once, all the modal voltage eigenvectors t k and all the propagation constants: Also, from modal analysis theory, 22 we know that which permits the evaluation of the characteristic wave impedance matrix Z w of the MTL,

| Line constants
At last, the per unit length impedance Z and shunt admittance matrix Y of the MTL system can be estimated. From modal analysis theory, 22,23 we know that Z and Y can be directly obtained from γ and T in Equation 25, and from Z w in Equation 29, through

| SIMULATION RESULTS AND VA LIDATION
The methodology proposed in this paper has been put to the proof. In this section, numerical simulation results are provided, validating the theoretical procedure and showing the effectiveness of the method. The problem of the estimation of the line constants and wave parameters of an MTL system based on measurements made at the external accessible ends of the system is an inverse engineering problem. To validate the successful solution of the inverse problem addressed in the preceding sections, one needs to know beforehand the targeted answer to it.

| Characterization of the reference MTL
Here, we describe in detail the physical constitution of the MTL system under analysis and solve the much simple direct problem, and find numerically its line constants, its wave parameters, and its ABCD matrix. The reference MTL system is a 400 kV (1.4 GVA) 3phase line, in flat configuration, 500-km long, with 2 shield conductors. The cross section of the bilaterally symmetric MTL structure is depicted in Figure 2.
The soil is assigned a typical resistivity value of 100 Ωm (average soil).
Each phase conductor is a bundle of 2 identical stranded cables separated of 40 cm. The shield conductors are grounded at each tower along the line length; the tower-to-tower distance varying typically in the range 300 to 400 m. For the frequency under analysis (f = 50 Hz), the grounding of the shield conductors can be seen as a continuous effect, and therefore, in the MTL propagation equations, one can set shield conductor voltages equal to 0. This transforms the 5-conductor system into a reduced order 3-conductor system (n = 3).
The corresponding p.u.l. series impedance matrix Z and p.u.l. shunt admittance matrix Y were computed, at 50 Hz, using the MNLA software (see Appendix A), yielding yielding the set of modal characteristic wave impedances: and yielding the modal transformation matrix: The above results show, as expected, that the propagation modes (traveling waves) are split into 2 categories: symmetric and antisymmetric modes. Mode 1 is the symmetric aerial mode, mode 2 is the symmetric ground mode (the slowest one, with the highest attenuation), and mode 3 is the antisymmetric aerial mode. Now, we use the preceding numerical information to obtain the characteristic wave impedance matrix Z w in phase coordinates, 22 The remaining D and C submatrices of the chain matrix are evaluated using D = A T ,C = B −1 (A 2 − 1), from Equation 9.
The MTL boundary conditions can be defined with the help of the following Thévenin equations where Z S and U S0 are the generator-side Thévenin matrix parameters (impedance and open-end voltage); likewise, Z R and U R0 are the load-side Thévenin matrix parameters. Once the particular MTL boundary conditions are enforced, all the MTL voltages and currents can be found. Combining Equations 1 and 33, we get

| Application and validation of the method
The procedure described in this work must be able to accurately estimate the matrix line constants of the MTL, no matter its particular boundary conditions at x = 0 and at x = l. The only requirement is that n (n = 3) independent measures of the MTL voltage and current phasors are available for processing. As has been referred, the MTL intrinsic parameters do not physically depend on the line terminations (boundary conditions); however, they are computed from the latter. From a numerical point of view (accuracy), some line terminations may be easier than others; in this paper, where the effectiveness of the proposed method is to be checked, we will focus attention on a particular situation perceived as numerically "difficult": the case when the sending end voltage phasors remain the same during measurements and the receiving end termination is close to perfect matching, that is, In fact, in this situation, we have In addition, when Z R is close to Z w , the absence of reflected waves makes the receiving end voltages to be almost proportional to those at the sending end, and therefore, with matrix U RR near to singular (recall that in a lossless matched MTL, the voltages U R and U S would be equal except for a rotation factor, dependent on the line length). According to Equation 20, the quasi-singularity of U RR could originate numerical problems in the computation of B.
Since we do not have access to real PMU data, we simulated them using the equation set in Equation 34, together with the data in Section 5.1. The 3 independent "measurements" were obtained considering in Equation 35 that where the delta parameter δ is small and will be chosen for discussion. For illustration purposes, we offer below the computed entries of U R , I R , and I S (in polar coordinates) that simulate the set of 3 measurements for the case of a variation of 0.1% in the load impedance, ie, δ = 10 −3 .
For k = 1, relative errors in Equation 38 are plotted against the delta parameter in Figure 3B.
Obtained results show that • admittance errors are larger than those in the impedance; • errors in the mutual elements are larger than those in the self, diagonal elements; • for δ larger than 10 −2 , the relative errors are negligibly small, of the order of 1 part per million. In addition, the errors are fairly constant and stable; • when δ decreases from 10 −2 to 10 −4 , the relative errors increase steadily. For δ = 10 −4 , the impedance error is about 10 −3 whereas the admittance error exceeds 10 −2 .
In addition, the errors are rather unstable.
The reason why admittance errors are larger than those in the impedance may be due to the estimated values of Y reveal the undesired presence of a contribution with real part, which should not exist, since the reference MTL does not have transverse losses.
The reason why the errors in the mutual elements exceed those in the self-elements may be due to the definition of the relative error in Equations 37 and 38 where, for the mutual elements, the denominator has a smaller absolute value, since the mutual elements are smaller than the diagonal elements.
The apparent instability of the computed errors, for δ < 10 −2 , is especially remarkable in the logarithmic representation that we used; it would not be visible in a linear scale. It is a mere reflex of the numerical noise occurring in the computation process.
In conclusion, even in the case of very small variations (0.01%) of the load condition, the proposed procedure can efficiently determine the matrix line constants of the MTL, with a relative error of 0.1% in the p.u.l. series impedance matrix and 1% in the shunt admittance matrix.
Finally, we need to refer to a recent work by Lowe,19 where the estimation of the line constants of an untransposed 3-phase line was presented. A rigorous comparison of results is unfortunately not possible. In this paper, we defined unambiguously the line load impedance for the set of three measurements, however, Lowe did not; he only says that "the impedance of the load is varied", with no specification. In any case, assuming that the measurements made by him and by us involve 3 clearly different loads, the agreement is obvious, Lowe reported 0% error while we report 1 part per million.

| CONCLUSION
A novel matrix procedure aimed at the estimation of the line constants and wave parameters of uniform MTLs based on PMU information was developed in this paper.
Contrary to other methods, the proposed procedure is rigorous, general, compact, and effective. Rigorous because it does not involve any simplifying approximations; it uses a distributed parameter approach and real traveling wave modes obtained via modal analysis techniques. General because it applies to uniform multiconductor lines with any number of conductors, any geometry (single and coupled double circuits) and any length. Compact because a single-matrix equation is used to estimate all the relevant entries of the chain matrix of the MTL from where the line constants are obtained. Effective because, for a MTL with n + 1 conductors, a minimum of n independent measurements is strictly required for the estimation of the line parameters with an error around 1 part per million (provided that the MTL boundary conditions change more than 1% between measurements).
As a continuation of this work, we suggest 2 paths, eventually linked. One is a sensitivity analysis considering the presence of noise in PMU data. Another is the consideration of oversampling and least squares methods for improved accuracy.