Performance Analysis of Conjugate Gradient Algorithms Applied to the Neuro-Fuzzy Feedback Linearization-Based Adaptive Control Paradigm for Multiple HVDC Links in AC/DC Power System

Abstract

The existing literature predominantly concentrates on the utilization of the gradient descent algorithm for control systems’ design in power systems for stability enhancement. In this paper, various flavors of the Conjugate Gradient (CG) algorithm have been employed to design the online neuro-fuzzy linearization-based adaptive control strategy for Line Commutated Converters’ (LCC) High Voltage Direct Current (HVDC) links embedded in a multi-machine test power system. The conjugate gradient algorithms are evaluated based on the damping of electro-mechanical oscillatory modes using MATLAB/Simulink. The results validate that all of the conjugate gradient algorithms have outperformed the gradient descent optimization scheme and other conventional and non-conventional control schemes.

1. Introduction

In recent years, the rapid growth in energy consumption and the implementation of the deregulated environment have resulted in the highly interconnected and stressed infrastructure of the power system. A change in the operating condition at any one location can have a magnified impact over a wide area as it propagates through the electrical network. Under the stressed operating conditions, even a small disturbance can threaten the stability of the power system and generate power system oscillations. Among these, Low-Frequency Oscillations (LFOs) would limit the power transfer capabilities between different interconnected areas and have a great impact on the most secure and economical operation of the power systems [1].

Rapid and effective damping of LFOs is an important issue as the sustained LFOs may result in almost instantaneous large-scale cascading failures. To avoid the economic loss and daily life interruption, it is highly necessary to dampen LFOs to continue the secure operation of the power system [2,3]. To improve the stability of the power system, a solution might be to utilize controllable power components within the system. One such component is the HVDC link [4].

Over the last few decades, rapid growth has been observed in the application of the HVDC technology in electric power transmission systems due to its economical and technical advantages. Mainly, HVDC transmission systems are employed to transmit electric bulk power over long distances and interconnect asynchronous AC systems having different frequencies [5]. Available technology options for HVDC transmission system includes Line Commuted Converters (LCC) and forced commutated converters. The LCC-HVDC systems have evolved as a much more mature and reliable technology [6]. The simply structured classic LCC-HVDC systems are available in high power ratings with lower power losses [7]. As HVDC system operation involves no inertia, it is possible to rapidly control the power flow through the HVDC link. The instant variation in power flow can be used to improve system operating flexibility, as well as maintain the stability of the power system during perturbed operating conditions. The provision of appropriate supplementary damping control of the classical HVDC system provides modulation of power flow through the HVDC link and effectively mitigates inter-area oscillations of interconnected power systems [8]. Recently, extensive research has been carried out to propose Power Oscillation Damping (POD) controllers applying classical control theory, optimal control schemes, robust control methods and the feedback linearization technique [9,10,11,12,13,14,15,16]. Classical control theory-based POD controllers are designed with a first-order approximation of the system dynamics at a certain operating condition. This method provides adequate damping for a system operating in the vicinity of equilibrium. However, their performance may decline due to the wide-ranging and rapidly varying operating conditions of the highly nonlinear and stochastic nature of the power system [17]. The linear optimal control schemes are based on the system model, which is linearized at a specific operating point and suffers the same drawbacks and limitations as classical control [18]. Robust control improves the performance of classical POD control. However, it requires an upper bound on uncertainties and is deficient in learning. Its performance deteriorates in a highly nonlinear power system with increased uncertainties [19].

For highly nonlinear systems, differential geometry-based Feedback Linearization Control (FBLC) is a better option to model a nonlinear control scheme [20]. However, it requires having exact knowledge of the system nonlinearities. Model uncertainty can result in degraded performance and hence limited applications in the power system [21]. The POD controller realization based on an adaptive feedback linearization becomes more suitable for applications in the highly nonlinear power systems with unknown or uncertain parameters [22]. Although model-based methods for identification and control of plants dynamics are well established and widely used, they may not provide satisfactory performance in today’s power system. Artificial intelligence-based control techniques have been applied for systems with both un-modeled dynamics and uncertainty in the parameter [23].

In [24], Artificial Neural Network (ANN)-based FBLC is employed to improve power system stability by mitigating LFOs. The opaque nature of ANNs includes learning capability, and they are capable of parallel processing. ANNs entail an ample amount of problem data for their training and are prone to overfitting [25]. ANNs’ drawbacks are tackled by the integration of their merits and the Fuzzy Logic System (FLS) that results in a more intelligent Adaptive Neuro-Fuzzy Inference System (ANFIS). The combination of FLS and ANNs acquires the capabilities of robustness, parallel processing and learning capability from ANNs and the imprecise knowledge-based modeling capability of the simple and natural structure FLS. ANFIS is dynamic and capable of learning from experience with sample data to estimate a function without any mathematical model. The online optimization of the parameters of ANFIS enhances its fast and accurate learning capabilities to achieve the desired output [26,27].

The gradient-based back-propagation algorithm has been widely used as a learning algorithm for ANFIS [28]. The Steepest Descent (SD) algorithm follows the negative of the gradient to approach the desired local minimum in a zigzag path. Newton methods have faster convergence than the SD method, but require computation of the Hessian and its inverse. A prominent iterative Conjugate Gradient (CG) method computes new search directions by using only one previous search vector at each time step and does not follow a pre-specified direction of the negative of the gradient. The speed of convergence is better than the SD method as the CG algorithm moves in the non-interfering directions. Computation of the Hessian is not a prerequisite in CG algorithms. It requires less memory to operate and less computation as compared to Hessian-based methods [29]. With simple implementation and fast convergence, CG methods become numerically competitive for applications in a highly nonlinear power system [30,31]. Furthermore, a real-time adaptation of ANFIS parameters has a swift response in comparison with batch mode adaptation and rapidly optimizes the parameters of ANFIS [32]. The performance of FBLC is further improved through on-line tuning of its coefficients. Optimized coefficients keep the control magnitude bounded and effective over a wide range of operating conditions [33]. Simple implementation of the Least Mean Square (LMS) algorithm makes it the most prevalent adaptive algorithm. Anyhow, the LMS algorithm is sensitive to the scaling of its input and exhibits data-dependent slow convergence behavior. Its variant, normalized Least Mean Square (nLMS) algorithm, is as simple as LMS and more robust. The nLMS algorithm is more stable on average than LMS and shows an improved balance between simplicity and performance in real-time applications [34].

A MIMO Adaptive Neuro-Fuzzy Feedback Linearization Control (ANFFBLC) for multiple HVDC links has been proposed in [35]. The ANFIS architecture is employed for real-time identification of unknown nonlinear dynamics of the power system. ANFIS parameters are optimized through a nonlinear CG algorithm proposed by Hager and Zhang. FBLC coefficients are on-line tuned to track the power system output for an extensive operating conditions. On-line self-tuning of FBLC coefficients and real-time adaptive neuro-fuzzy identification aspects of ANFFBLC make it an appropriate POD control for HVDC systems. The main contributions of this research include:

• Implementation of the SD algorithm and six CG methods for online optimization of ANFFBLC parameters.
• Evaluation and comparison of damping performance of the SD algorithm and six CG methods along with conventional and non-conventional control schemes.

Different case studies are performed in a multi-area AC/DC test power system to validate the damping performance of ANFFBLC with different CG algorithms. Results show a promising response with ANFFBLC and improved stability of the power system as compared with conventional and non-conventional control schemes.

The rest of the paper is organized as follows. Section 2 describes the modeling and control of the AC/DC power system. The proposed closed-loop control system design is explained in Section 3. Results are presented and discussed in Section 4. Section 5 concludes the findings of this research work.

2. AC/DC Power System Model Description

The major components of test AC/DC power system comprise the synchronous generators with control, LCC-HVDC links, transmission lines and transformers. Synchronous machines and control dynamics are expressed by differential equations, while algebraic equations describe the load flow and network model.

