Computationally Efficient Optimal Control for Unstable Power System Models

Uddin, Mahtab; Uddin, M. Monir; Hakim Khan, Md. Abdul

doi:https://doi.org/10.1155/2021/8071869

Mathematical Problems in Engineering

On this page

Abstract Introduction Preliminaries Conclusion Data Availability Disclosure Conflicts of Interest Acknowledgments References Copyright Related Articles

Research Article | Open Access

Volume 2021 | Article ID 8071869 | https://doi.org/10.1155/2021/8071869

Computationally Efficient Optimal Control for Unstable Power System Models

Mahtab Uddin,^1,2M. Monir Uddin,³and Md. Abdul Hakim Khan²

Academic Editor: Filippo Cacace

Received03 Jun 2021

Revised14 Aug 2021

Accepted08 Sept 2021

Published19 Oct 2021

Abstract

In this article, the focus is mainly on gaining the optimal control for the unstable power system models and stabilizing them through the Riccati-based feedback stabilization process with sparsity-preserving techniques. We are to find the solution of the Continuous-time Algebraic Riccati Equations (CAREs) governed from the unstable power system models derived from the Brazilian Inter-Connected Power System (BIPS) models, which are large-scale sparse index-1 descriptor systems. We propose the projection-based Rational Krylov Subspace Method (RKSM) for the iterative computation of the solution of the CAREs. The novelties of RKSM are sparsity-preserving computations and the implementation of time-convenient adaptive shift parameters. We modify the Low-Rank Cholesky-Factor integrated Alternating Direction Implicit (LRCF-ADI) technique-based nested iterative Kleinman–Newton (KN) method to a sparse form and adjust this to solve the desired CAREs. We compare the results achieved by the Kleinman–Newton method with that of using the RKSM. The applicability and adaptability of the proposed techniques are justified numerically with MATLAB simulations. Transient behaviors of the target models are investigated for comparative analysis through the tabular and graphical approaches.

1. Introduction

For the practical purposes, optimal control is a vital part of the engineering interest, e.g., industrial control systems, system defense strategy, voltage stability of the power systems, and signal processing [1, 2]. Multitasking systems having various components arise in many fields of engineering applications, such as microelectronics, microelectromechanical systems, cybersecurity, computer control of industrial processes, communication systems, and the like. These systems are composed of branches of subsystems and are functioned by large mathematical models utilizing the interrelated inner mathematical system of higher dimensions.

Power system models are one of the prime branches of the application of optimal controls. For the multiconnected power systems, avoiding cyber-attacks, ensuring the stability of the power connections, and maintaining the compatible frequency level are essential. For those inevitable issues, optimal controls have exigent roles [3].

The dynamic of a large-scale power system model can be described by the Differential Algebraic Equations (DAEs) aswhere, is the vector with differential variables, is the vector with algebraic variables, and is the vector of parameters with [4, 5]. In the state-space representation, the dynamic state variables and instantaneous variables are defined for the specific system, where the parameter defines the configuration and the operation condition.

The state variables in are time-dependent generator voltages, and the parameter is composed of the system parameters. The control devices together form , and the power flow balance forms . In case of voltage stability, some equations of the will not be considered. For a fixed parameter , linearizing system (1) around the equilibrium point will provide the following Linear Time-Invariant (LTI) continuous-time system with the input-output equations in the sparse form with the block matrices aswhere , , , and with very large and represent differential coefficient matrix, state matrix, control multiplier matrix, state multiplier matrix, and direct transmission map, respectively [6]. In system (2), is the state vector and is control (input), while is the output vector and considering as the initial state. In most of the state-space representations, the direct transmission remains absent and because of that . Because is singular (i.e., ), system (2) is called the descriptor system [7, 8].

Here, , with are state vectors, and other submatrices are sparse in appropriate dimensions. If and have full rank, and is nonsingular (i.e., ), the system is called the index-1 descriptor system [9]. In the current work, we will focus on the stabilization of index-1 descriptor system only. By proper substitution and elimination, descriptor system (2) can be converted to the generalized LTI continuous-time systemwhere we have considered the following relations:

Lemma 1 (equivalence of transfer functions [10]). Assume for , the transfer functions and are obtained from index-1 descriptor system (2) and converted generalized system (3), respectively. Then, the transfer functions and are identical, and hence, those systems are equivalent.

Though systems (2) and (3) are equivalent, system (2) is sparse and system (3) is dense. This explicit conversion is contradictory with the aim of the work, and it will be bypassed by a more efficient way.

LTI continuous-time systems are the pivot ingredient of the present control theory and many other areas of science and engineering [11]. Continuous-time Algebraic Riccati Equation (CARE) appears in many branches of engineering applications; especially in electrical fields [12, 13]. The CARE connected to system (3) is defined as

