An algorithm for ﬁtting passive equivalent circuits for lumped parameter frequency dependent transmission line models

Accurately ﬁtting rational functions to the frequency response of modal impedances is crucial for including frequency dependency in lumped parameter models of transmission lines. Vector ﬁtting is widely used for ﬁtting rational functions to the frequency response of modal impedances and then an R-L equivalent circuit is obtained from the rational functions. A single-step method based on the properties of Foster equivalent circuit is proposed to directly ﬁt an R-L equivalent circuit to the frequency response of modal impedances. The positive peaks, negative peaks and the positive zero crossings of the slope change plot of the frequency response are found to provide a good approximation of zeros, poles and the ﬂat region locations of the frequency response. An enhanced ﬁtting algorithm is proposed based on these observations. A close enough ﬁtting is achieved using the proposed method with less number of passive elements. Using the proposed model in EMTP-RV, 400, 765, 1200 kV transmission lines and a 11-bus 500 kV network are simulated for switching transients. The results are compared with the constant parameter cascaded 𝜋 -model and the Frequency-dependent line model in EMTP-RV. The switching transient results of the proposed model are found to be


INTRODUCTION
High penetration of renewables through power converters has reduced the gap between the electro-mechanical dynamic transients and the electromagnetic transients to a negligible extent [1]. The shrinking difference between these two types of transients necessitates accurate modelling of individual components in power system which are valid for a wide frequency range (dc to several megahertz) [2]. In case of transmission lines, frequency dependent models gained importance due to intensification of harmonics and considerable distortion in wave shape of voltages and currents with the constant parameter (CP) models [3,4]. Especially frequency dependent models are extensively used in real time simulation [5] to deal with switching operations and system non-linearities for considerably long time [6]. Much efforts have been devoted for the development of frequency dependent models of transmission lines for digital computer simulations [3][4][5][6][7]. Among these models, Marti's model [6] [8,9] by representing line as a cascaded connection of multiple lumped-parameter -sections. This model is adopted in [10] to develop a lumped parameter frequency dependent transmission line (LPFD). Procedure followed in [10] is as follows, • The longitudinal parameters of each section are approximated by a rational function F ( ). • A lumped parameter R-L network is fitted to F ( ).
• The fitted longitudinal parameters are inserted in each section.
Even though [10] rightly represented frequency dependent line, it resulted in higher number of passive elements in fitted circuit which increased the order of state space matrix. Fitting an equivalent circuit for a rational function is addressed by several researchers in the literature [11][12][13][14][15]. Fitting the frequency dependent impedance to a higher order R-L circuit by vector fitting and Pade approximation via the Lanczos (PVL) process is proposed in [11,12]. An R-L circuit is fitted to the frequency response of a resistance in [13]. The functions characterized by smooth behaviour in frequency domain are fitted to a model, which gives a rational function with only real negative poles and real positive residues in [14]. Passive electrical equivalent circuits are derived from these rational functions.
A single step methodology is proposed in [15], which drastically reduces the number of passive elements required for fitting the frequency response of the characteristic impedance. This method is shown to produce same accuracy as that of vector fitting used in Marti's model of EMTP-RV. In [16], similar approach is used to fit an R-L equivalent circuit to the modal impedances of LPFD line proposed in [10]. This approach depends on the local minimization of fitting error by splitting the entire frequency range into multiple windows. Even though this approach achieves reduced number of elements, some amount of trial and error is involved in the initial values of the network elements and the number of windows (order of network). Fitting approach presented in [16] is used to obtain LPFD transmission line models of WECC 3 machine 9 bus system in [17]. This paper is an extended version of [16]. The algorithm is enhanced in this paper with a systematic selection procedure for the initial R-L values and order of the network based on the fundamental properties of the Foster R-L circuit [18]. It is found that the slope change plot of the frequency response provides a good approximation to zeros, poles and flat regions of the frequency response. Initial values of the circuit elements are decided based on the slope changes in the magnitude response. The proposed enhancements avoid the window approach and corresponding exhaustive search needed in [16]. The enhanced algorithm has been applied to single circuit 440 kV line with delta conductor configuration [19], single and double cirucit 765 kV line with delta and flat conductor configuration [20] and single circuit 1200 kV line with flat conductor configuration [20] transmission lines and the LPFD models are simulated in EMTP-RV. Switching transients in a 11-bus 4-generator 500 kV network are also simulated by replacing all the lines in the network by the proposed model [21]. The results are compared with the constant parameter cascaded -model and the standard frequency dependent line model available in EMTP-RV. The proposed LPFD line simulation results are found to be comparable to the EMTP-RV's Marti's model.