2.1. Power System Components Modeling

The synchronous machines with their associated control constitutethe major components of the modern power system. This research work implements a sixth order synchronous machine, while considering the dynamics of the stator, field and damper windings. The sixth order model of the i-th synchronous machine is given by the following differential equations [36].

$$\begin{array}{ccc}\hfill T_{do}^{\prime }\frac{d}{dt}{E}_q^{\prime }& =& E_{fd}-E_q^{\prime }-(x_d-x_d^{\prime })\left\{I_d-\frac{x_d^{\prime }-x_d^{″}}{{(x_d^{\prime }-x_{ls})}^2}({\psi }_{1d}+(x_d^{\prime }-x_{ls})I_d-E_q^{\prime })\right\}\hfill \end{array}$$

(1)

$$\begin{array}{ccc}\hfill T_{do}^{″}\frac{d}{dt}{\psi }_{1d}& =& -{\psi }_{1d}+E_q^{\prime }+(x_d^{\prime }-x_{ls})I_d\hfill \end{array}$$

(2)

$$\begin{array}{ccc}\hfill T_{qo}^{\prime }\frac{d}{dt}{E}_d^{\prime }& =& -E_d^{\prime }+(x_q-x_q^{\prime })\left\{I_q-\frac{x_q^{\prime }-x_q^{″}}{{(x_q^{\prime }-x_{ls})}^2}({\psi }_{2q}+(x_q^{\prime }-x_{ls})I_q+E_d^{\prime })\right\}\hfill \end{array}$$

(3)

$$\begin{array}{ccc}\hfill T_{qo}^{″}\frac{d}{dt}{\psi }_{2q}& =& -{\psi }_{2q}-E_d^{\prime }+(x_q^{\prime }-x_{ls})I_q\hfill \end{array}$$

(4)

$$\begin{array}{ccc}\hfill \frac{d}{dt}\delta & =& {\omega }_r-{\omega }_s\hfill \end{array}$$

(5)

$$\begin{array}{ccc}\hfill \frac{d}{dt}{\omega }_r& =& \frac{1}{2H}\left(P_m-P_e-D\Delta {\omega }_r\right)\hfill \end{array}$$

(6)

where:

• $E_{fd}$ is the field voltage.
• $E_{d/q}^{\prime }$ is the transient voltage in the d/q-axis.
• $I_{d/q}$ is the armature current in the d/q-axis.
• $x_{d/q},x_{d/q}^{\prime },x_{d/q}^{″}$ are synchronous, transient and sub-transient reactances in the d/q-axis, respectively.
• $x_{ls}$ is the armature leakage reactance.
• $T_{do/qo}^{\prime }$ is the transient time constants of the d/q-axis.
• $T_{do/qo}^{″}$ is the sub-transient time constants of the d/q-axis.
• ${\psi }_{1d}$ is the d-axis flux linkage of the damper winding.
• ${\psi }_{2q}$ is the q-axis flux linkage of the damper winding.
• H is the inertia constant.
• $P_e$ is the generator electrical power.
• $P_m$ is the generator mechanical power.
• $\delta$ is the generator rotor angle.
• ${\omega }_r$ is the generator rotor speed.

The dynamics in the d-axis are narrated by Equations (1) and (2), while Equations (3) and (4) define the dynamics in the q-axis. The swing equation of the synchronous machine is represented by Equations (5) and (6). With negligible stator’s resistance, the stator voltage and machine’s real power are given as:

$$\begin{array}{ccc}\hfill E_d& =& \frac{x_q^{″}-x_{ls}}{x_q^{\prime }-x_{ls}}{E}_d^{\prime }-\frac{x_q^{\prime }-x_q^{″}}{x_q^{\prime }-x_{ls}}{\psi }_{2q}+x_q^{″}{I}_q\hfill \end{array}$$

(7)

$$\begin{array}{ccc}\hfill E_q& =& \frac{x_d^{″}-x_{ls}}{x_d^{\prime }-x_{ls}}{E}_q^{\prime }+\frac{x_d^{\prime }-x_d^{″}}{x_d^{\prime }-x_{ls}}{\psi }_{1d}+x_d^{″}{I}_d\hfill \end{array}$$

(8)

$$\begin{array}{ccc}\hfill E_t& =& \sqrt{E_d^2+E_q^2}\hfill \end{array}$$

(9)

$$\begin{array}{ccc}\hfill P_e& =& E_d{I}_d+E_q{I}_q\hfill \end{array}$$

(10)

where $E_t$ is the generator terminal voltage. The complete generator model includes a detailed nonlinear hydraulic turbine with a governor system. The excitation system is implemented with the IEEE Type 1 voltage regulator combined with an exciter. The long transmission line is modeled with the distributed parameter, and distribution lines are modeled with lumped parameters and $\pi$-representation. All transmission lines, transformers and loads are combined to model the power system network, and the nodal equations describe the relationships between the system voltages and currents.

2.2. HVDC Dynamics

LCC-HVDC is comprised of the commutation transformer, thyristor bridges, the AC filter system to limit the harmonic content in the AC current, the DC smoothing reactor, the mono-polar DC transmission line and the control system.

2.2.1. LCC-HVDC Converter Modeling

The LCC-HVDC system represented by an average model is shown in Figure 1. The variables are shown with appropriate subscripts I and R for the inverter and rectifier poles, respectively. The rectifier and inverter poles are average models of 12-pulse converters with equivalent DC voltage and AC current sources. The voltage source generates converter voltage at the DC side, and the current source injects the fundamental component of the current into the AC network. The linear model of the commutation transformer is also included in the average model representation. The model does not represent high-frequency switching harmonics, which are unnecessary for dynamic analysis and control dynamics. However, dynamics resulting from the control system and the power system interaction are preserved [37]. The rectifier and the inverter are interconnected through a DC transmission line with smoothing reactors for each pole.

The average DC voltage with firing delay angle, $\alpha$ is defined as:

$$\begin{array}{c}\hfill V_{dc}=\frac{6}{\pi }{\int }_{-\pi /6+\alpha }^{\pi /6+\alpha }\sqrt{2}a{V}_{t}cos(\omega t)d(\omega t)=\frac{6\sqrt{2}}{\pi }a{V}_{t}cos\alpha \end{array}$$

(11)

where $V_t$ is the converter AC bus voltage, a is the transformer’s off-nominal tap ratio and $\frac{6\sqrt{2}}{\pi }a{V}_t=V_{dco}$ is the ideal no-load DC voltage; where:

$$\begin{array}{c}\hfill \frac{6\sqrt{2}}{\pi }a{V}_t=V_{dco}\end{array}$$

(12)

The voltage source in the average model is powered by Equation (11). The figure also models the thyristor on-state voltage drop and losses in the converter and transformer. The diode represents the unidirectional current flow without any losses associated with it. The DC voltage at the rectifier and inverter, with the effect of commutation overlap, is given as:

$$\begin{array}{ccc}\hfill V_{dcR}& =& V_{dco}cos\alpha -2\frac{3X_c}{\pi }{I}_{dcR}\hfill \end{array}$$

(13)

$$\begin{array}{ccc}\hfill V_{dcI}& =& V_{dco}cos\gamma -2\frac{3X_c}{\pi }{I}_{dcI}\hfill \end{array}$$

(14)

In the above equations, $I_{dc}$ is the converter direct current, $X_c$ is the commutation reactance and $\gamma =\pi -\alpha -\mu$ is defined as the extinction advance angle. Here, $\mu$ is the commutation interval. The fundamental primary phase current injected into the AC line is computed in the $dq$-axis frame as:

$$\begin{array}{ccc}\hfill I_d& =& \frac{V_{dco}(cos2\alpha -cos(2\alpha +2\mu ))}{4\sqrt{6}a{X}_c}\hfill \end{array}$$