If the Hamiltonian matrix corresponding to system (3) has no pure imaginary eigenvalues, then the solution of the CARE (5) exists and unique [14]. The solution of (5) is symmetric positive definite and called stabilizing for the stable closed-loop matrix . Riccati-based feedback matrix has a prime role in the stabilization approaches for unstable systems [15, 16]. To find an optimal feedback matrix , the Linear Quadratic Regulator (LQR) problem technique can be applied, where the cost functional is defined as

Cost functional (6) can be optimized as by applying the optimal control generated by the optimal feedback matrix associated with the solution matrix of the CARE (6). Using the optimal feedback matrix , an unstable LTI continuous-time system can be optimally stabilized by replacing with . The stabilized system can be written as

The eligibility of Rational Krylov Subspace Method (RKSM) for the large-scale LTI continuous-time systems has discussed by Simoncini [17], and the application of adaptive RKSM to solve large-scale CARE for finding optimal control of the LTI systems has narrated by Druskin and Simoncini [18]. Analysis of the basic properties of RKSM for solving large-scale CAREs subject to LTI systems was investigated by Simoncini [19], where the author briefed a new concept of shift parameters that is efficient for the perturbed systems. Very detailed discussion on the numerical solution of large-scale CAREs and LQR-based optimal control problems are given by Benner et al. [20], where the authors narrated the extensions Newton method by means of Alternating Direction Implicit (ADI) technique with convergence properties and the comparative analysis with the RKSM approach. For solving large-scale matrix equations, the Kleinman–Newton method based on Low-Rank Cholesky-Factor ADI (LRCF-ADI) technique is discussed by Kürschner [21].

Because there is a scarcity of efficient computational solvers or feasible simulation tools for large-scale CAREs governed from the unstable power system models, the concentration of this work is to develop sparsity-preserving efficient techniques to find the solution of those CAREs. We are proposing a sparsity-preserving and rapid convergent form of the RKSM algorithm for finding solutions of CAREs associated with the unstable power system models and implementation of Riccati-based feedback stabilization. Also, a modified sparse form of LRCF-ADI-based Kleinman–Newton method for solving CAREs subject to unstable power system models and corresponding stabilization approach is proposed. The proposed techniques are applied for the stabilization of the transient behaviors of unstable Brazilian Inter-Connected Power System (BIPS) models [22]. Moreover, the comparison of the computational results will be provided in both tabular and graphical methods.

2. Preliminaries

In this section, we discuss the background of the proposed techniques, which are derived for generalized LTI continuous-time systems. It includes the derivation of the Rational Krylov Subspace Method (RKSM) technique for solving the Riccati equation, basic structure of the Kleinman–Newton method, and derivation of Low-Rank Cholesky-Factor Alternative Direction Implicit (LRCF-ADI) method for solving Lyapunov equations.

2.1. Generalized RKSM Technique for Solving Riccati Equations

In a study [19], Simoncini applied RKSM approach for solving the CARE in the formassociated with the LTI continuous-time system

If the eigenvalues of the matrix pair satisfy that ensures that the solution of the CARE (8) exists and unique. The orthogonal projector spanned by the -dimensional rational Krylov subspace for a set of shift parameters is defined as

If are the eigenvalues of approximates the mirror eigenspace of and is its border, the shifts are computed from

According to the Galerkin condition and after simplification by matrix algebra, a low-rank CARE can be obtained aswhere , , , , and . Equation (12) is an approximated low-rank CARE and can be solved by any conventional method or MATLAB care command. Here, is taken as low-rank approximation of , corresponding to the low-rank CARE (12). Then, residual of the -th iteration iswhere denotes the Frobenius norm and is a block upper triangular matrix in the QR factorization of the matrix defined aswhere is a block upper Hessenberg matrix, and is the matrix formed by the last columns of the -order identity matrix. For such that , the relative residual can be estimated as

Through RKSM, the low-rank factor of the approximate solution of the CARE (8) needs to be estimated, such that . The low-rank solution is symmetric positive definite, and the original solution can be approximated as . By the eigenvalue decomposition to the approximate solution and truncating the negligible eigenvalues, the possible lowest order factor of can be estimated as

Here, consists the negligible eigenvalues. Summary of the above process is given in Algorithm 1.

	Input:, (number of iterations) and (initial shifts).
	Output: Low-rank factored solution such that .
(1)	Compute (QR factorization).
(2)	Choose .
(3)	while not converged or do
(4)	Solve .
(5)	Compute shift for the next iteration.
(6)	Using Arnoldi algorithm, orthogonalize against to obtain , such that .
(7)	Assuming , , , and , for , solve the reduced-order Riccati equation (12).
(8)	Compute for convergence.
(9)	end while
(10)	Compute eigenvalue decomposition .
(11)	For negligible eigenvalues, truncate and compute .

2.2. Modified LRCF-ADI Method for Solving Lyapunov Equations

