Effect of Shift Parameters in Rational Krylov Subspace Method for Solving Riccati Equations Arise from Power System Models

The simulation of continuous-time algebraic Riccati equations (CARE) derived from the very large power system models is a highly laborious task and most cases infeasible due to the sophisticated structural ingredients. The computation is very time costly and the rate of convergence can be severely affected in the direct solvers. To overcome those adversities, an iterative approach Rational Krylov Subspace Method (RKSM) is introduced to deal with those large-scale CAREs. The solutions of those CAREs and hence the optimal feedback matrices can be efficiently explored by the RKSM approach to stabilize the power system models of unstable category. In this approach, shift parameters play a vital role in the convergence of the computations and size of the solution spaces. The goal of the work is to investigate the effect of different types of shift parameters on the stabilization process. To attain the mentioned objective, a modified version of the iterative RKSM algorithm is proposed by employing two types of shift parameters, namely, the adaptive ADI shifts and heuristic shifts. Qualitative discussions for the outcomes for those shift parameters are narrated by tabular and figurative methods.

In the system (1), is the differential coefficient matrix and , , , are the state submatrices, respectively, whereas , are the control multiplier sub-matrices, , are the state multiplier sub-matrices, and is the direct gain for with . In the power system models, the matrix is zero or absent. Also, output. The mentioned sub-matrices and are invertible (Hossain and Uddin, 2019).
Over the years, computational techniques are modified for the feasibility of massive systems. But computational intricacy and memory extravagant simulation tolls keep the approaches impractical. Thus, large dimensional empirical systems are needed to be reduced to smaller dimensional systems. The strategy of finding a reduced-order model (ROM) of a given system is called model order reduction (MOR), currently earning vital attention to the researchers. ROM-based techniques are extensively adopted in the technological fields of the modern scientific world (Flagg et al., 2010). An adjustable ROM size and keeping the structural design invariant of the test system is the earnest sought (Rahman et al., 2020). The mechanism of the MOR techniques needs to essentially robust maintaining minimized norm of the error.
To duly synchronize the iterative steps of the manipulation, the index-1 descriptor system (1) is to be restructured by implementing the well-known Schur-complements as , .
Then, the index-1 descriptor system (1) can be converted to an equivalent generalized LTI continuous-time system as The essence of LTI continuous-time systems is undeniable in the phenomena of the real-world engineering models, such as mechatronics, aero-nautics, and system and control theory (Rahman et al., 2021).
In the analysis of system stability and its appliances, the Continuous-time Algebraic Riccati Equation Riccati-based feedback stabilization is obligate for an unstable type of system (3), and the optimal feedback matrix is the premier ingredient in this process (Chen and Qui, 2015). To stabilize the target system optimally, the matrix is to be replaced by . Then, the stabilized system can be formed as Computationally proficient solvers or analytical tools for large-scale CARE arising from the index-1 descriptors systems are still not technically cheap. Some Kleinman-Newton methods exist, which are very complicated, time laborious, and predefined structures are required (Mena and Saak, 2008). Alternating direction implicit (ADI) based approaches is computationally extravagant for the prerequisite large solution space, preconditioned initial parameters, and time-costly matrix factorization (Hasan and Uddin, 2018).
The linearization ability and enforcement of initial priories boost up the convergence rate of the simulations by the projection-based approach RKSM, which is efficiently applicable in systems with perturbations (Simoncini et al., 2013). We are extending the techniques of the RKSM approach for the standard systems that is discussed in (Simoncini, 2016). The extension is including the modification of the existing work in such a way that it can apply for index-1 descriptor systems. We will investigate the effect of the shift parameters in the RKSM algorithm and compare this for the adaptive ADI shifts and heuristic shifts, respectively.
Also, to optimally stabilize the unstable index-1 descriptor system, we will apply the Riccati-based feedback matrix and analyze the stabilized system through eigenspace and transient behaviors.