(15)

$$\begin{array}{ccc}\hfill I_q& =& -\frac{V_{dco}(2\mu +sin2\alpha -sin(2\alpha +2\mu ))}{4\sqrt{6}a{X}_c}\hfill \end{array}$$

(16)

$$\begin{array}{ccc}\hfill \mu & =& {cos}^{-1}\left(cos\alpha -\frac{2X_c{I}_d}{V_{dco}}\right)-\alpha \hfill \end{array}$$

(17)

The active and reactive power exchanged with the AC system is calculated as:

$$\begin{array}{c}\hfill P_{ac}=2\frac{\pi }{6}a{V}_{dco}{I}_d\end{array}$$

(18)

$$\begin{array}{c}\hfill Q_{ac}=2\frac{\pi }{6}a{V}_{dco}{I}_q\end{array}$$

(19)

2.2.2. DC Transmission System

The DC transmission system model includes a DC line with smoothing reactors at both pole terminals. The T-model for the DC line dynamics of the LCC-HVDC system is described by two reactor currents and a capacitor voltage. With DC voltages, $V_{dcR}$ and $V_{dcI}$, available at the terminals of DC system, the dynamics of DC transmission line are given as [16]:

$$\begin{array}{ccc}\hfill L_d\dot{I_{dcR}}& =& V_{dcR}-V_C-I_{dcR}\frac{R_d}{2}\hfill \end{array}$$

(20)

$$\begin{array}{ccc}\hfill L_d\dot{I_{dcI}}& =& V_C-V_{dI}-I_{dcI}\frac{R_d}{2}\hfill \end{array}$$

(21)

$$\begin{array}{ccc}\hfill C_d\dot{V_c}& =& (I_{dcR}-I_{dcI})\hfill \end{array}$$

(22)

where $R_d$, $C_d$ and $V_c$ are resistance, capacitance and voltage across the capacitance of the DC line, respectively.

2.2.3. LCC-HVDC Control

The control architecture of an LCC-HVDC system is implemented with internal control and external POD control. The rapid control of the desired power flow through HVDC link is realized through swift valve gate control. In the internal control, master control determines the appropriate current order $I_{ref}$ and transmits the same to the rectifier and inverter controls. The current order is calculated as:

$$\begin{array}{c}\hfill I_{ref}=\frac{P_{ref}}{V_{dcR}}+I_{damp}\end{array}$$

(23)

where $P_{ref}$ is the desired HVDC power and the damping signal $I_{damp}$ is derived through ANFFBLC. During normal operation, the power order determines the current order, and $I_{damp}$ is zero. Under disturbed operating conditions, the ANFFBLC generates the appropriate damping signal that adjusts the $I_{ref}$ and, hence, modulates HVDC power. The conventional Proportional-Integral (PI)-based pole control generates the respective ignition delay angle $\alpha$ subject to current reference and the Voltage-Dependent Current Order Limiter (VDCOL). In the case of low commutating voltage, the VDCOL redefines the reduced current reference according to reduced DC voltage once it becomes less than the predefined threshold [38].

In the normal operation, the rectifier maintains a constant current through the DC transmission system. The control mode is shifted to constant ignition angle with minimum delay angle ${\alpha }_{min}$ during disturbed operating conditions. The rectifier current control is defined as:

$$\begin{array}{c}\hfill {\dot{\alpha }}_R=K_P{\dot{I}}_{dcR}+K_I{I}_{eR}\end{array}$$

(24)

where the current error is $I_{eR}=I_{dcR}-I_{ref}$ and ${\alpha }_{R\_min}\le {\alpha }_R\le {\alpha }_{R\_max}$. The inverter is operated in a constant extinction angle control mode that maintains the DC system voltage at the inverter. For the reduced commutation voltage in distressed operating conditions, the inverter control switches to constant current control mode. The dynamics of the inverter current control, voltage control and gamma control are described respectively as:

$$\begin{array}{ccc}\hfill {\dot{\alpha }}_I^I& =& K_P{\dot{I}}_{dcI}+K_I{I}_{eI}\hfill \end{array}$$

(25)

$$\begin{array}{ccc}\hfill {\dot{\alpha }}_I^V& =& K_P{\dot{V}}_{dcI}+K_I{V}_{eI}\hfill \end{array}$$

(26)

$$\begin{array}{ccc}\hfill {\dot{\alpha }}_I^{\gamma }& =& K_P\dot{\gamma }+K_I{\gamma }_{eI}\hfill \end{array}$$

(27)

$$\begin{array}{ccc}\hfill {\alpha }_I{|}_{{\alpha }_{I\_min}}^{{\alpha }_{I\_max}}& =& min\{{\alpha }_I^I,{\alpha }_I^V,{\alpha }_I^{\gamma }\}\hfill \end{array}$$

(28)

where current, voltage and $\gamma$ error, respectively, are $I_{eI}=I_{dcI}-I_{ref}+I_m$, $V_{eI}=V_{dcI}-V_{dcref}-\left(V_{dcm}\frac{I_{ref}-I_{dcI}}{Im}\right)$, ${\gamma }_{eI}=\gamma -{\gamma }_{ref}-\left({\gamma }_m\frac{I_{ref}-I_{dcI}}{Im}\right)$. The inverter control also provides current support for the situation when the DC current drops below a threshold equal to $I_{ref}-I_m$. Furthermore, $\alpha$ for the inverter is maintained to ensure $\gamma \ge {\gamma }_{min}$ [39]. $I_m$, $V_m$ and ${\gamma }_m$ are current, voltage and extinction advance angle margins at the inverter pole, respectively.

3. Closed-Loop Control System Design

The proposed ANFFBLC scheme is a model-free indirect control strategy and requires minimal knowledge of the power system. As shown in Figure 2, the dynamics of the AC/DC power system are identified through an Adaptive Neuro-Fuzzy Identifier (ANFI) using the Wide-Area Measurement System (WAMS)-based measured speed signals of generators. The CG algorithm on-line optimizes the ANFI parameters to minimize the identification error. The coefficients of the FBLC are also tuned through the nLMS algorithm to maintain its performance over a wide range of operating conditions. Based on the identified power system model, the ANFFBLC generates the appropriate damping signal for each HVDC link. The following section explains the implementation and function of the ANFFBLC to improve the damping of the power system oscillations.

3.1. Feedback Linearization Control

Feedback linearization is a well-known technique used to transform a nonlinear system into a linear one by a state transformation and nonlinear feedback. A well developed linear control design now can be applied to a linearized system [40]. Consider an-th-order multi-variable nonlinear system [41]:

$$\begin{array}{c}\hfill \begin{array}{c}{y}_1^n=f_1(x)+g_{11}(x)u_1+\cdots +g_{1m}(x)u_m\hfill \\ \begin{array}{ccc}& ⋮& \end{array}\hfill \\ y_m^n=f_m(x)+g_{m1}(x)u_1+\cdots +g_{mm}(x)u_m\hfill \end{array}\end{array}$$

(29)

where $x={[x_1,\dots ,x_1^{n-1},\dots ,x_m,x_1^{n-1}]}^T\in {\Re }^{mn}$ is the state variable available for measurement, $u={[u_1,u_2,\dots ,u_m]}^T\in {\Re }^m$ represents the control input and $y={[y_1,y_2,\dots ,y_m]}^T\in {\Re }^m$ is the output of MIMO nonlinear system. $f_{\tau }(x)\in {\Re }^m$ and $g_{\tau h}(x)\in {\Re }^{m\times m}$ represent unknown smooth nonlinearities of the dynamic system with $\tau ,h=1,2,\dots ,m$. The system in Equation (29) can also be expressed in companion form as:

$$\begin{array}{c}\hfill \mathbf{y}=\mathbf{f}(\mathbf{x})+\mathbf{G}(\mathbf{x})\mathbf{u}\end{array}$$