The generalized Continuous-time Algebraic Lyapunov Equation (CALE) associated with LTI continuous-time system (9) is

The shift parameters are allowed, and the initial iteration is taken as . Assume as the low-rank Cholesky factor of such that [23]. The ADI algorithm can be formed in terms of Cholesky factor of , and there will be no need to estimate or store at each iteration as only is required [24].

Considering , the conventional low-rank Cholesky factor ADI iterations yield the form as

Assume a set of adjustable shift parameters, for two subsequent block iterates and of the ADI technique related to the pair of complex conjugated shifts , it holdswhere indicates the complex conjugate with and iterates associated to real shifts are always purely real [25].

Then, the following matrix for the basis extension can be obtained:

Then, for a pair of complex conjugate shifts at any iteration, the low-rank factor can be computed as

The modified techniques discussed above is summarized in Algorithm 2.

	Input:, (number of iterations) and shift parameters .
	Output: Low-rank Cholesky-factor such that .
(1)	Consider .
(2)	for to do
(3)	ifthen
(4)	Solve .
(5)	else
(6)	Compute .
(7)	end if
(8)	ifthen
(9)	Update .
(10)	else
(11)	Assume and .
(12)	Update .
(13)	Compute .
(14)
(15)	end if
(16)	end for

2.3. Basic Kleinman–Newton Method

For system (9), the residual of the CARE (8) is

From the Fréchet derivative at , we have

Consider the Newton iteration and apply (23), we have

Now, put in (24). Then, after simplification, we get

Then, assume and , then equation (25) reduces to a generalized Lyapunov equation such as

Generalized CALE (26) can be solved for by any conventional method, such as Low-Rank Cholesky-Factor Alternative Direction Implicit (LRCF-ADI) method and the corresponding feedback matrix can be estimated. The whole mechanism is called the Kleinman–Newton method for solving generalized CARE [26]. The summary of the method is given in Algorithm 3.

	Input:, (number of iterations) and (initial assumption).
	Output: Approximate solution and feedback matrix .
(1)	whiledo
(2)	Compute and .
(3)	For solve .
(4)	Compute .
(5)	end while

3. Solving Riccati Equations Arising from the Index-1 Descriptor Systems

In this section, we discuss the updated RKSM techniques for solving the Riccati equation derived from index-1 descriptor systems, it includes the stopping criteria, sparsity preservation, and estimation of the optimal feedback matrix. Then, the LRCF-ADI-based Kleinman–Newton method with the adjustment for index-1 descriptor systems and finding the optimal feedback matrix are discussed. Finally, the stabilized system is narrated accordingly.

3.1. Updated Rational Krylov Subspace Method

Let us consider LTI continuous-time system (3) and introduce an orthogonal projector spanned by the -dimensional rational Krylov subspace for a set of shift parameters that is defined as

To find the shift parameters , the eigenvalues of the matrix pencil are evaluated and sorted numerically. The required number of dominant eigenvalues are to store for the upcoming iterations. If some of the dominant eigenvalues have the negative real part, they are to shift to the positive complex plane by the mirror imaging process.

Again, consider CARE (5) and apply the Galerkin condition on it. Then, after the simplification by matrix algebra, a low-rank CARE can be achieved aswhere , , , , and . Equation (28) is a low-rank CARE and can be solved by MATLAB care command or any existing methods, such as Schur-decomposition method.

For the quick and smooth convergence of the proposed algorithm, adjustable shift selection is crucial, and we are adopting the adaptive shift approach for index-1 descriptor systems [27]. This process required to be recursive, and in each step, the subspace to all the projector generated with the current set of shifts will be extended.

3.1.1. Stopping Criteria and Related Theorem

Arnoldi relation is a very essential tool for the computation of residual of the RKSM iterations. To avoid extra matrix-vector multiplies per iteration, the computation of projected matrix can be performed more efficiently than the explicit product . To find the stopping criteria, we have to consider the following lemmas.

Lemma 2 (Arnoldi relation [28]). Let be the rational Krylov subspace with the shift parameters . Then, satisfies the Arnoldi relation as follows:where is the QR decomposition of the right-hand side matrix with .

Lemma 3 (building the projected matrix [29]). Let the column vectors of be an orthonormal basis of the rational Krylov subspace with the block diagonal matrix , where is the set of shift parameters used in the algorithm. Then, for the projected matrix , the following relation holds:

Theorem 1 (residual of the RKSM iterations). Let be the orthogonal projector spanned by the rational Krylov subspace and is the solution of the CARE (5) using the low-rank solution . Then, the residual of -th iteration can be computed aswhere denotes the Frobenius norm, and is a block upper triangular matrix in the QR factorization of the matrix defined as