Rational Krylov subspace method for solving index-1 descriptor systems
If the eigenvalues of the matrix pair ( ) satisfy ̅ for the all that ensures the solution of the CARE (4) exists and unique. For a given set of competent shift parameters , the requisite dimensional rational Krylov subspace can be generated to pursue the orthogonal projector . The desired rational Krylov subspace can be constructed as Again, consider the CARE (4) and apply the Galerkin condition to it. After the simplification by matrix algebra, a low-rank CARE can be achieved as Here, ̂ , ̂ , ̂ , ̂ , and ̂ . The equation (6) is a low-rank CARE and can be solved by MATLAB care command or any existing methods, such as the Schur-decomposition method.
For the quick and smooth convergence of the proposed approach, adjustable shift selection is crucial (Druskin et al., 2010). There are many schemes to find suitable shift parameters, several approaches are given in (Moret and Popolizio, 2014;Benner et al., 2010) and references therein. In the present work, we are adopting the adaptive ADI shift and the heuristic shift parameters approach for index-1 descriptor systems.
In the iterative simulation process, to generate the feasible orthogonal projector the existing set of shift parameters must be extended within the solution space. By the solution ̂ of the lowrank CARE (6), the solution of the CARE (4) will be approximated as ̂ .
To have a stopping condition, a suitable approximation of the residual is required. The residual of th iteration can be executed in terms of the Frobenius norm ‖ ‖ as Here the block upper triangular matrix is found from the QR-factorization of the matrix is defined as The factor needs to be computed from the QR-factorization . Then, the relative residual can be estimated as The discussion about the terms and notations of the equation (8) are given in (Uddin, 2020) and the references therein.

Sparsity preservation techniques
In system (3), the matrix is dense and for this reason, the rate of convergence of their structured system will be disturbed and accuracy will be affected. To avoid those inconveniences, instead o f the explicit conversion of , a more convenient approach will be done. To find , at each of the iterations, a simple linear system will be solved as Here is the truncated term. Due to sparse system structure, the higher dimensional linear system (9) is conveniently solvable by the conventional direct solvers (Uddin, 2020). For the improvement of the consistency of the RKSM approach, explicit forms of the reduced-order matrices will not be used to construct the reduced-order system. The reduced-order matrices in sparse form can be attained in the following way The low-rank solution ̂ is symmetric and positive-definite and can be factorized as ̂ . Consider the desired low-rank factored solution of the CARE (4) as . Thus, the solution of the (4) can be attained as such that For the future steps, the factored solution will be stored and hence 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-1.
Finally, applying , the optimally stabilized LTI continuous-time system can be written as (5). To preserve the structure of the system, it needs to back to the original structure (1), and for this, the target system needs to be re-written as

Numerical results
The eligibility and efficiency of the proposed techniques for desired shift parameters are investigated by subjecting a Brazilian Interconnected Power System (BIPS) model, namely BIPS-606 with 606 differential and 6529 algebraic variables, which is of the type unstable index-1 descriptor system (Rommes et al., 2019;Freitas et al., 2008). The proposed techniques are applied with a tolerance level to find low-rank approximate solutions for both adaptive ADI and heuristic shifts.
To stabilize the target model the optimal feedback matrix is used. The comparisons of several computational aspects are given in Table 1. From the above table, it is observed that by using adaptive ADI shifts a smaller solution space can be achieved, which is feasible for the memory allocation of the computational tool. In contrast, heuristic shifts are better for minimizing computational time.
The comparison of the relative residuals and the eigenvalue stabilization for the RKSM approach exerting the adaptive ADI and heuristic shifts are shown in Fig 1 and Fig 2.   Fig 1: Comparison of the relative residuals. From the above figures, it is evident that that adaptive ADI shifts provide better relative residuals than heuristic shifts. But the stabilized eigenvalues for both of the shifts are very identical. From the above sub-figures, it can be said that by utilizing both of the adaptive and heuristic shifts the unstable index-1 descriptor systems can be efficiently stabilized.

CONCLUSION:
In this work, the techniques to optimally stabilize the unstable index-1 descriptor system by the optimal feedback matrix attained by the solution of the Riccati equation, are introduced and embedded. To serve the purpose, the CARE corresponding to a model derived from the Brazilian Interconnected Power System has been efficiently solved by the iterative RKSM approach utilizing adaptive ADI shift and heuristic shift parameters. The robustness and time-dealing efficiency of the computation is justified by applying the techniques to the target model. From the tabular and graphical comparisons, it can be concluded that the RKSM approach utilized both the adaptive ADI and heuristic shifts, stabilized the target model with the desired efficiency. But adaptive ADI shifts are feasible for the number of iterations and space dimensions, whereas heuristic shifts are suitable for the numerical rank and computation time.
Thus, adaptive ADI shifts are applicable for the feasibility of memory allocation and heuristic shifts are comfortable for quick convergence.

ACKNOWLEDGEMENT:
We would like to make manifest our profound gratitude to Md. Motlubar Rahman, Associate Professor in Mathematics, Directorate of Secondary and Higher Education, Ministry of Education, Bangladesh for his invaluable suggestions and support in writing this article successfully.

CONFLICTS OF INTEREST:
The authors declared that there is no conflict of interest in this article.