(30)

where ${\mathbf{y}}^n={\left[y_1^n,y_2^n,\dots ,y_m^n\right]}^T\in {\Re }^m$, $\mathbf{f}(\mathbf{x})={\left[f_1(x),f_2(x),\dots ,f_m(x)\right]}^T\in {\Re }^m$ and $\mathbf{G}(\mathbf{x})=\left[\begin{array}{ccc}{g}_{11}& \dots & g_{1m}\\ ⋮& \ddots & ⋮\\ g_{m1}& \dots & g_{mm}\end{array}\right]\in {\Re }^{m\times m}$.

The control objective is to derive control law $\mathbf{u}={[u_1,u_2,\dots ,u_m]}^T$ that forces the output $\mathbf{y}(t)={[y_1(t),y_2(t),\dots ,y_m(t)]}^T$ to follow the a desired trajectory ${\mathbf{y}}_d(t)={[y_{1d}(t),y_{2d}(t),\dots ,y_{md}(t)]}^T$. The error matrix $\underline{\mathbf{e}}$ is defined as:

$$\begin{array}{c}\hfill \underline{\mathbf{e}}=\left[\begin{array}{cccc}{y}_1-y_{1d}& y_2-y_{2d}& \cdots & y_m-y_{md}\\ {\dot{y}}_1-{\dot{y}}_{1d}& {\dot{y}}_2-{\dot{y}}_{2d}& \cdots & {\dot{y}}_m-{\dot{y}}_{md}\\ ⋮& ⋮& ⋮& ⋮\\ y_1^{(n-1)}-y_{1d}^{(n-1)}& y_2^{(n-1)}-y_{2d}^{(n-1)}& \cdots & y_m^{(n-1)}-y_{md}^{(n-1)}\end{array}\right]\end{array}$$

(31)

The filter tracking error $\mathrm{{\rm Y}}$ is defined as:

$$\begin{array}{ccc}\hfill {\mathrm{{\rm Y}}}^T& =& {\left[\begin{array}{ccccc}{\theta }_1& {\theta }_2& \cdots & {\theta }_{n-1}& 1\end{array}\right]}^T\underline{\mathbf{e}}=\left[\begin{array}{cc}{\mathrm{\Theta }}^T& 1\end{array}\right]\underline{\mathbf{e}}\hfill \end{array}$$

(32)

where $\mathrm{\Theta }=\left[\begin{array}{cccc}{\theta }_1& {\theta }_2& \cdots & {\theta }_{n-1}\end{array}\right]$ is a vector of FBLC coefficients, and its appropriate choice ensures that $z^{n-1}+{\theta }_{n-1}{z}^{n-1}+\cdots +{\theta }_1$ is asymptotically stable, i.e., $\underline{\mathbf{e}}(k)\to 0$ as $\mathrm{{\rm Y}}(k)\to 0$. The constant value of $\mathrm{\Theta }$ does not provide desired performance over a wide range of operating conditions. In this research work, $\mathrm{\Theta }$ is updated at each time step using the nLMS-based self-tuning scheme. The derivative of Equation (32) results in:

$$\begin{array}{c}\hfill {\dot{\mathrm{{\rm Y}}}}^T={\left[\begin{array}{cccc}{\theta }_1& {\theta }_2& \cdots & \begin{array}{cc}{\theta }_{n-1}& 1\end{array}\end{array}\right]}^T\left[\begin{array}{cccc}{\dot{y}}_1-{\dot{y}}_{1d}& {\dot{y}}_2-{\dot{y}}_{2d}& \cdots & {\dot{y}}_m-{\dot{y}}_{md}\\ {\ddot{y}}_1-{\ddot{y}}_{1d}& {\ddot{y}}_2-{\ddot{y}}_{2d}& \cdots & {\ddot{y}}_m-{\ddot{y}}_{md}\\ ⋮& ⋮& ⋮& ⋮\\ y_1^n-y_{1d}^n& y_2^n-y_{2d}^n& \cdots & y_m^n-y_{md}^n\end{array}\right]\end{array}$$

(33)

$$\begin{array}{c}\hfill {\dot{\mathrm{{\rm Y}}}}^T={\left[\begin{array}{c}{y}_1^n\\ y_2^n\\ ⋮\\ y_m^n\end{array}\right]}^T+{\left[\begin{array}{c}-y_{1d}^n\\ -y_{2d}^n\\ ⋮\\ -y_{md}^n\end{array}\right]}^T+\left[\begin{array}{cc}0& {\mathrm{\Theta }}^T\end{array}\right]\underline{\mathbf{e}}=\left[\begin{array}{c}{y}_1^n\\ y_2^n\\ ⋮\\ y_m^n\end{array}\right]+{\mathbf{y}}_D^T\in {\Re }^{1\times m}\end{array}$$

(34)

Using Equations (30) and (34), we have:

$$\begin{array}{c}\hfill \dot{\mathrm{{\rm Y}}}={\mathbf{y}}^n+{\mathbf{Y}}_D=\mathbf{f}(\mathbf{x})+\mathbf{G}(\mathbf{x})\mathbf{u}+{\mathbf{y}}_D\end{array}$$

(35)

Using $\mathrm{{\rm Y}}=exp(-K_{\mathrm{\Gamma }}t)\Rightarrow \dot{\mathrm{{\rm Y}}}=-K_{\mathrm{\Gamma }}\mathrm{{\rm Y}}$ and $K_{\mathrm{\Gamma }}>0$ ensures that $\underline{\mathbf{e}}\to 0$ as $\mathrm{{\rm Y}}\to 0$. Therefore, the above equation can be written as:

$$\begin{array}{c}\hfill -K_{\mathrm{\Gamma }}\mathrm{{\rm Y}}=\mathbf{f}(\mathbf{x})+\mathbf{G}(\mathbf{x})\mathbf{u}+{\mathbf{y}}_D\end{array}$$

(36)

For the known nonlinear functions $\mathbf{f}(\mathbf{x})$ and $\mathbf{G}(\mathbf{x})$, the following control law is designed according to the feedback linearization method, which would cancel the nonlinearities of the system and would bring $\mathrm{{\rm Y}}$ to zero:

$$\begin{array}{c}\hfill \mathbf{u}={\mathbf{G}}^{-1}(\mathbf{x})[-\mathbf{f}(\mathbf{x})-K_{\mathrm{\Gamma }}\mathrm{{\rm Y}}-{\mathbf{y}}_D]\end{array}$$

(37)

However, it is impractical to implement the control law if the nonlinear functions $\mathbf{f}$ and $\mathbf{G}$ are uncertain or unknown. Therefore, this study proposes to identify the nonlinear functions $\widehat{\mathbf{f}}$ and $\widehat{\mathbf{G}}$ at each time-step through an Adaptive Neuro-Fuzzy Identifier (ANFI) and to online optimize the design parameter $\mathrm{\Theta }$. The control law is now given as:

$$\begin{array}{c}\hfill \mathbf{u}={\widehat{\mathbf{G}}}^{-1}(\mathbf{x})[-\widehat{\mathbf{f}}(\mathbf{x})-\mathrm{\Gamma }]\end{array}$$

(38)

where $\mathrm{\Gamma }=-K_{\mathrm{\Gamma }}\mathrm{{\rm Y}}-{\mathbf{y}}_D$.

Assumption 1.

As the nonlinear MIMO system is controllable, $\mathit{G}(\mathit{x})\in {\Re }^{m\times m}$ is invertible for all $x\in U_c\subset {\Re }^{mn}$.

The following section explains the identification of nonlinear functions and online estimation of the design parameter.

3.2. Adaptive Neuro-Fuzzy Identification

