A New Algorithm for Solving Large-Scale Generalized Eigenvalue Problem Based on Projection Methods

In this paper, we consider four methods for determining certain eigenvalues and corresponding eigenvectors of large-scale generalized eigenvalue problems which are located in a certain region. In these methods, a small pencil that contains only the desired eigenvalue is derived using moments that have obtained via numerical integration. Our purpose is to improve the numerical stability of the moment-based method and compare its stability with three other methods. Numerical examples show that the block version of the moment-based (SS) method with the Rayleigh–Ritz procedure has higher numerical stability than respect to other methods.


Introduction
Many problems arising in different fields of science and engineering can be reduced to the generalized eigenvalue problem [1][2][3]: where A, B are n × n real or complex, large, sparse, and only a few of the eigenvalues are desired. Also, when is B � I (identity matrix), we have a standard eigenvalue problem. Computing eigenpairs (λ, x) of the generalized and standard eigenvalue problems is one of the important problems in many scientific applications [4][5][6][7]. ere are several methods for solving such eigenvalue problems [8]. Among these methods, the iterative methods are used to generate a subspace that contains the desired eigenvectors. Techniques based on the Krylov subspaces are powerful tools for building desired subspaces for large-scale eigenvalue problems [9][10][11]. Expressed methods in this article find all of the zeros that lie in a circle using numerical integration. In this paper, we briefly describe moment-based method in Section 2, Rayleigh-Ritz with contour integral method in Section 3, block version of the Sakurai-Sugiura method in Section 4, and block version of the SS method with Rayleigh-Ritz procedure in Section 5 for solving generalized eigenvalue problem (1). In Section 6, we provide four numerical tests for comparing four methods, and in Section 6, we apply the BSSRR method with selected matrices from different application areas, and finally, we draw some conclusions in Section 7.

Moment-Based Method (SS Method)
For solving (1), we consider computing entire poles of a rational function: ose are eigenvalues λ of equation (1) and lie in a circle using numerical integration. Let Γ be positively oriented closed Jordan curve [12] in the complex plane and λ 1 , . . . , λ n be distinct eigenvalues that lie in the Γ. Let  (4) en, we have the following theorem.
Proof. In [13], by approximating the integral of equation (3) via the N-point trapezoidal rule, we obtain Let ξ 1 , . . . , ξ m be the eigenvalues of pencil H < m − λH m . We regard λ j � c + ξ j , 1 ≤ j ≤ m as the approximations for λ 1 , . . . , λ m and Also, let V m be the Vandermonde matrix for ξ 1 , . . . , ξ m . en, the approximations for the eigenvectors are obtained by

Rayleigh-Ritz with Contour Integral Method (CIRR Method)
We consider (1), let A, B ∈ R n×n be symmetric and let B be positive definite and (λ j , x j ), 1 ≤ j ≤ n be eigenpairs of the matrix pencil (A, B). We apply a Rayleigh-Ritz procedure with an orthonormal basis Q ∈ R n×m . e projected matrices are given by A � Q T AQ and B � Q T BQ. Q ∈ R n×m is used to generate a sequence of subspace containing approximations to the desired eigenvector. e Ritz values of the projected pencil (A, tB) are taken as approximate eigenvalues for original pencil (A, B) with corresponding Ritz vectors. In this method, by applying the Rayleigh-Ritz procedure moments are not explicitly used [17]. e algorithm is as follows.
Rayleigh-Ritz procedure (1) Construct an orthonormal basis Q Theorem 2. Let s k be defined by (4). Suppose that ]is expanded by the eigenvectors en, Proof. It follows from (4) and (8) that Since (λ j , x j ) is an eigenpair of the matrix pencil (A, B), By the residue theorem, we obtain the result. We define the m × m Vandermonde matrix with From the equation (9), we have Proof. Since λ 1 , . . . , λ m are mutually distinct and α j ≠ 0 for 1 ≤ j ≤ m. V and D are nonsingular. erefore, it follows from (13) that Since the vectors q 1 , . . . , q m are orthonormal basis of span s 0 , . . . , s m− 1 , equation (14) holds.

□
For nonzero vector v ∈ R n , we define the moments: where c is located inside Γ. Also, we obtain the following approximations via the N-point trapezoidal rule:

Block Sakurai-Sugiura Method (BSS Method)
In this method, for solving (1), we reformulate the SS method in the context of the resolvent theory. is method has the potential to resolve degenerated eigenvalues.

Theorem 4.
Let zB − A be a regular pencil of order N. en, there exist nonsingular matrices P, Q ∈ C N×N such that where J i , N i ∈ C k i ×k i are Jordan blocks, N i is nilpotent, and I k denotes the identity matrix of order k. Proof. In [12].
Here, because P, Q are the regular matrices, we can define P � P − 1 and Q � Q − 1 . According to the (18), we will partition row vectors in P and Q into P i , Q i ∈ C k i ×N , and column vectors in P and Q into P i , Q i ∈ C N×k i , respectively, for i � 1, . . . , r.
where α i is an eigenvalue of Jordan block J i .
Using the resolvent of the Jordan block, And