LUMPED PARAMETER FREQUENCY DEPENDENT TRANSMISSION LINE MODEL
In this section, mathematical formulation of the lumped parameter frequency dependent transmission line model proposed in [10] is described. The mathematical model of a single phase transmission line shown in Figure 1 with per unit length parameters is given by,  (1) and (2), one can get the voltage and current propagation equations for the line, In (5) and (6), the vectors V m ( ) and I m ( ) are transversal voltages and longitudinal currents in modal domain. Matrices Z m ( ) and Y m ( ) are per unit length longitudinal impedance and shunt admittance matrices respectively. The relationship between phase and modal quantities is given as follows, From (7) and (8), it is evident that the mode is completely decoupled from the other two modes. Whereas, and zero modes are mutually coupled with one another through Z 0 and Y 0 . These non-zero off-diagonal elements become insignificant when compared to the diagonal elements depending on the line characteristics [22]. Therefore, mutual terms Z 0 and Y 0 are discarded in (7) and (8). Then all the modes are declared as exact modes and transmission line representation in modal domain is obtained as shown in Figure 2.
Each mode is represented by the cascaded connection of multiple -circuits as shown in Figure 3 where Z ms represents the corresponding modal impedance (Z or Z or Z 0 ) [10,[23][24][25]. Each modal impedance is combination of resistance and inductance. Frequency response for these modal resistances and inductances is shown in Figure 4. It can be observed that after a certain frequency the inductance dominates the resistance. The resistance and inductance remain constant during low frequencies. Variation in the inductance of Modeand Modeis not significant. Only 2-5% variation of inductance from its initial value can be observed for all the voltage levels in these modes. Inductance variation in case of Mode-0 is quite significant and one can observe 8-90% variation from its initial value for all the voltage levels for this mode. However, resistance variation in all the modes is significant. Resistance decreases considerably with increase in frequency in all the modes. However, resistance variation in all the modes is significant. Resistance decreases considerably with increase in frequency in all the modes. For given a given AC input at one end of the transmission line, the effective impedance of the line seen by the source dictates the terminal characteristics at the other end of the line. Hence the modal impedances need to be fitted accurately to reproduce the terminal characteristics. Frequency response of the modal impedances is shown in Figure 5. Frequency dependency of modal impedances is incorporated by fitting a Foster R-L equivalent circuit shown in Figure 6 to the frequency response of the modal impedances.

ENHANCED ALGORITHM FOR FITTING PASSIVE EQUIVALENT CIRCUIT
A single-step algorithm is proposed in [16] to fit modal equivalent circuits to the frequency response of the modal impedances of the LPFD line model. In this section the fundamental ideas presented in [16] are repeated for the sake of completeness  Figure 6 is Z ms ( j ) in (10) can be represented as a sum of multiple first order transfer functions as follows: where, . G n (s) = sL n R n R n + sL n . A second-order equivalent circuit with R 0 = 0.5 Ω, R 1 = 5 Ω, R 2 = 50 Ω, L 0 = 5 mH , L 1 = 10 mH and L 2 = 0.5 mH is analyzed to study the impact of each element on the frequency response of Z ms . The values of passive elements are chosen such that their frequency response falls in the frequency range represented by windows in Figure 7. Impact of each element on the frequency response is shown in Figures 8-11 by varying the corresponding element while keeping other elements constant.
Frequency response of the Z ms is shown in Figure 8 by varying R 0 from 0.3Ω to 0.7Ω. From Figure 8, it is evident that the initial magnitude of Z ms depends only on the value of R 0 and it decreases/increases with decrease/increase in the value of R 0 [16]. Frequency response of the Z ms is shown in Figure 9 by varying L 1 from 8 mH to 12 mH. From Figure 8, it is evident that the slope of the frequency response in window 1 is a function of L 1 only and it decreases/increases with decrease/increase in the value of L 1 . Frequency response of the Z ms is shown in Figure 10 for the variation of R 1 from 3Ω to 7Ω. From Figure 10, it is evident that the slope of the frequency response in window 1 is a function of R 1 only and it decreases/increases with decrease/increase in the value of R 1 . Figure 11 shows the impact of L 2 & R 2 on the frequency response of Z ms , it can be observed that these passive elements corresponds to window 2 . Slope and height of frequency response in window 2 depends on L 2 and R 2 respectively. At higher frequencies, the magnitude of Z ms is increasing and it is not constant. It indicates that there is no resistance in parallel with an inductance in this window. So the inductance affecting the frequency response at higher frequencies is considered as L 0 and its corresponding window is labelled as window 0 . Analysis of Z ms magnitude variation concerning to each element in the circuit shows that the height of every flat zone depends on the value of the resistance and the slope of the magnitude variation depends on the inductance as shown in Figures 14 and 15. Based on the above observations, an algorithm is developed in [16] to ensure that at any frequency, the mismatch between the actual frequency response of a modal impedance (Z ms = Z or Z or Z 0 ) and the fitted Z fit ( ) remains within a user-defined error bound. This approach depends on the local minimization of fitting error by splitting the entire frequency range into multiple windows. Even though this approach gives reduced number of elements, when applied to transmission line models of different voltages levels, it is found that some amount of trial and error is required in the selection of initial values of the network elements and the number of windows (order of the network).
The following properties of the Foster R-L equivalent circuit [18] play an important role in the fitting process.
• All poles and zeros of the impedance are simple and located on negative real axis of s-plane. • Zeros and poles interlace. The lowest critical frequency is a zero which may be at s = 0. • Slope of the magnitude response is positive.
The magnitude response, its Slope and the change in slope (ΔSlope) of Z ms calculated using the following equations at the points on log scale are plotted in Figure 16. .
It can be found from Figure 16 that the slope is always positive and the value of the Slope, ΔSlope become zero during the flat regions. The positive zero crossings of the ΔSlope represent these flat regions and the magnitudes at those points become functions of R 0 , R 1 , R 2 . It can be observed that the ΔSlope plot at its negative peaks is very close to the pole locations of the Z ms , i.e.
respectively. Also, one can observe that the ΔSlope positive peaks are very close to the zero locations of the Z ms .
The positive zero crossings of the ΔSlope are obtained and then the corresponding magnitudes are approximated to the corresponding resistance values. After that the negative peaks of ΔSlope are identified and they are taken as the initial values of the pole locations. Around the vicinity of those initial locations, the corner frequency is chosen by minimizing the root mean square error (RMSE) between the actual and the fitted frequency responses. The order of the network is decided by the number of positive zero crossings/negative peaks. Absence of these indicate that the response can be fitted with just the series branch R 0 , L 0 . The series branch R 0 is obtained from the initial flat region. The positive peaks of the ΔSlope are obtained to get approximate zero locations of the frequency response. The value of L 0 is chosen to minimize the RMSE error in the frequency range between the smallest and the largest zero locations. It is found that these modifications significantly improved the convergence, avoids windowing and significant search efforts for optimal values in the algorithm of [16]. The enhanced fitting algorithm is described in Algorithm 1.

Analysis of modal impedances
The frequency response of modal impedances of 440, 500, 765 and 1200 kV transmission lines with different tower configurations [19,20,[26][27][28][29] Figure 18. It can be observed that the frequency at which the slope changes occur is different for different lines but the number of positive and negative peaks are same for all the lines. This indicates that the number of poles and zeros needed to fit the responses will be same for all the lines tested, but the locations of poles and zeros would be different. Figure 18a shows that the ΔSlope of Modefor all the lines has one positive peak and this indicates that the Modeimpedance has a single zero and no poles in its equivalent circuit. ΔSlope of Mode-0 shown in Figure 18b has two positive peaks and one negative peak which indicates that the Mode-0 ALGORITHM 1 Enhanced R-L circuit fitting algorithm 1: Calculate R 0 as R 0 = Z ms (1)

6:
Define frequencies corresponding to negative peaks as the corner frequencies of poles p1 , p2 , … ., pk 7: procedure Initialize G 1 (s), G 2 (s), .., G k (s) 8:  impedance has two zeros and one pole in its equivalent circuit. It is also evident from these responses that the modal impedances of all the transmission lines irrespective of their voltage levels and tower configurations can be fitted using simple Foster R-L equivalent circuits using the proposed algorithm.

Simulation of switching transients
The performance of the proposed algorithm is used for simulation of the switching transients of 440, 765 and 1200 kV transmission lines. Technical specifications of the transmission lines used for the study are summarized in Tables 1-3. The modal impedances and admittances of all the lines are obtained using (7) and (8) for a frequency range of 1 mHz to 100 MHz with 10 points/decade. Mutual terms in these matrices are neglected since they are insignificant [10]. In order to decide the number of data points required for the modal impedance   Figure 19 shows the Δslope variation with frequency for Mode-, Modeand Mode-0 for a 400kV line. Figure 19 shows that the change in slope is different in case of 5 points/decade and it is similar in all other cases. Fitting accuracy with 10 points/decade is improved in the proposed algorithm by minimizing the RMSE between the actual and the fitted frequency responses around the vicinity of those identified pole zero location. Foster R-L equivalents are fitted to the frequency response of the individual modal impedances using the proposed algo-   sufficient to fit these modal responses. Also, one can observe the presence of two positive peaks (one positive zero crossing) and one negative peak for all the lines in ΔSlope of Mode-0 responses in Figures 21c, 22 and 23c. This indicates that the Mode-0 impedance has two zeros and one pole. So, a first order system having two zeros (G 0 (s) and G 1 (s)) with R 0 , R 1 , L 0 and L 1 is fitted to the Mode-0 impedance of all the lines. The number of poles and zeros required for fitting the longitudinal parameters of the propagation modes are summarized in Table 4. The fitting results show that Modeand Modeequivalent circuits are similar to the equivalent circuits but the Mode-0 equivalent circuit has an additional R 1 and L 1 . The lumped parameter frequency dependent line model incorporate frequency dependency of the line parameters, a natural extension to the constant (CP) section model. In the following simulations, 10 km sections of modal impedances are considered for an improved simulation accuracy [30]. Fitted equivalent circuit parameters of each mode for all the lines are summarized in Table 5. Figure 20 shows the simulation setup modelled in EMTP-RV for analyzing the accuracy of the proposed algorithm. In case of switching transient studies, source is modelled as a voltage behind an equivalent impedance [31]. Assuming a short circuit ratio of 10 and X ∕R ratio of 20, source impedance of (0.005 +     Short circuit at the line-end After energization of the line at 0.5 ms, switch S 2 is closed at 100 ms to create a short circuit at the receiving-end. Switching transients in the current at the receiving-end are shown in Figure 25. Current at the receiving-end starts increasing immediately since the disturbance is at the receiving-end itself. Magnitude of this current increases in steps due to the short circuit transients. It can be observed from the figures that the CP-model reports very low short circuit currents compared to other two models. However, the FD-Proposed and FD-EMTP report more or less close short circuit currents. Voltage transients after short circuit removal Short circuit at the receiving-end is removed by opening switch S 2 at 200 ms and the results for voltage transients at the receiving-end are shown in Figure 26. First peak of the receiving-end voltage after the removal of the short circuit is very high. It is evident from Figure 26 that the peak value and the respective transients of the receiving-end voltage with the FD-Proposed closely follow the FD-EMTP line results. However, significant deviations from FD-EMTP can be observed in the CP-model responses.

Simulation of an interconnected network
A 11-bus 4-generator 500 kV interconnected power network shown in Figure 27 is used for testing the proposed algorithm [21]. This test system consists of seven 500 kV transmission lines of same type with the double circuit conductor configuration and the complete data of the test system is available in [21]. Foster R-L equivalent circuits are fitted to the modal  impedances of 500 kV line and fitting results are presented in Figure 28. It can be observed from Figure 28 that the Modeand Modeof 500 kV line also have a single zero in its equivalent circuit and Mode-0 has two zeros and one pole in its equivalent circuit. ΔSlope of all the modes of 500 kV line also follows the same trend as the modal impedances analyzed in the previous section. Fitted equivalent circuit parameters of each mode of the line are given in Table 6.  whereas the CP line voltage has significant deviation from the FD-EMTP line. After this, CB B is closed to integrate G4 to the remaining test system. A shot circuit is created at LD6 bus to evaluate the proposed line performance in short circuit conditions and the corresponding fault current and voltage at that bus after removal of short circuit are presented in Figure 29b and Figure 29c respectively. Magnitude of the fault current in case of FD-Proposed line is same as the FD-EMTP line as shown in Figure 29b. But, the CP-line model reports very low short circuit currents compared to the other two models. Similar to the voltage transients observed after the removal of short circuit in case of simple lines, first peak value and magnitude variation in the FD-Proposed line closely follow the FD-EMTP line as shown in Figure 29c. However, CP-model shows significant deviation from the FD-EMTP model.

CONCLUSION
This paper presents a simplified algorithm to fit Foster R-L equivalent circuit to the modal impedances of a lumped parameter frequency dependent (LPFD) transmission line. Unlike the previous approaches used for fitting where an equivalent circuit is obtained from the transfer function of the frequency response, the proposed approach directly fits an equivalent circuit to the modal impedance frequency responses. Number of passive elements required for each modal impedance is obtained from the negative peaks or positive zero crossings of the frequency response slope change plots. The series branch parameter R 0 is approximated from the initial flat region of the frequency response. The initial resistance parameters of the passive elements are approximated from the positive zero crossings and the initial inductances are approximated from the negative peaks. The root mean square error (RMSE) between the fitted and the actual responses is used as a measurement index in the selection of optimal values of the inductances.