The five-layered ANFI architecture was implemented to identify the nonlinear dynamics of the power system. The inputs to the ANFI are relative speed deviations of different machines w.r.t. the swing machine, and the outputs are the identified nonlinear dynamics of the power system. In the forward pass, the whole network works in a layer fashion, and in the backward pass, error signals are propagated to update the antecedent and consequent parts’ parameters according to the chain rule.

Assumption 2.

The system relative degree n is known, and y, $\dot{y}$, …, $y^n$ are bounded output and available for measurement.

The ANFI architecture of Figure 3 is implemented with $m^2+m$ neuro-fuzzy subsystem blocks that realizes a fuzzy IF-THEN rule in the following form:

• Rule $j_f$: If $x_1$ is $F_{\tau 1}^j$ and $x_2$ is $F_{\tau 2}^j$ and ⋯ and $x_m$ is $F_{\tau m}^j$, Then $\widehat{f_{\tau }}$ is $w_{f_{\tau }}^j$
• Rule $j_g$: If $x_1$ is $G_{\tau h1}^j$ and $x_2$ is $G_{\tau h2}^j$ and ⋯ and $x_m$ is $G_{\tau km}^j$, Then $\widehat{g_{\tau h}}$ is $w_{g_{\tau h}}^j$

where $F_{\tau i}^j$ and $G_{\tau hi}^j$ represent fuzzy sets of the antecedent part for the i-th input, $i=1,2,\cdots ,m$, and the j-th rule, $j=1,2,\cdots ,p$, and $w_{f_{\tau }}^j$ and $w_{g_{\tau h}}^j$, $\tau ,h=1,2,\cdots ,m$, are the consequent weights. The network is implemented with Gaussian membership functions in the antecedent part, a product inference engine and singleton variables in the consequent part.

The function of the nodes in each layer is defined as:

Layer 1: Each node in this layer only transmits input values to the next layer.

Layer 2: The adaptive nodes represent the input fuzzy sets described by the Gaussian membership function. The output of each node is calculated as:

$$\begin{array}{c}\hfill A_i^j=exp\left(-\frac{1}{2}{\left(\frac{x_i-c_i^j}{{\sigma }_i^j}\right)}^2\right)\end{array}$$

(39)

where $c_i^j$ and ${\sigma }_i^j$, receptively, are the mean and standard deviation of the Gaussian membership function for the j-th term of the i-th input variable $x_i$.

Layer 3: Rule nodes constitute the antecedent of the fuzzy rule base. The firing strength of each rule is calculated using the product T-norm. The output of this layer for the j-th rule is defined as:

$$\begin{array}{c}\hfill A^j={\displaystyle \prod _{i=1}^n}{A}_i^j(x_i)\end{array}$$

(40)

Layer 4: Nodes in Layer 4 are called consequent nodes and represent the possible THEN-part of the fuzzy rules. The weights of the adaptive links in this layer represent the singleton $w_{f_{\tau }}^j$ and $w_{g_{\tau h}}^j$. The nodes in this layer compute the weighted firing strengths of each rule from the previous layer.

Layer 5: Nodes in this layer calculate the normalized sum of weighted firing strengths over the sum of firing strengths. The computed defuzzified crisp outputs are identified ${\widehat{f}}_{\tau }(x)$ and ${\widehat{g}}_{\tau h}(x)$ and given as:

$$\begin{array}{c}\hfill {\widehat{f}}_{\tau }=\frac{{\sum }_{j=1}^p{w}_{f_{\tau }}^j{A}^j}{{\sum }_{j=1}^p{A}^j}={\displaystyle \sum _{j=1}^p}{\xi }^j{w}_{f_{\tau }}^j\end{array}$$

(41)

$$\begin{array}{c}\hfill {\widehat{g}}_{\tau h}=\frac{{\sum }_{j=1}^p{w}_{g_{\tau h}}^j{A}^j}{{\sum }_{j=1}^p{A}^j}={\displaystyle \sum _{j=1}^p}{\xi }^j{w}_{g_{\tau h}}^j\end{array}$$

(42)

where ${\xi }^j=\frac{A^j}{{\sum }_{j=1}^p{A}^j}$ is the j-th fuzzy basis function. The MIMO nonlinear system is identified as:

$$\begin{array}{ccc}\hfill \widehat{\mathbf{y}}& =& \widehat{\mathbf{f}}(\mathbf{x})+\widehat{\mathbf{G}}(\mathbf{x})\mathbf{u}\hfill \end{array}$$

(43)

where $\widehat{\mathbf{f}}(\mathbf{x})={\left[\begin{array}{cccc}{\widehat{f}}_1(\mathbf{x})& {\widehat{f}}_2(x)& \cdots & {\widehat{f}}_m(\mathbf{x})\end{array}\right]}^T\mathrm{and}\widehat{\mathbf{G}}(\mathbf{x})=\left[\begin{array}{ccc}{\widehat{g}}_{11}(x)& \cdots & {\widehat{g}}_{1m}(x)\\ ⋮& \ddots & ⋮\\ {\widehat{g}}_{m1}(x)& \cdots & {\widehat{g}}_{mm}(\mathbf{x})\end{array}\right]$.

The identification error, $e_I$, is defined as the error between identified and real dynamics and given as:

$$\begin{array}{c}\hfill e_I(k)=\frac{1}{2}{\epsilon }^T\epsilon =\frac{1}{2}{\left[\widehat{\mathbf{y}}(k)-\mathbf{y}(k)\right]}^2\end{array}$$

(44)

The objective is to accurately identify the dynamics of the power system by minimizing the identification error $e_I(k)$. The conjugate gradient-based parameter learning algorithm is developed to minimize the error through optimization of the parameters of the membership functions and the link weights in the consequent part. The following section explains the CG-based optimization algorithm for the adaptive parameters of ANFI.

Conjugate Gradient Algorithm for Parameter Optimization

CG algorithms comprise a class of optimization algorithms that have strong local and global convergence properties. The nonlinear CG methods are relatively robust in numerical computation and provide sufficient descent independence of line search [42]. Our problem is to minimize the identification error by adapting the parameters of the ANFFBLC.

$$\rho =\{c_{i{f}_{\tau }}^j,{\sigma }_{i{f}_{\tau }}^j,c_{i{g}_{\tau h}}^j,{\sigma }_{i{g}_{\tau h}}^j,w_{f_{\tau }}^j,w_{g_{\tau h}}^j\}$$

where $\rho$ represents the adaptation parameters of the ANFIS. A nonlinear CG algorithm computes a sequence ${\rho }^k,k\ge 1$ with initial guess ${\rho }^0$, using the recurrence as:

$$\begin{array}{c}\hfill {\rho }^{(k+1)}={\rho }^{(k)}+{\chi }_i{\mathbf{d}}^{(k)}\end{array}$$

(45)

where the $\chi >0$ is a design parameter and ${\mathbf{d}}^{(k)}$ is the conjugate search direction, computed as:

$$\begin{array}{c}\hfill {\mathbf{d}}^{(k)}=\left\{\begin{array}{cc}-{\mathbf{g}}^{(k)},\hfill & k=0\\ -{\mathbf{g}}^{(k)}+{\beta }_k{\mathbf{d}}^{(k-1)},\hfill & k\ge 1\end{array}\right.\end{array}$$

(46)

Here, ${\mathbf{g}}^{(k)}=\nabla e_I(\rho )$ is the gradient of the error function at the k-th instant, and ${\beta }_k$ is the CG coefficient, chosen in such a way that the new search direction is orthogonal to the subspace generated by previous search directions. In this research work, the following different CG algorithms are opted to corresponding different choices for the computation of ${\beta }_k$.

• Fletcher and Reeves (FR) presented the first nonlinear CG algorithm with update choice for the CG coefficient, ${\beta }_k^{FR}$, given as [43]:

  $$\begin{array}{c}\hfill {\beta }_k^{FR}=\frac{{∥{\mathbf{g}}^{(k+1)}∥}^2}{{∥{\mathbf{g}}^{(k)}∥}^2}\end{array}$$

  (47)
• Polak-Ribière and Polayk (PRP) proposed a CG method that updated the CG coefficient, ${\beta }_k^{PRP}$, using the following formula [44].

  $$\begin{array}{c}\hfill {\beta }_k^{PRP}=\frac{{\left[{\mathbf{g}}^{(k+1)}\right]}^T\left[{\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)}\right]}{{\left[{\mathbf{g}}^{(k)}\right]}^T{\mathbf{g}}^{(k)}}\end{array}$$

  (48)