Proof. Assume and consider the reduced QR factorization . Then, by putting the relations in equation (29) of Lemma 2, the relation can be written asThe residual of CARE (5) can be written asConsider the approximate solution using the low-rank solution as and equation (30) in Lemma 3, then applying (33) in (34), we getThus, the proof follows from the Frobenius norm of equation (35).
Then, the relative residual with respect to can be estimated as follows:

3.1.2. Sparsity Preservation

The matrix in (3) is in dense form, which is contradictory to the aim of the work and the rate of convergence of the converted system is slow enough. So, to bypass these drawbacks at each iteration, a shifted linear system needs to be solved for as

Here, is the truncated term. Linear system (37) is higher dimensional but sparse and can be solved by the conventional sparse-direct solvers very efficiently [30]. To improve the consistency of the RKSM approach, explicit form of the reduced-order matrices must not be used to construct reduced-order system. The sparsity-preserving reduced-order matrices can be attained by the following way:

3.1.3. Treatment for the Unstable Systems

If the system is unstable, a Bernoulli stabilization is required through an initial-feedback matrix to estimate , and the matrix needs to be replaced [31]. Then, system (3) and the CARE (5) need to be redefined as

Then, for every iteration, needs to be updated by the solution , which is the low-rank approximant of the solution of (40). For this, the rational Krylov subspace for the projector needs to be redefined as

For the stabilized system using the initial-feedback matrix , expressions (37) can be written as

To evaluate the shifted linear systems, explicit inversion of should be avoided in practice, instead the Sherman–Morrison–Woodbury formula needs to be used as follows:

3.1.4. Estimation of the Optimal Feedback Matrix

The low-rank solution is symmetric and positive definite and can be factorized as . The original solution can be approximated as . Finally, the desired low-rank factored solution of the CARE (5) will be stored, and the optimal feedback matrix can be estimated. This process is iterative and will continue until the desired convergence is achieved. The whole process is summarized in Algorithm 4.

	Input: (initial feedback matrix), (number of iterations), and (initial shifts).
	Output: Low-rank factored solution such that optimal feedback matrix .
(1)	Compute (QR factorization).
(2)	Choose .
(3)	Choose .
(4)	while not converged or do
(5)	Solve the linear system (42) for .
(6)	Compute adaptive shift parameters for the next iterations (if store is empty).
(7)	Using Arnoldi algorithm, orthogonalize against to obtain , such that .
(8)	Assuming , , , and , for solve the reduced-order Riccati equation (28).
(9)	Update .
(10)	Compute for convergence.
(11)	end while
(12)	Compute eigenvalue decomposition .
(13)	For negligible eigenvalues, truncate and compute .
(14)	Compute the optimal feedback matrix .

3.2. Updated LRCF-ADI-Based Modified Kleinman–Newton Method

In each iteration of Algorithm 3, the generalized CALE needs to be solved for once, and there are several techniques available to do it. In his Ph.D. thesis, Kuerschner discussed low-rank ADI approaches (Algorithm 3.2 in Chapter 3 and Algorithm 6.2 in Chapter 6) for solving generalized CALE derived from generalized CARE in the iterative loops of Kleinman–Newton algorithm [21]. Now, we need to implement above mechanisms for the index-1 descriptor system in the sparse form. For the adjustment, some modifications are required as given below.

Here, the initial stage of the computation of the shift parameters is the same as narrated in previous section. But due to the structure of the subspace for the ADI technique, the dominant eigenvalues with the negative real part are to store, and those having the positive real part are to shift to the negative complex plane by the mirror imaging process.

3.2.1. Convergence Criteria and Recurrence Relations

The computation of residuals of the ADI iterations can be achieved using the simplified technique implementing the adjustable shift parameters. This approach will be efficient for time dealing and memory allocation. Following lemmas represent some important properties of the ADI iterates.

Lemma 4 (convergence of ADI iterates [32]). Let be the solution at the -th iteration of the CALE as follows:

Consider and are the successive iterates of the ADI method. Then, for all , the relation holdswhere .

Lemma 5 (residual of the ADI iterates [32]). Let an iterate of the CALE (44) by the ADI method. Then, considering for all , the residuals at have the form

Theorem 2 (residual factor of the ADI iterations). The residual for CALE (44) at -th iteration of ADI method is of rank at most and is given bywith , for all .

Proof. Consider is the solution of CALE (44) and is its -th iterate by ADI method. Then, the residual of the -th iteration in terms of can be written asFor all , the Stein’s equation is equivalent to CALE (44), and it can be written asUsing the Stein’s equation (49) in Lemma 4, we can find some initial iterates for which the residuals . Then, according to Lemma 5, residual (48) can be defined aswith .
Thus, for , the ADI iterations in the LRCF-ADI algorithm needs to be stopped, where is a given margin of tolerance. Then, with the residual factor relation for can be derived asThen, using (51) the residual factor can be derived in a recursive form asIn case of real setting, needs to be considered to find the following form:The summary of above techniques is given in Algorithm 5.

	Input:, (tolerance), (number of iterations), and shift parameters .
	Output: Low-rank Cholesky-factor such that .