Theorem 6. e localized moment matrix is written as
Proof. In [18].
Proof. In [18].  Proof. By choosing row vectors of C H and vectors of D to be As for the rank of Q Γ D, we consider that column vectors of Q Γ D from the Krylov series of J Γ starting from Q Γ v. Because J Γ has not degenerated, and elements of Q Γ v are nonzero, these column vectors are linearly independent, and thus the rank of Q Γ D is f k Γ .

Block Version of the SS Method with the Rayleigh-Ritz Procedure (BSSRR Method)
We suggest a new algorithm for computing all poles of analytic function (2) with the use of the algorithm in [19]. As the eigenpairs (λ i , x i ) of equation (1) can be obtained from (25) en, the block version of the SS method with Rayliegh-Ritz procedure [20] constructs the LM-dimensional subspace: For the Rayleigh-Ritz procedure, the subspace δ M contains all eigenvectors of (1) δ M � span x 1 , x 2 , . . . , x m for m ≤ LM. With using N-point trapezoidal rule for equation (25), we have where z j is the quadrature point and ω j is the corresponding weight [21]. Based on eorem 4, we analyse the relationship between the contour integral spectral projection and the Krylov subspace.

An Arnoldi-Based Interpretation of the Contour Integral
Spectral Projection. Since the matrices P, Q are nonsingular, we define P � P − 1 , Q � Q − 1 . According to the Jordan block structure of (18), we partition row vectors in P, Q into P i , Q i ∈ C n i ×n and column vectors in P, Q into P i , Q i ∈ C n×n i , respectively, for i � 1, 2, . . . , r. en, we can derive the following lemma and theorem. where Proof. From eorem 6 and the binomial theorem we have the following relation: Here, since erefore, Lemma 1 is proved.
en, the subspace κ * m (A, V) is defined by the sum of the Krylov subspaces, i.e., Proof. From the definition of S k (25) and Lemma 1, we have erefore, eorem 10 is proved.

□
Remark 1. eorem 10 shows that the block version of the SS method with the Rayleigh-Ritz procedure can be regarded as the Rayleigh-Ritz procedure based on the block Krylov subspace κ * M (S Γ , P Γ V). Here, we note that, in the block version of the SS method with the Rayleigh-Ritz procedure, the basis vectors of κ * M (S Γ , P Γ V) are explicitly computed by (25) and the QR decomposition of S (Algorithm 1).

Numerical Experiments
In  Table 1.

Mathematical Problems in Engineering
Example 2. We let that A, B were complex, random matrices, and B was positive definite. After applying Algorithms 1-4, we obtained numerical results that have been shown in Table 2. Also, the relative residual for described methods has been drawn in Figure 1 for n � 1000.
Example 3. In this example, A, B were taken sparse, symmetric, and random, and A was positive definite. After applying Algorithms 1-4, we obtained numerical results that have been shown in Table 3. Also, the relative residual for described methods has been drawn in Figure 2 for n � 1000.
Example 4. We consider matrices: After applying Algorithms 1-4, we obtained numerical results that have been shown in Table 4. Also, the relative residual for described methods has been drawn in Figure 3 for n � 1000.  Example 5. In this example, we selected seventeen matrices from the UF sparse matrix collection. Two major requirements were used in the selection procedure: matrices with different parameters and matrices arising in different application areas were chosen. We consider the following symptoms parameters. e order N, the number of nonzero elements NZ, and the condition number CON. e application areas of the selected matrices are listed in Table 5. Matrices from the different areas were selected and thus obtained results by running the matrices will be typical in several scientific fields. We applied the BSSRR method for the calculation of relative residual generalized eigenvalue problem (1) when A is one of the selected matrices in Table 5 and B is identity matrix. As the dimension A, B is equal. We computed relative residual respect to − 2 In Table 6 and computed respect to other norms in Table 7, too, the number of eigenpairs for each matrix was sixteen in Table 7.

General Conclusions and Plans for Future Work
Several specific conclusions were drawn in connection with the numerical results presented in the previous section. Some general conclusions are given as follows: (1) All numerical experiments indicate that CIRR, BSS, and BSSRR methods have higher stability than respect to the SS method (2) BSSRR method has less relative residual respect to SS, CIRR, and BSS methods (3) If − 2 is used for calculation of relative residual in the BSSRR method, then we have higher accuracy and less consuming time Designing quadrature points with higher performance and a more precise error analysis of the BSSRR method is a part of our future work.

Data Availability
No data were used to support this study.

Conflicts of Interest
e authors declare that they have no conflicts of interest. Mathematical Problems in Engineering 9