• The CG algorithm suggested by Fletcher (CD) has a strong convergence property, and the coefficient, ${\beta }_k^{CD}$, is updated using [45]:

  $$\begin{array}{c}\hfill {\beta }_k^{CD}=\frac{{∥{\mathbf{g}}^{(k+1)}∥}^2}{-{\left[{\mathbf{d}}^{(k)}\right]}^T{\mathbf{g}}^{(k)}}\end{array}$$

  (49)
• The effect of inexact linear search was considered by Liu and Storey (LS) to develop a generalized CG scheme with coefficient ${\beta }_k^{LS}$ updated as [46]:

  $$\begin{array}{c}\hfill {\beta }_k^{LS}=\frac{{\left[{\mathbf{g}}^{(k+1)}\right]}^T\left[{\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)}\right]}{-{\left[{\mathbf{d}}^{(k)}\right]}^T{\mathbf{g}}^{(k)}}\end{array}$$

  (50)
• In [47], a CG algorithm was presented by Dai and Yuan (DY) with a sufficient descent property. The CG coefficient, ${\beta }_k^{DY}$, is updated as:

  $$\begin{array}{c}\hfill {\beta }_k^{DY}=\frac{{∥{\mathbf{g}}^{(k+1)}∥}^2}{{\left[{\mathbf{d}}^{(k)}\right]}^T\left[{\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)}\right]}\end{array}$$

  (51)
• A modified CG algorithm is proposed by Hager and Zhang (HZ) that has a more complicated update formula [48]. The CG coefficient, ${\beta }_k^N$, is updated as:

  $$\begin{array}{c}\hfill {\beta }_k^N={\left[({\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)})-2{\mathbf{d}}^{(k)}\frac{\parallel {\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)}{\parallel }^2}{{\left[{\mathbf{d}}^{(k)}\right]}^T\left[{\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)}\right]}\right]}^T\frac{{\mathbf{g}}^{(k+1)}}{{\left[{\mathbf{d}}^{(k)}\right]}^T\left[{\mathbf{g}}^{(k+1)}-{\mathbf{g}}^{(k)}\right]}\end{array}$$

  (52)

All of these methods do not require the computation of the Hessian matrix that is usually preferred over the methods requiring the Hessian in each iteration. The gradient of identification error is computed in each iteration of the CG algorithm. Using the chain rule of differentiation, the gradient of Equation (44) w.r.t. the respective consequent part parameter is given as:

$$\frac{\partial e_I}{\partial w_{f_{\tau }}^j}=\frac{\partial e_I}{\partial \widehat{y}}\frac{\partial \widehat{y}}{\partial \widehat{f}}\frac{\partial {\widehat{f}}_{\tau }}{\partial w_{f_{\tau }}^j}=\epsilon {\xi }_{f_{\tau }}^j$$

(53)

$$\frac{\partial e_I}{\partial w_{g_{\tau h}}^j}=\frac{\partial e_I}{\partial \widehat{y}}\frac{\partial \widehat{y}}{\partial {\widehat{g}}_{\tau h}}\frac{\partial {\widehat{g}}_{\tau h}}{\partial w_{g_{\tau h}}^j}=\epsilon {\xi }_{g_{\tau h}}^{j}u$$

(54)

where $\frac{\partial e_I}{\partial \widehat{y}}=\epsilon$, $\frac{\partial \widehat{y}}{\partial {\widehat{f}}_{\tau }}=1$, $\frac{\partial \widehat{y}}{\partial {\widehat{g}}_{\tau h}}=u$, $\frac{\partial {\widehat{f}}_{\tau }}{\partial w_{f_{\tau }}^j}={\xi }_{f_{\tau }}^j$ and $\frac{\partial {\widehat{g}}_{\tau h}}{\partial w_{g_{\tau h}}^j}={\xi }_{g_{\tau h}}^j$. For the identification of ${\widehat{f}}_{\tau }$, the gradient of the identification error w.r.t. the parameters of the Gaussian membership function is calculated as:

$$\begin{array}{ccc}\hfill \frac{\partial e_I}{\partial c_{i{f}_{\tau }}^j}& =& \frac{\partial e_I}{\partial \widehat{y}}\frac{\partial \widehat{y}}{\partial \widehat{f_{\tau }}}\frac{\partial {\widehat{f}}_{\tau }}{\partial A_{f_{\tau }}^j}\frac{\partial A_{f_{\tau }}^j}{\partial c_{i{f}_{\tau }}^j}=\epsilon \left(\frac{w_{f_{\tau }}^j-{\widehat{f}}_{\tau }}{{\sum }_j{A}_{f_{\tau }}^j}\right)A_{f_{\tau }}^j\left(\frac{x_i-c_{i{f}_{\tau }}^j}{{({\sigma }_{i{f}_{\tau }}^j)}^2}\right)\hfill \end{array}$$

(55)

$$\begin{array}{ccc}\hfill \frac{\partial e_I}{\partial {\sigma }_{i{f}_{\tau }}^j}& =& \frac{\partial e_I}{\partial \widehat{y}}\frac{\partial \widehat{y}}{\partial {\widehat{f}}_{\tau }}\frac{\partial {\widehat{f}}_{\tau }}{\partial A_{f_{\tau }}^j}\frac{\partial A_{f_{\tau }}^j}{\partial {\sigma }_{i{f}_{\tau }}^j}=\epsilon \left(\frac{w_{f_{\tau }}^j-{\widehat{f}}_{\tau }}{{\sum }_j{A}_{f_{\tau }}^j}\right)A_{f_{\tau }}^j\left(\frac{{(x_i-c_{i{f}_{\tau }}^j)}^2}{{({\sigma }_{i{f}_{\tau }}^j)}^3}\right)\hfill \end{array}$$

(56)

and similarly for the identification of ${\widehat{g}}_{\tau h}$:

$$\begin{array}{ccc}\hfill \frac{\partial e_I}{\partial c_{i{g}_{\tau h}}^j}& =& \frac{\partial e_I}{\partial \widehat{y}}\frac{\partial \widehat{y}}{\partial {\widehat{g}}_{\tau h}}\frac{\partial {\widehat{g}}_{\tau h}}{\partial A_{g_{\tau h}}^j}\frac{\partial A_{g_{\tau h}}^j}{\partial c_{i{g}_{\tau h}}^j}=\epsilon \left(\frac{a_{\tau h}^j-{\widehat{g}}_{\tau h}}{{\sum }_j{A}_{g_{\tau h}}^j}\right)A_{g_{\tau h}}^j\left(\frac{x_i-c_{i{g}_{\tau h}}^j}{{({\sigma }_{i{g}_{\tau h}}^j)}^2}\right)u\hfill \end{array}$$

(57)

$$\begin{array}{ccc}\hfill \frac{\partial e_I}{\partial {\sigma }_{i{g}_{\tau h}}^j}& =& \frac{\partial e_I}{\partial \widehat{y}}\frac{\partial \widehat{y}}{\partial {\widehat{g}}_{\tau h}}\frac{\partial {\widehat{g}}_{\tau h}}{\partial A_{g_{\tau h}}^j}\frac{\partial A_{g_{\tau h}}^j}{\partial {\sigma }_{i{g}_{\tau h}}^j}=\epsilon \left(\frac{a_{\tau h}^j-{\widehat{g}}_{\tau h}}{{\sum }_j{A}_{g_{\tau h}}^j}\right)A_{g_{\tau h}}^j\left(\frac{{(x_i-c_{i{g}_{\tau h}}^j)}^2}{{({\sigma }_{i{g}_{\tau h}}^j)}^3}\right)u\hfill \end{array}$$