(1)	Consider , and .
(2)	while or do
(3)	Solve.
(4)	ifthen
(5)	Update .
(6)	Compute
(7)	else
(8)	Assume and .
(9)	Update .
(10)	Compute .
(11)
(12)	end if
(13)
(14)	end for

3.2.2. Adjustment for the Unstable Sparse Systems

For the unstable index-1 descriptor system, initial feedback matrix needs to be introduced, and instead of , corresponding shift parameters are needed to be computed from eigenpair . The sparse form of the eigenpair can be structured as

To find , in each ADI (inner) iteration, a shifted linear system needs to be solved as

Thus, can be obtained from the sparse form of the shifted linear system structured as

3.2.3. Estimation of the Optimal Feedback Matrix

The feedback matrix needs to be computed in each ADI (inner) iteration, and the optimal feedback matrix needs to be stored after the final Newton (outer) iteration. The summary of the modified method is given in Algorithm 6.

Input: (initial feedback matrix), (tolerance), and (number of iterations).
Output: Low-rank Cholesky-factor such that and optimal feedback matrix .
(1)	for to do
(2)	Choose , and .
(3)	Assume .
(4)	Compute adaptive shift parameters from the eigenpair defined in (54).
(5)	whiledo
(6)
(7)	Solve the linear system (56) for .
(8)	ifthen
(9)	Update .
(10)	Compute .
(11)	Compute .
(12)	else
(13)	Assume , , and .
(14)	Compute .
(15)	Update .
(16)	Compute .
(17)	Compute .
(18)
(19)	end if
(20)	end while
(21)	Update and
(22)	end for

3.3. Computation of the Optimally Stabilized System and the Optimal Control

The optimal feedback matrix can be achieved by the feasible solution of Riccati equation (5). Then, applying , optimally stabilized LTI continuous-time system can be written as (7). To preserve the structure of the system, it needs to back to the original structure (2), and for this, the submatrices , are replaced by , , respectively.

Finally, the desired optimal control of the targeted power system models can be computed.

4. Numerical Results

The stability of the target models is investigated, and the unstable models are stabilized through the Riccati-based feedback stabilization process. The proposed methods are employed to find the solution of the Riccati equation arising from the BIPS models, and corresponding feedback matrices are generated for system stabilization. Also, initial Bernoulli feedback stabilization is implemented for the convenient rate of convergence. All the results have been achieved using the MATLAB 8.5.0 (R2015a) on a Windows machine having Intel-Xeon Silver 4114 CPU 2.20 GHz clock speed, 2 cores each, and 64 GB of total RAM.

4.1. Brazilian Interconnected Power System Models

Power system models are an essential part of engineering fields that consists of simulations based on power generations and grid networks. The computation required to analyze electrical power systems employing mathematical models utilizing real-time data. There are several applications of the power system model, i.e., electric power generation, utility transmission and distribution, railway power systems, and industrial power generation [33].

The Brazilian Interconnected Power Systems (BIPS) is one of the most convenient examples of the power system models with various test systems [34]. Table 1 provides the details about the models. The detailed structure of those will be found at the webpage of Joost Rommes, where all of them are index-1 descriptor system. The models mod-606, mod-1998, mod-2476, and mod-3078 have the unstable eigenvalues, whereas the models mod-1142, mod-1450, and mod-1693 have stable eigenvalues [35]. Here, the names of the models are considered according to their number of states.

4.2. Comparison of the Results Found by RKSM and KN-LRCF-ADI

The CAREs arising from the models mod-606, mod-1998, and mod-2476 are efficiently solved and stabilized the corresponding models by both RKSM and KN-LRCF-ADI techniques. As model mod-3078 is semistable, the computation of CARE derived from this model is not possible by LRCF-ADI techniques, but by the RKSM approach, model mod-3078 successfully stabilized, and the numerical result for model mod-3078 is investigated for RKSM only.

Table 2 depicts the numerical results of the stabilization process via RKSM for the unstable BIPS models and various properties of the stabilized systems are illustrated, whereas Table 3 displays the several modes of ADI techniques in KN-LRCF-ADI method for stabilizing the unstable BIPS models, including characteristics of the stabilized models. In Tables 2 and 3, we have analyzed the same features of the stabilized BIPS models, we can easily compare the efficiency and robustness of the proposed methods.

From Tables 2 and 3, it can be said that the proposed RKSM approach has quick convergence ability and occupies very small solution space to provide the efficient solution of the CAREs. In contrast, LRCF-ADI-based Kleinman–Newton has several approaches for finding the solution of CAREs, whereas most of the approaches required higher computation time. Also, there are deviations in the numerical ranks of the factored solution of CAREs in the Kleinman–Newton approaches, and in all of the cases, RKSM provides significantly better result.