(58)

The ANFIS parameters $\rho =\{c_{i{f}_{\tau }}^j,{\sigma }_{i{f}_{\tau }}^j,c_{i{g}_{\tau h}}^j,{\sigma }_{i{g}_{\tau h}}^j,w_{f_{\tau }}^j,w_{g_{\tau h}}^j\}$ are updated through CG algorithms using the respective gradients from Equations (53)–(58) and learning rates ${\chi }_i$ in Equations (45) and (46). The optimized ANFI parameters enable the real-time identification of the power system dynamics to design the control law of Equation (38) with the minimized identification error from Equation (44).

3.3. nLMS-Based Self-Tuning of FBLC Coefficients

The objective of the ANFFBLC is to derive a control law $\mathbf{u}\in {\Re }^m$ such that the system defined in Equation (30) follows the desired trajectory ${\mathbf{y}}_d\in {\Re }^m$ with an acceptable accuracy and bounded states and controls. The control law defined in Equation (38) is based on the identified system model and input $\mathrm{\Gamma }=-K_{\mathrm{\Gamma }}-{\mathbf{y}}_D\in {\Re }^m$. In Equation (32), $\mathrm{\Theta }$ is estimated at each time step using the nLMS algorithm to ensure the appropriate performance of FBLC over a wide range of operating conditions [24].

Assumption 3.

The desired trajectory $y_d(t)$ is continuous, bounded and available for online control computation.

Assumption 4.

There are no zero dynamics.

Equation (32) can be re-written as:

$$\begin{array}{c}\hfill \mathrm{{\rm Y}}=\left[\begin{array}{cc}{\mathrm{\Theta }}^t& 1\end{array}\right]\left[\begin{array}{c}e\\ \dot{e}\\ ⋮\\ e^{(n-1)}\end{array}\right]=\left[\begin{array}{cc}{\mathrm{\Theta }}^t& 1\end{array}\right]\left[\begin{array}{c}\varphi \\ \mathrm{\Phi }\end{array}\right]\end{array}$$

(59)

where ${\varphi }^T=\left[\begin{array}{cccc}e& \dot{e}& \cdots & e^{(n-2)}\end{array}\right]$ and $\mathrm{\Phi }=e^{(n-1)}$. If $\widehat{\mathrm{\Theta }}$ is the estimate of $\mathrm{\Theta }$, the nLMS-based update law can be written as:

$$\begin{array}{c}\hfill {\widehat{\mathrm{\Theta }}}^{(k)}={\widehat{\mathrm{\Theta }}}^{(k-1)}+{\chi }_s\frac{{\phi }^{(k)}}{\phi {}^{(k)T}{\phi }^{(k)}+{\chi }_z}\left\{{\mathrm{{\rm Y}}}^{(k)}-{\phi }^{(k)}{\widehat{\mathrm{\Theta }}}^{(k-1)}-{\varrho }^{(k)}\right\}\end{array}$$

(60)

where $0<{\chi }_s<2$ is the step size and ${\chi }_z>0$ is a design parameter introduced to avoid the division by zero. The online tuning of FBLC coefficients ensures the appropriate choice of $\mathrm{\Theta }$ that results in the desired damping performance of ANFFBLC over a wide range of operating conditions.

3.4. Computational Steps for Closed-Loop Control System

The following steps summarize the realization of the algorithm for damping LFOs by ANFFBLC-based multiple HVDC links in the AC/DC power system.

• The wide area measurement system transmits the actual speed of generators, ${\omega }_1,{\omega }_2,\cdots ,{\omega }_7$, to HVDC control.
• Speed deviations of the generators w.r.t. the swing machine are computed as plant output, ${\mathbf{y}}^T=\left[\begin{array}{cc}{y}_1& y_2\end{array}\right]$, where $y_1={\sum }_{i=2}^5({\omega }_i-{\omega }_1)$ and $y_2={\sum }_{i=6}^7({\omega }_i-{\omega }_1)$.
• ANFFBLC captures the nonlinear dynamics, $\widehat{\mathbf{f}}$ and $\widehat{\mathbf{G}}$, of the power system using $\mathbf{y}$ as explained in Section 3.2. The parameters of ANFIS of ANFFBLC are instantaneously optimized through the CG algorithm to minimize the identification error defined in Equation (44).
• At the same instant, the online estimation generates the appropriate input $\mathrm{\Gamma }$. The coefficients, $\mathrm{\Theta }$, are optimized through the nLMS algorithm to minimize the tracking error of Equation (31).
• ANFFBLC generates the appropriate control law, ${\mathbf{u}}^T=\left[\begin{array}{cc}{u}_1& u_2\end{array}\right]$ that is based on identified functions $\widehat{\mathbf{f}}$ and $\widehat{\mathbf{G}}$ and optimized $\mathrm{\Gamma }$ as given in Equation (38).
• The ANFFBLC output constitutes $I_{damp1}$ and $I_{damp2}$ for master controls of HVDC Link 1 and Link 2, respectively. The $I_{damp}$ modulates the current order for HVDC link using Equation (23).
• At each pole control, the current order $I_{ord}$ generates the corresponding ignition delay angle, $\alpha$, that controls the power flow through the HVDC system
• During normal operation, $I_{damp}=0$, and power flow through the HVDC system is set to a pre-specified value. During perturbed operating conditions, the generators oscillate against each other at a speed different from the set value. The speed deviations are detected by ANFFBLC, and the appropriate damping signal is generated for each HVDC system. The precise power flow control through the HVDC link improves the damping of power oscillations present in the system and improves its stability.

4. Simulation Results and Discussion

The test power system comprises the simplified model of the highly interconnected Hydro-Québec (HQ) power system with two HVDC links implemented in SIMULINK/SimPowerSystem [49]. Seven power plants and loads in different areas are connected through a huge 735-kV transmission network. Major load buses are MTLand QUE with 15,500-MW and 6000-MW loads, respectively. The load includes dynamic load, constant power PQ , constant impedance (Z) and induction motor load. Two HVDC links, each 1000 MW at 500 kV, are installed in the transmission network connecting buses LG27 to MTL and CHU to QUE for bulk power transmission. Sequential AC/DC load flow is used to initialize the system with machine ${\mathrm{G}}_1$ taken as a swing bus.

The test power system is subjected to various contingencies with different operating conditions, locations and duration. The system response is evaluated through the observation of the following parameters.

• Low-frequency oscillations: Observation of over-shoot and settling time for modes with oscillations of generators ${\mathrm{G}}_2$ to ${\mathrm{G}}_6$ against the reference generator ${\mathrm{G}}_1$ i.e., ${\omega }_i-{\omega }_1,i=2,3,\cdots ,7$.
• Performance indexes: The two used indexes for quantitative measures are Integral of Time-weighted Squared Error (ITSE) and Integral of Time-weighted Absolute Error (ITAE) [50]. For performance indexes, the error is calculated as $e_p={\sum }_{i=2}^7({\omega }_1-{\omega }_i)$.

In each scenario, for the comparison, to verify the effectiveness of damping performance, simulations were carried out with conventional PID control, Adaptive PID (AdapPID) [51], Direct Intelligent Neuro-Fuzzy Control (DirINF) [28], SD-based ANFFBLC and six CG methods.

4.1. Scenario # 1

At t = 1 s, a three-phase to ground fault occurs at ${\mathrm{HVDC}}_1$ inverter bus MTL. The fault is cleared after 100 ms with a 1000-MW reduction in load at the MTL bus. Targeted oscillation modes are shown in Figure 4. Well-tuned conventional and non-conventional control schemes have competitive damping performance as compared to ANFFBLC with SD and different CG optimization techniques. During the first cycle, the overshoot is reduced by 6–26% with ANFFBLC damping schemes as compared to benchmark controls. For all modes, LFOs are damped completely with all machines settling to steady-state speed. The damping performance of HVDC links with damping control is also assessed through performance indexes. Two performance indexes ITSE and ITAE are shown in Figure 5. The ITSE plot represents time-weighted speed deviation error during the transient-state. The best improvement of 55% in the minimization of the transient-state error is achieved by ANFFBLC-HZ in comparison with PID control, while the ITSE plot of ANFFBLC-FR shows 50% improvement over PID. The ITAE plot represents time-weighted errors persisting for a long time and depicts the damping efforts during the steady-state conditions. The ITAE index for ANFFBLC-HZ shows 36% improvement in damping persistent oscillations as compared to the conventional control. ANFFBLC with SD and CG optimization depicts damping performance in close competition; however, the ANFFBLC-HZ index shows the most effective damping of LFOs. ANFFBLC shows 31–36%, 18–24% and 13–19% improvement in steady-state stability as compared to PID, AdapPID and DirINF, respectively.

4.2. Scenario # 2

The CHM-JCA transmission line is exposed to a three-phase to ground fault at time $t=1$ s. The fault is cleared after 120 ms with the permanent outage of the transmission line. Oscillation modes, shown in Figure 6, reveal the more severe nature of the disturbance as post-fault operating conditions have changed. The plots show excessive overshoots with benchmark control schemes, while the ANFFBLC has rapidly dampedthe oscillations, and system settleto steady-state. In some cases, the competitive performance of conventional control is observed during steady-state; however, in the transient-state, 8–16% reduced overshoots are observed with ANFFBLC.

The ITSE and ITAE indexes shown in Figure 7 represent the damping performance comparison. HVDC links with MIMO ANFFBLC have lesser ITSE index values at any time and more flat profiles as compared to the ITSE plot of the benchmark control systems. As compared to conventional PID control, the ITSE plot for ANFFBLC-HZ depicts 63% improvement in damping the transient-state oscillations. The plots for the ANFFBLC schemes get flattened prior to the PID plot that shows the rapid minimization of error by ANFFBLC in the transient-state. Time-weighted index ITAE for ANFFBLC-HZ shows 49% improvement in the steady-state stability of the power system as compared to PID control. ITSE improvement of 46–58%, 38–52% and 30–46% is observed with ANFFBLC as compared to the PID, AdapPID and DirINF control schemes, respectively. While in the steady-state, improvement with ANFFBLC is 28–37%, 20–30% and 10–21%, as compared to the PID, AdapPID and DirINF controllers schemes, respectively. The ANFFBLC-HZ scheme has well contained the LFOs in the transient-state and steady-state conditions and maintains the best performance among different CG methods.

4.3. Scenario # 3

The test power system is subjected to a sequence of faults to investigate the effectiveness and robust nature of the proposed control. A 100-ms three-phase to ground fault occurs at $t=1$ s at bus QUE. Another fault arises at bus SAG at time $t=5$ s for a duration of 110 ms. Both faults are self-cleared without any structural change in the power system. The response of the system is observed through LFOs modes as shown in Figure 8. In terms of overshoot, the performances of ANFFBLC and conventional control are competitive after the first fault. However, after the second fault, the performance of the ANFFBLC is much more improved than benchmark controllers in terms of overshoot and settling time. The oscillatory behavior of the power system with PID control is observed for different oscillation modes. A better picture of damping performance is described by the performance indexes shown in Figure 9. As perceived by the ITSE plots, transient-state errors are minimized with ANFFBLC and benchmark controllers for both faults. However, ITSE indexes for ANFFBLC schemes have a lesser value and rapidly maintain constant index values during the steady-state condition after the second fault. The ITSE index for ANFFBLC schemes shows overall 39–61%, 34–57% and 15–45% improvement in the minimization of the speed deviation error during transient-state conditions as compared to PID, AdapPID and DirINF controls, respectively. Similar performance trends are observed in ITAE index with 23–35%, 21–34% and 12–26% improvement with ANFFBLC schemes. In Scenario 3, ANFFBLC-HZ has repeated the best damping performance like previous scenarios and shows the most improved transient-state and steady-state conditions.

4.4. Performance Comparison of CG Algorithms

In an application of different CG algorithms in the power system, all methods performed very close to each other and better than conventional control and non-conventional controllers. The damping performance of ANFFBLC with different CG schemes and benchmark controllers is summarized in

Table 1. Performance improvement of ANFFBLC as compared to conventional control.

Scenario	Performance	%Age Improvement AdapPID, DirINF and ANFFBLC w.r.t. PID
	Index	AdapPID	DirINF	SD	FR	CD	DY	PRP	LS	HZ
1	ITAE	16	21	31	33	33	35	35	35	36
	ITSE	24	33	47	50	51	51	53	54	55
2	ITAE	9	20	28	28	30	31	32	33	37
	ITSE	12	23	46	47	50	52	53	54	58
3	ITAE	2	13	23	26	31	32	32	33	35
	ITSE	8	29	39	44	54	54	56	56	61

. It shows the performance improvement attained by ANFFBLC as compared to conventional and non-conventional control in transient-state and steady-state conditions. The best performance is observed with the ANFFBLC-HZ scheme, while the FR method only performed better than the benchmark controls. The FR scheme suffers a jamming problem in some cases [52], and a modification was proposed in the Polak-Ribière scheme that enables a restart feature to avoid jamming and more rapid convergence as compared to the FR method. The proposed CD method with relaxed constraints is closely related to the FR scheme, and the search directions satisfy the sufficient descent condition under strong Wolfe line search. However, in some particular cases, the CD method may converge to a point where $g_k\ne 0$ [53]. The DY method always generates descent directions and converges globally with the line search satisfying the standard Wolfe condition. The DY method is also susceptible to jamming. The LS method exhibits similar performance as the PRP scheme with exact line search and has an in-built restart feature to avoid the jamming problem. Computation performance of the PRP and LS methods is better than the FR, PRP and DY methods. The HZ method modified the linear CG method of Hestenes and Stiefel that provides sufficient descent with a relaxed accuracy of the line search. In [48], the global convergence is proven to show that for any function and any line search, the sufficient descent condition is satisfied, and the jamming is avoided. That made the ANFFBLC-HZ scheme more robust in computation and provided better optimization of ANFIS parameters to ensure minimized identification error and the best damping performance over the other classical CG algorithms [54].

5. Conclusions

The article presented an MIMO POD controller for multiple HVDC transmission systems. The ANFFBLC schemes effectively modulate the real power flow through the HVDC system to enhance its damping assistance during perturbed conditions. Excited power oscillations under any disturbance are apprehended by the ANFFBLC on the basis of the measured speed signal of generators. SD and six CG algorithms are employed to optimize the parameters of ANFIS that minimize the identification error and captured the updated plant dynamics without a priori knowledge of the system model. The damping performance is investigated through the multi-machine AC/DC power system exposed to the disturbances of different severity. Results obtained for a wide range of operating conditions indicate the improved damping performance of different proposed CG-based ANFFBLC schemes as compared to ANFFBLC-SD, conventional PID control, AdapPID and DirINF. Among CG algorithms, the HZ method shows the best optimization capabilities with strong convergence to optimal minima of the optimization function. Rapid identification of the plant model with the HZ method enables ANFFBLC-HZ to derive the desired control output that effectively damped LFO in the power system. The future work includes the implementation of MIMO nonlinear POD control for multiple HVDC links and FACTS controllers, as well as the investigation of the effect on damping the low-frequency oscillations.