4.3. Stabilization of Eigenvalues

Figures 1–3 illustrate the eigenvalue stabilization of the models mod-606, mod-1998, and mod-2476 by both the RKSM and LRCF-ADI-based Kleinman–Newton techniques.

(a)

(b)

(a)

(b)

(a)

(b)

From Figures 1–3, it has been observed that the eigenvalues of the models mod-606 and mod-2476 are stabilized very well, but it is marginal for the model mod-1998. Thus, it can be concluded that both the RKSM and LRCF-ADI-based Kleinman–Newton techniques have adequate efficiency to stabilize the unstable descriptor systems by closed-loop structure via Riccati-based feedback stabilization.

Figure 4 displays the applicability of the RKSM for the semistable descriptor system, whereas LRCF-ADI-based methods are ineffective in this case. Because the eigenvalues of the semistable model mod-3078 are very close to the imaginary axis, it is not possible to realize in the normal view. So, for the simulation tool capacity and visual convenience, we have considered the a very magnified view of the eigenspaces.

4.4. Stabilization of Step Responses

The investigation of Figures 5 to 12 consist the step responses for some dominant input/output relations to compare the RKSM and the LRCF-ADI-based Kleinman–Newton approaches via the system stabilization. Because the power system models are of input-output relations, there are 16 step responses for each model. For the effective comparison, we have investigated only some graphically significant step responses.

From Figures 5 to 10, it is evident that both RKSM and LRCF-ADI-integrated Kleinman–Newton techniques are applicable for the Riccati-based feedback stabilization of unstable power system models. The graphical comparisons indicate by RKSM is suitably robust. On the other hand, though sometimes the Kleinman–Newton approach provides very good accuracy, it has some scattered behaviors.

Moreover, Figures 11 and 12 signify the applicability of the RKSM technique for the Riccati-based feedback stabilization for the semistable index-1 descriptor systems.

5. Conclusion

From the tabular and graphical comparisons of the results of numerical computations, we have observed that by both RKSM and KN-LRCF-ADI techniques, the CAREs arising from the unstable index-1 descriptor systems are efficiently solved, and the corresponding models are properly stabilized. The semistable index-1 descriptor system is successfully stabilized through Riccati-based feedback stabilization by RKSM, whereas KN-LRCF-ADI is still not suitable for it. There are deviations of the numerical ranks of the factored solutions of CAREs in the Kleinman–Newton approaches, whereas RKSM provides significantly better results for all the cases. RKSM approach has quick convergence ability and occupies very small solution spaces to provide efficient solutions to the CAREs. In contrast, LRCF-ADI-based Kleinman–Newton method has several approaches for finding the solutions of CAREs, where almost all of the approaches required higher computation time. Riccati-based feedback stabilization for the index-1 descriptor systems by the RKSM approach is very effective and robust. Contrariwise, LRCF-ADI-based Kleinman–Newton method is slightly scattered in case of the stabilization of step responses. Thus, it can be concluded that the RKSM is suitably applicable to the unstable index-1 descriptor systems for Riccati-based feedback stabilization, and this method is more preferable to the Kleinman–Newton method in the sense of computation time and memory allocation.

In this work, the MATLAB library command care is used in RKSM to find the solution of the CAREs governed from the reduced-order models and used the direct backward inversion technique for solving shifted sparse linear systems. In the future, we will try to find the self-sufficient RKSM algorithm for solving CAREs by applying the “Restarted Hessenberg Method” with a quick rate of convergence and implement the preconditioning techniques for solving shifted sparse linear systems (see [36–38] and the references therein). Because the LRCF-ADI techniques are inefficient for the semistable systems, we will work for it as well.

Data Availability

The data are available from the following link: https://sites.google.com/site/rommes/software.

Disclosure

This work is under the project "Computation of Optimal Control for Differential-Algebraic Equations (DAE) with Engineering Applications." An earlier version of this manuscript has been presented as a preprint in “Riccati-based feedback stabilization for unstable power system models” and can be found at the link https://arxiv.org/abs/2006.14210.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This project was funded by United International University, Dhaka, Bangladesh. It starts from October 01, 2019, and the reference is IAR/01/19/SE/18.

References

G. M. Huang, L. Zhao, and X. Song, “A new bifurcation analysis for power system dynamic voltage stability studies,” in Proceedings of the IEEE Power Engineering Society Winter Meeting, vol. 2, pp. 882–887, New York, NY, USA, 2002.
View at: Google Scholar
M. G. Angle, S. Madnick, J. L. Kirtley, and S. Khan, “Identifying and anticipating cyberattacks that could cause physical damage to industrial control systems,” IEEE Power and Energy Technology Systems Journal, vol. 6, no. 4, pp. 172–182, 2019.
View at: Publisher Site | Google Scholar
M. Cui and J. Wang, “Deeply hidden moving-target-defense for cybersecure unbalanced distribution systems considering voltage stability,” IEEE Transactions on Power Systems, vol. 36, no. 3, pp. 1961–1972, 2020.
View at: Publisher Site | Google Scholar
I. Kose, Introduction to State-Space Control Theory, Department of Mechanical Engineering, Bogazici University, Istanbul, Turkey, 2003.
D. Groß, C. Arghir, and F. Dörfler, “On the steady-state behavior of a nonlinear power system model,” Automatica, vol. 90, pp. 248–254, 2018.
View at: Publisher Site | Google Scholar
M. A. Mahmud, H. R. Pota, and M. J. Hossain, “Full-order nonlinear observer-based excitation controller design for interconnected power systems via exact linearization approach,” International Journal of Electrical Power & Energy Systems, vol. 41, no. 1, pp. 54–62, 2012.
View at: Publisher Site | Google Scholar
M. Sumon Hossain and M. Monir Uddin, “Iterative methods for solving large sparse Lyapunov equations and application to model reduction of index 1 differential-algebraic-equations,” Numerical Algebra, Control and Optimization, vol. 9, no. 2, pp. 173–186, 2019.
View at: Publisher Site | Google Scholar
P. Benner, J. Saak, and M. M. Uddin, “Reduced-order modeling of index-1 vibrational systems using interpolatory projections,” in Proceedings of the 19th International Conference on Computer and Information Technology (ICCIT), pp. 134–138, IEEE, Dhaka, Bangladesh, 2016.
View at: Publisher Site | Google Scholar
M. M. Uddin, “Computational methods for model reduction of large-scale sparse structured descriptor systems,” Otto-von-Guericke University, Magdeburg, Germany, 2015, Ph.D. dissertation.
View at: Google Scholar
S. Gugercin, T. Stykel, and S. Wyatt, “Model reduction of descriptor systems by interpolatory projection methods,” SIAM Journal on Scientific Computing, vol. 35, no. 5, pp. B1010–B1033, 2013.
View at: Publisher Site | Google Scholar
P. Benner, J. Saak, and M. M. Uddin, “Structure preserving model order reduction of large sparse second-order index-1 systems and application to a mechatronics model,” Mathematical and Computer Modelling of Dynamical Systems, vol. 22, no. 6, pp. 509–523, 2016.
View at: Publisher Site | Google Scholar
E. K. W. Chu, “Solving large-scale algebraic Riccati equations by doubling,” in Proceedings of the Seventeenth Conference of the International Linear Algebra Society, vol. 22, Braunschweig, Germany, 2011.
View at: Google Scholar
W. Chen and L. Qiu, “Linear quadratic optimal control of continuous-time LTI systems with random input gains,” IEEE Transactions on Automatic Control, vol. 61, no. 7, pp. 2008–2013, 2016.
View at: Publisher Site | Google Scholar
H. Abou-Kandil, G. Freiling, V. Ionescu, and G. Jank, Matrix Riccati Equations in Control and Systems Theory, Springer, Berlin, Germany, 2012.
M. Uddin, M. A. H. Khan, and M. M. Uddin, “Riccati based optimal control for linear quadratic regulator problems,” in Proceedings of the 5th International Conference on Advances in Electrical Engineering (ICAEE), pp. 290–295, IEEE, Dhaka, Bangladesh, 2019.
View at: Google Scholar
M. Uddin, M. A. H. Khan, and M. M. Uddin, “Efficient computation of Riccati based optimal control for power system models,” in Proceedings of the 22nd International Conference of Computer and Information Technology (ICCIT), pp. 260–265, IEEE, Dhaka, Bangladesh, 2019.
View at: Google Scholar
V. Simoncini, “A new iterative method for solving large-scale Lyapunov matrix equations,” SIAM Journal on Scientific Computing, vol. 29, no. 3, pp. 1268–1288, 2007.
View at: Publisher Site | Google Scholar
V. Druskin and V. Simoncini, “Adaptive rational Krylov subspaces for large-scale dynamical systems,” Systems & Control Letters, vol. 60, no. 8, pp. 546–560, 2011.
View at: Publisher Site | Google Scholar
V. Simoncini, “Analysis of the rational Krylov subspace projection method for large-scale algebraic Riccati equations,” SIAM Journal on Matrix Analysis and Applications, vol. 37, no. 4, pp. 1655–1674, 2016.
View at: Publisher Site | Google Scholar
P. Benner, J. R. Li, and T. Penzl, “Numerical solution of large-scale Lyapunov equations, Riccati equations, and linear-quadratic optimal control problems,” Numerical Linear Algebra with Applications, vol. 15, pp. 1–23, 2008.
View at: Publisher Site | Google Scholar
P. Kürschner, “Efficient low-rank solution of large-scale matrix equations,” Shaker Verlag Aachen, Düren, Germany, 2016, Ph.D. dissertation.
View at: Google Scholar
S. L. Varricchio, F. D. Freitas, and N. Martins, “Hybrid modal‐balanced truncation method based on power system transfer function energy concepts,” IET Generation, Transmission & Distribution, vol. 9, no. 11, pp. 1186–1194, 2015.
View at: Publisher Site | Google Scholar
J. R. Li, “Model reduction of large linear systems via low rank system gramians,” Massachusetts Institute of Technology, Cambridge, MA, USA, 2000, Ph.D. dissertation.
View at: Google Scholar
T. Penzl, “A cyclic low-rank smith method for large sparse Lyapunov equations,” SIAM Journal on Scientific Computing, vol. 21, no. 4, pp. 1401–1418, 1999.
View at: Publisher Site | Google Scholar
P. Benner, P. Kürschner, and J. Saak, “Efficient handling of complex shift parameters in the low-rank Cholesky factor ADI method,” Numerical Algorithms, vol. 62, no. 2, pp. 225–251, 2013.
View at: Publisher Site | Google Scholar
P. Benner and T. Stykel, “Numerical solution of projected algebraic Riccati equations,” SIAM Journal on Numerical Analysis, vol. 52, no. 2, pp. 581–600, 2014.
View at: Publisher Site | Google Scholar
P. Benner, P. K¨urschner, and J. Saak, “Self-generating and efficient shift parameters in ADI methods for large Lyapunov and sylvester equations,” Electronic Transactions on Numerical Analysis, vol. 43, pp. 142–162, 2014.
View at: Google Scholar
O. Abidi, M. Hached, and K. Jbilou, “Adaptive rational block Arnoldi methods for model reductions in large-scale MIMO dynamical systems,” New Trends in Mathematical Science, vol. 4, no. 2, p. 227, 2016.
View at: Publisher Site | Google Scholar
Y. Lin and V. Simoncini, “Minimal residual methods for large scale Lyapunov equations,” Applied Numerical Mathematics, vol. 72, pp. 52–71, 2013.
View at: Publisher Site | Google Scholar
M. M. Uddin, Computational Methods for Approximation of Large-Scale Dynamical Systems, Chapman and Hall/CRC Press, Boca Raton, FL, USA, 2019.
P. Benner and T. Stykel, “Model order reduction for differential-algebraic equations: a survey,” Surveys in Differential-Algebraic Equations IV, Springer, Berlin, Germany, 2017.
View at: Publisher Site | Google Scholar
F. Feitzinger, T. Hylla, and E. W. Sachs, “Inexact kleinman-Newton method for Riccati equations,” SIAM Journal on Matrix Analysis and Applications, vol. 31, no. 2, pp. 272–288, 2009.
View at: Publisher Site | Google Scholar
N. Martins, P. C. Pellanda, and J. Rommes, “Computation of transfer function dominant zeros with applications to oscillation damping control of large power systems,” IEEE Transactions on Power Systems, vol. 22, no. 4, pp. 1657–1664, 2007.
View at: Publisher Site | Google Scholar
F. D. Freitas and A. Simoes Costa, “Computationally efficient optimal control methods applied to power systems,” IEEE Transactions on Power Systems, vol. 14, no. 3, pp. 1036–1045, 1999.
View at: Publisher Site | Google Scholar
R. B. Leandro, A. S. e Silva, I. C. Decker, and M. N. Agostini, “Identification of the oscillation modes of a large power system using ambient data,” Journal of Control, Automation and Electrical Systems, vol. 26, no. 4, pp. 441–453, 2015.
View at: Publisher Site | Google Scholar
W.-H. Luo, T.-Z. Huang, L. Li, Y. Zhang, and X.-M. Gu, “Efficient preconditioner updates for unsymmetric shifted linear systems,” Computers & Mathematics with Applications, vol. 67, no. 9, pp. 1643–1655, 2014.
View at: Publisher Site | Google Scholar
X.-M. Gu, T.-Z. Huang, G. Yin, B. Carpentieri, C. Wen, and L. Du, “Restarted Hessenberg method for solving shifted nonsymmetric linear systems,” Journal of Computational and Applied Mathematics, vol. 331, pp. 166–177, 2018.
View at: Publisher Site | Google Scholar
X. M. Gu, T. Z. Huang, B. Carpentieri, A. Imakura, K. Zhang, and L. Du, “Efficient variants of the CMRH method for solving a sequence of multi-shifted non-Hermitian linear systems simultaneously,” Journal of Computational and Applied Mathematics, vol. 375, Article ID 112788, 2020.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2021 Mahtab Uddin et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

297

Downloads

427

Citations