{'fraction_non_alphanumeric': 0.04815457318896931, 'fraction_numerical': 0.020338577495962195, 'mean_word_length': 4.202271666407344, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
An Eigenwise Parallel Augmented Subspace Method for Eigenvalue Problems

A type of parallel augmented subspace scheme for eigenvalue problems is proposed by using coarse space in the multigrid method. With the help of coarse space in multigrid method, solving the eigenvalue problem in the finest space is decomposed into solving the standard linear boundary value problems and very low dimensional eigenvalue problems. The computational efficiency can be improved since there is no direct eigenvalue solving in the finest space and the multigrid method can act as the solver for the deduced linear boundary value equations. Furthermore, for different eigenvalues, the corresponding boundary value problem and low dimensional eigenvalue problem can be solved in the parallel way since they are independent of each other and there exists no data exchanging. This property means that we do not need to do the orthogonalization in the highest dimensional spaces. This is the main aim of this paper since avoiding orthogonalization can improve the scalability of the proposed numerical method. Some numerical examples are provided to validate the proposed parallel augmented subspace method.


Introduction
Solving large scale eigenvalue problems is one of fundamental problems in modern science and engineering society. It is always a very difficult task to solve high-dimensional eigenvalue problems which come from practical physical and chemical sciences. Different from the case of boundary value problems, there is no many efficient numerical methods for solving eigenvalue problems with optimal complexity. The large scale eigenvalue problems pose significant challenges for scientific computing. In order to solve these large scale sparse eigenvalue problems, Krylov subspace type methods (Implicitly Restarted Lanczos/Arnoldi Method (IRLM/IRAM) [20]), the Preconditioned INVerse ITeration (PINVIT) method [21,22,23], the Locally Optimal Block Preconditioned Conjugate Gradient (LOBPCG) method [24,25], and the Jacobi-Davidson-type techniques [26] have been developed. All these popular methods include the orthogonalization step which is a bottleneck for designing efficient parallel schemes. Recently, a type of multilevel correction method is proposed for solving eigenvalue problems in [15,29,30]. In this multilevel correction scheme, an augmented subspace is constructed with the help of coarse space from the multigrid method. The application of this augmented subspace leads to that the solution of eigenvalue problem on the final level mesh can be reduced to a series of solutions of boundary value problems on the multilevel meshes and a series of solutions of the eigenvalue problem on the low dimensional augmented subspace. The multilevel correction method gives a way to construct the multigrid method for eigenvalue problems.
It is well known that the multigrid method [8,19,31] provides an optimal numerical method for linear elliptic boundary value problems. The error bounds of the approximate solution obtained from these efficient numerical algorithms are comparable to the theoretical bounds determined by the finite element discretization, while the amount of computational work involved is only proportional to the number of unknowns in the discretized equations. For more details of the multigrid method, please refer to [4,5,6,7,8,13,14,16,18,19,31,32] and the references cited therein.
This paper aims to design a type of eigenvalue-wise parallel method for eigenvalue problems with the help of the coarse space from the multigrid method. It is well known that there exist many work considering the applications of multigrid method for eigenvalue problems. For example, there have existed applications of the multigrid method to the PINVIT and LOBPCG methods. But, in these applications, the multigrid method only acts as the precondition for the included linear equations. This means that the multigrid method only improves the efficiency of the inner iteration and does not change the outer iteration. Unfortunately, in these state-of-the-art, the applications of multigrid method do not deduce a new eigensolver. The idea of designing the eigenvalue-wise parallel augmented subspace method for eigenvalue problems is based on the combination of the multilevel correction method [15,29,30] and parallel computing technique. With the help of coarse space in multigrid method, the eigenvalue problem solving is transformed into a series of solutions of the corresponding linear boundary value problems on the sequence of finite element spaces and eigenvalue problems on a very low dimensional augmented space. Further, in order to improve the parallel scalability, the correction process for different eigenpair is executed independently. This property means there is no orthogonalization in the highest dimensional space which account for a large portion of wall time in the parallel computation.
An outline of the paper goes as follows. In Section 2, we introduce the finite element method for the eigenvalue problem and the corresponding basic error estimates. A type of parallel augmented subspace method for solving the eigenvalue problem by finite element method is given in Section 3, and the corresponding computational work estimate are given in Section 4. In Section 5, four numerical examples are presented to validate our theoretical analysis. Some concluding remarks are given in the last section.

Finite element method for eigenvalue problem
This section is devoted to introducing some notation and the standard finite element method for the eigenvalue problem. In this paper, we shall use the standard notation for Sobolev spaces W s,p (Ω) and their associated norms and semi-norms (cf. [1]). For p = 2, we denote H s (Ω) = W s,2 (Ω) and H 1 0 (Ω) = {v ∈ H 1 (Ω) : v| ∂Ω = 0}, where v| Ω = 0 is in the sense of trace, · s,Ω = · s,2,Ω . In some places, · s,2,Ω should be viewed as piecewise defined if it is necessary. The letter C (with or without subscripts) denotes a generic positive constant which may be different at its different occurrences through the paper.
For simplicity, we consider the following model problem to illustrate the main idea: Find (λ, u) such that −∇ · (A∇u) = λu, in Ω, u = 0, on ∂Ω, where A is a symmetric and positive definite matrix with suitable regularity, Ω ⊂ R d (d = 2, 3) is a bounded domain with Lipschitz boundary ∂Ω.
In order to use the finite element method to solve the eigenvalue problem (2.1), we need to define the corresponding variational form as follows: Find (λ, u) ∈ R × V such that a(u, u) = 1 and a(u, v) = λb(u, v), ∀v ∈ V, (2.2) where V := H 1 0 (Ω) and The norms · a and · b are defined by It is well known that the eigenvalue problem (2.2) has an eigenvalue sequence {λ j } (cf. [2,9]): and associated eigenfunctions where a(u i , u j ) = δ ij (δ ij denotes the Kronecker function). In the sequence {λ j }, the λ j are repeated according to their geometric multiplicity. For our analysis, recall the following definition for the smallest eigenvalue (see [3,9]) Now, let us define the finite element approximations of the problem (2.2). First we generate a shape-regular triangulation T h of the computing domain Ω ⊂ R d (d = 2, 3) into triangles or rectangles for d = 2 (tetrahedrons or hexahedrons for d = 3). The diameter of a cell K ∈ T h is denoted by h K and the mesh size h describes the maximum diameter of all cells K ∈ T h . Based on the mesh T h , we can construct a finite element space denoted by V h ⊂ V . For simplicity, we set V h as the linear finite element space which is defined as follows where P 1 denotes the linear function space.
The standard finite element scheme for eigenvalue problem (2.2) is: From [2,3], the discrete eigenvalue problem (2.6) has eigenvalues: and corresponding eigenfunctionsū . From the min-max principle [2,3], we have the following upper bound result Let M (λ i ) denote the eigenspace corresponding to the eigenvalue λ i which is defined by w is an eigenfunction of (2.2) corresponding to λ i , (2.8) and define Let us define the following quantity: where T : In order to understand the method more clearly, we provide a detailed process to derive the error estimate for the eigenpair approximation by the finite element method. For this aim, we define the finite element projection P h as follows It is obvious that Before stating error estimates of the subspace projection method, we introduce a lemma which comes from [27]. For completeness, a proof is also provided here.
The following Rayleigh quotient expansion of the eigenvalue error is the tool to obtain the error estimates of the eigenvalue approximations.
is an eigenpair of the eigenvalue problem (2.2). Then for any w ∈ V \{0}, the following expansion holds: The following theorem is similar to the corresponding results in [3,9]. In order to analyze and understand the proposed numerical algorithms in this paper, the following error estimates include the explicit constants. Since this paper considers the parallel scheme for different eigenpair, we state the error estimates for any single eigenpair.
Theorem 2.1. Let (λ, u) denote an exact eigenpair of the eigenvalue problem (2.2). Assume the eigenpair approximation (λ i,h ,ū i,h ) has the property thatμ i,h = 1/λ i,h is closest to µ = 1/λ. The corresponding spectral projectors E i,h : V → span{ū i,h } and E : V → span{u} are defined as follows and a(Ew, u i ) = a(w, u i ), ∀w ∈ V.
Then the following error estimate holds

14)
where η(V h ) is define in (2.10) and δ λ,h is defined as follows Furthermore, the eigenvector approximationū i,h has following error estimate in L 2 -norm Proof. Similarly to the duality argument in the finite element method, the following inequality holds Since (I − E i,h )P h u ∈ V h and a((I − E i,h )P h u,ū i,h ) = 0, the following orthogonal expansion holds where α j = a(P h u,ū j,h ). From Lemma 2.1, we have From the property of eigenvectorsū 1,h , · · · ,ū N h ,h , the following equalities hold which leads to the following property From (2.6) and definitions of eigenvectorsū 1,h , · · · ,ū N h ,h , we have following equalities for normal orthogonal basis functions Then from (2.18), (2.19), (2.20) and (2.21), the following estimates hold This is the desired result (2.14).
Similarly, from (2.18), (2.19), (2.20) and (2.21), the following estimates hold Combining (2.17) and (2.23) leads to following inequalities From (2.17), (2.24) and the triangle inequality, we have the following error estimate for the eigenvector approximation in the L 2 -norm This is the second desired result (2.16) and the proof is complete.
In order to analyze the method which will be given in this paper, we state some error estimates in the following theorem.
Theorem 2.2. Under the conditions of Theorem 2.1, the following error estimates hold

29)
and Based on the error estimates in Theorem 2.1, the property u a = ū i,h a = 1 and (2.31), we have following estimations Then the desired result (2.26) can be deduced by the combination of (2.14) and (2.32).
We now investigate the distance of P h u fromū i,h . First, the following estimate holds Combining (2.37) and the triangle inequality leads to following inequalities which in turn implies that This completes the proof of the desired result (2.28).

Parallel augmented subspace method
In this section, we will propose the eigenvalue-wise parallel augmented subspace method for eigenvalue problems based on the multilevel correction scheme [15,29,30]. With the help of the coarse space in multigrid method, the method can transform the solution of the eigenvalue problem into a series of solutions of the corresponding linear boundary value problems on the sequence of finite element spaces and eigenvalue problems on a very low dimensional augmented space. For different eigenpairs, we can do the correction process independently and it is not necessary to do orthogonalization in the finest level of fine finite element space. Thus the proposed algorithm has a good scalability. Since the eigenvalue problems are only solved in a low dimensional space, the numerical solution in this new version of augmented subspace method is not significantly more expensive than the solution of the corresponding linear boundary value problems.
In order to describe the parallel augmented subspace method clearly, we first introduce the sequence of finite element spaces. We generate a coarse mesh T H with the mesh size H and the coarse linear finite element space V H is defined on the mesh T H . Then we define a sequence of triangulations T h k of Ω ⊂ R d as follows. Suppose that T h1 (produced from T H by some regular refinements) is given and let T h k be obtained from T h k−1 via one regular refinement step (produce β d subelements) such that where positive number β > 1 denotes the refinement index. Based on this sequence of meshes, we construct the corresponding nested linear finite element spaces such that The sequence of finite element spaces V h1 ⊂ V h2 ⊂ · · · ⊂ V hn and the finite element space V H have the following relations of approximation accuracy Proposition 3.1. For simplicity of theoretical analysis, we assume the domain Ω is convex in this paper. The standard error estimates [8,11,27] for the linear finite element method implies where C is the constant independent from the mesh size and eigenpair (λ i , u i ) of (2.2).

One correction step and efficient implementation
In order to design the eigenvalue-wise parallel augmented subspace method, we first introduce an one correction step in this subsection.
Assume we have obtained an eigenpair approximations (λ for a certain exact eigenpair. Now we introduce a type of iteration step as follows to improve the accuracy of the given eigenpair approximation (λ For simplicity of notation, we assume that the eigenvalue gap δ λ,h has a uniform lower bound which is denoted by δ λ (which can be seen as the "true" separation of the eigenvalue λ from others) in the following parts of this paper. This assumption is reasonable when the mesh size H is small enough. We refer to [17,Theorem 4.6]and Theorem 2.1 in this paper for details on the dependence Algorithm 1: One Correction Step 1. Define the following linear boundary value problem: Find u Solve (3.6) by some multigrid steps to obtain a new eigenfuction u } and solve the following eigenvalue Summarize the above two steps by defining of error estimates on the eigenvalue gap. Furthermore, we also assume the concerned eigenpair approximation (λ h k , u h k ) is closet to the exact eigenpair (λ h k ,ū h k ) of (2.6) and (λ, u) of (2.2) in this section.
for some constant C 1 . The multigrid iteration for the linear equation (3.6) has the following uniform contraction rate with θ < 1 independent from k and ℓ.
Then the eigenpair approximation (λ where the constants γ,C λ andD λ are defined as follows (3.14) Using (3.10) and (3.16), we deduce that The eigenvalue problem (3.8) can be seen as a low dimensional subspace approximation of the eigenvalue problem (2.6). Using (2.28), Theorems 2.1, 2.2, and their proof, we obtain that and Then we have the desired results (3.11) and (3.12) and conclude the proof.
Remark 3.1. Definition (3.13), Theorems 2.1 and 2.2 imply that γ is less than 1 when η(V H ) is small enough. If λ is large or the spectral gap δ λ is small, then we need to use a smaller η(V H ) or H. Furthermore, we can increase the multigrid steps to reduce θ and then γ. These theoretical restrictions do not limit practical applications where (in numerical implementations), H is simply chosen (just) small enough so that the number of elements of corresponding coarsest space (just) exceeds the required number of eigenpairs (H and the coarsest space are adapted to the number of eigenpairs to be computed).
We would like to point out that the given eigenpair (λ is not necessary to be the one corresponding the smallest eigenvalue. So when we need to solve more than one eigenpairs, the one correction step can be carried out independently for every eigenpair and there exists no data exchanging. This property means that we can avoid doing the time-consuming orthogonalization in the high dimensional space V h k . Now, let us give details for the second step of Algorithm 1. Solving the eigenvalue problem (3.8) provides several eigepairs. Since the desired eigenvalue maybe not the first (smallest) one, we should choose the suitable or the desired eigenpair from the ones of (3.8). Now, let us consider the details to choose the desired eigenpair which has the best accuracy among all the eigenpairs of eigenvalue problem (3.8). For this aim, we come to consider the matrix version of the small scaled eigenvalue problem (3.8). Let N H and {φ k,H } 1≤k≤NH denote the dimension and Lagrange basis functions for the coarse finite element space V H . The function in V H,h k can be denoted by u H,h k = u H + α k u h k . Solving eigenvalue problem (3.8) is to obtain the function u H ∈ V H and the value α k ∈ R. Let u H = NH k=1 u k φ k,H and define the vector u H as u H = [u 1 , · · · , u NH ] T . Based on the structure of the space V H,h k , the matrix version of the eigenvalue problem (3.8) can be written as follows where u H ∈ R NH , α k ∈ R, column vectors b H,h k and c H,h k , scalars β k and ζ k are defined as follows In the practical calculation, the desired eigenpair (λ ) may not be the eigenpair corresponding to the smallest eigenvalue and we will produce a series of (u H , α k ) by solving eigenvalue problem (3.8). In this case, we need to choose the approximate solution which has the largest component in the direction span{ u h k } which is the desired eigenpair in the one correction step defined by Algorithm 1.
Since there holds After obtaining the approximate solutions (u H , α k ), we only need to calculate u H · b H,h k + α k ζ k for every eigenpair which are obtained by solving (3.8) numerically and then choose the one with the largest absolute value as the desired solution.
In the second step of Algorithm 1, we can use the shift-inverse technique since an approximate solution has been obtained in the previous step. Furthermore, we can use the different level of space to act as the coarse space V H in the one correction step for different eigenvalue.

Parallel augmented subspace method
In this subsection, we introduce an eigenvalue-wise parallel augmented subspace method based on the one correction step defined in Algorithm 1.
Here, we design the parallel method to compute m eigenpair approximations of (2.2). For simplicity, we denote the desired eigenpairs by (λ 1 , u 1 ), · · · , (λ m , u m ) and assume there exist m processes denoted by {P 1 , · · · , P m } for the parallel computing. When the number of processes is not equal to the number of desired eigenparis, in order to improve the parallel efficiency, the distribution of desired eigenparis onto the processes should be equal as far as possible to arrive the load balancing. About this point, we refer to the concerned papers for load balancing. The corresponding parallel augmented subspace algorithm is described in Algorithm 2. From Algorithm 2, the computation for m eigenpairs is decomposed into m processes. In order to make the initial eigenfunction approximation u 1,h1 , · · · , u m,h1 be orthogonal each other, in the first step of Algorithm 2, the eigenvalue problem is solved in the first process. We adopt this strategy since (3.22) is a low dimensional eigenvalue problem corresponding to the one in the finest space. Similar to the idea in the full mulgrid method for boundary value problems, the step 2. (A) in Algorithm 2 is used to give an initial eigenpair approximation in the finest space V hn .
2. For i = 1, · · · , m, do the following multilevel correction steps on the process P i in the parallel way (A). For k = 1, · · · , n − 2, do the following iteration: . For ℓ = 0, · · · , ̟ − 1, do the following one correction steps (B). Do the following iterations on the finest level space V hn : . For ℓ = 0, · · · , ̟ n − 1, do the following one correction steps i,hn ) as the output in the n-th level space V hn .
Algorithm 2 shows the idea to design the parallel method for different eigenpairs. In each process, the main computation in the one correction step defined by Algorithm 1 is to solve the linear equation (3.6) in the fine space V h k . It is an easy and direct idea to use the parallel scheme to solve this linear equation based on the mesh distribution on different processes. This type of parallel method is well-developed and there exist many mature software packages such as Parallel Hierarchy Grid (PHG). But we would like to say this is another sense of parallel scheme and this paper is concerned with the eigenvalue-wise parallel method. These discussion means we can design a two level parallel scheme for the eigenvalue problem solving. After implementing Algorithm 2, the resulting eigenpair approximation (λ i,hn , u i,hn ) has the following error estimates Proof. Define e i,k :=ū i,h k − u i,h k . From step 1 in Algorithm 2, it is obvious e i,1 = 0. Then the assumption (3.9) in Theorem 3.1 is satisfied for k = 1. From the definitions of Algorithms 1 and 2, Theorem 3.1 and recursive argument, the assumption (3.9) holds for each level of space V h k (k = 1, · · · , n) with C 1 =C λ in (3.14). Then the convergence rate (3.11) is valid for all k = 1, · · · , n and ℓ = 0, · · · , ̟ − 1.
For k = 2, · · · , n − 1, by Theorem 3.1 and recursive argument, we have By iterating inequality (3.27), the following inequalities hold Then the combination of Theorem 3.1, (3.28) and Algorithm 2 leads to the following error estimates For such choices of ̟ and ̟ n , we arrive at the desired result (3.24).
Remark 3.2. The proof of Theorem 3.2 implies that the assumption (3.9) in Theorem 3.1 holds for C 1 =C λ in each level of space V h k (k = 1, · · · , n). The structure of Algorithm 2, implies that C λ does not change as the algorithm progresses from the initial space V h1 to the finest one V hn . From the estimate (3.24), it can be observed that the final algebraic accuracy depends strongly on γ ̟n . Furthermore, increasing ̟ on the coarse levels spaces V h2 , · · · , V hn−1 can not improve the final algebraic accuracy. For this reason, we set ̟ = 1 on the coarse level spaces V h2 , · · · , V hn−1 . Now we briefly analyze the orthogonality of different eigenfunction obtained by Algorithm 2. Suppose u i,hn =ū i,hn + r i and u j,hn =ū j,hn + r j corresponding toλ i,hn =λ j,hn . By Theorem 3.2, we have the error estimates for r i and r j . Furthermore, the orthogonality of u i,hn and u j,hn has following estimate b(u i,hn , u j,hn ) = b(ū i,hn + r i ,ū j,hn + r j ) = b(r i ,ū j,hn ) + b(ū i,hn , r j ) + b(r i , r j ) So Algorithm 2 can keep the orthogonality for different eigenfunction when we do enough correction steps (̟ n is enough large) such that the algebraic accuracy is enough small in the finest space V hn . Theorem 3.3. After implementing Algorithm 2, there exists an eigenpair (λ i , u i ) of (2.2) such that the eigenpair approximation (λ i,hn , u i,hn ) has the following error estimates for i = 1, · · · , m Proof. From Theorems 2.2 and 3.2, (2.26) and (3.23), we have following estimates This is the desired result (3.29).
From (2.26), (2.33), (3.24), (3.25) and (3.29), u i − u i,hn b has the following estimates This is the desired result (3.30). From (2.13) and (3.29), the error estimate for |λ i − λ i,hn | can be deduced as follows Then the desired result (3.31) is obtained and the proof is complete.
From Proposition 3.1 and Theorem 3.3, it is easy to deduce following explicit error estimates for the eigenpair approximation (λ i,hn , u i,hn ) by Algorithm 2.
Corollary 3.1. After implementing Algorithm 2, there exists an eigenpair (λ i , u i ) of (2.2) such that the eigenpair approximation (λ i,hn , u i,hn ) has the following error estimates for i = 1, · · · , m u i − u i,hn a ≤ 2 1 + Ch 2 n + 2βγ ̟n λ i h n , (3.32) where the constant C depends on λ 1 , spectral gap δ λ , γ and β but independent of the mesh size h n .
Remark 3.3. When m = 1, Algorithm 2 becomes a serial algorithm. Even in this case, we can deal with different eigenpair individually, which always has a better efficiency than traditional algorithm when the number of desired eigenpairs is large enough.

Work estimates of parallel augmented subspace method
Now we turn our attention to the estimate of computational work for the parallel augmented subspace scheme defined by Algorithm 2. First, we define the dimension of each level finite element space as N k := dimV h k . Then we have N n , k = 1, 2, · · · , n. Iterating (4.2) and using the fact (4.1), we obtain Remark 4.1. Since there exists no data transfer between different processes, the total computational work of Algorithm 2 is equal to that of one process in Theorem 4.1 for one eigenpair. Further, since γ has a uniform bound from 1 (γ < 1), then we do not need to do many correction steps in each level of finite element space. As in Remark 3.1, we choose ̟ = 1 for k = 2, · · · , n − 2 and ̟ n is dependent on the algebraic accuracy ε. Then the final computational work in each processor should be O(N n | log ε|) and the included constant is independent of the number m of the desired eigenpairs.

Numerical results
In this section, we provide four numerical examples to validate the proposed numerical method in this paper.

The model eigenvalue problem
In this subsection, we use Algorithm 2 to solve the following model eigenvalue problem: Find (λ, u) ∈ R × V such that ∇u 0 = 1 and where Ω = (0, 1) × (0, 1) × (0, 1). In order to check the parallel property of Algorithm 2, we compute the first 200 eigenpairs of (5.1). In this example, we choose H = 1/16, β = 2, ̟ = ̟ n = 1 and n = 5. In the first step of one correction step defined by Algorithm 1, 1 multigrid step with 2 Conjugate-Gradient (CG) steps for pre-and post-smoothing is adopted to solve the linear problem (3.6). Figure 1 shows the corresponding error estimates of |λ i −λ i,hn | for i = 1, · · · , 200 and the CPU time for each eigenpair, respectively. From Figure 1, we can find that Algorithm 2 has the optimal error estimate, and We also test the algebraic errors |λ i,hn − λ i,hn | between the numerical approximation by Algorithm 2 and the exact finite element approximation for the first 20 eigenvalues on the finest level of mesh. The corresponding results are presented in Figure 2 which shows the algebraic accuracy improves with the growth of number of correction steps ̟ n . The performance of Algorithm 2 for computing the first 1000 eigenpairs is also investigated. Figure 3 shows the error estimate and CPU time for each eigenvalue. From Figure 3, we can also find the eigenvalue-wise parallel method has optimal convergence order even for the first 1000 eigenpairs, which shows the efficiency of Algorithm 2 and validity of Theorem 3.3 and Corollary 3.1. In order to check the orthogonality of approximate eigenfunctions by Algorithm 2, we investigate the inner products of eigenfunctions corresponding to different eigenvalues. We compute inner products for the first 100 approximate eigenfunctions on the finest level of mesh by Algorithm 2. Figure 4 shows the biggest values of inner product of eigenfunctions according to different eigenvalue along with the growth of correction steps on the finest level of mesh. The results in Figure 4 shows that Algorithm 2 can keep the orthogonality when the algebraic accuracy is small enough. These results show the validity of Remark 3.2.

A more general eigenvalue problem
In this example, we consider the following second order elliptic eigenvalue problem: Find (λ, u) ∈ R × V such that ∇u 0 = 1 and    −∇ · (A∇u) + ϕu = λu in Ω, u = 0 on ∂Ω, u 0,Ω = 1, and Ω = (0, 1) × (0, 1) × (0, 1). In order to check the parallel property of Algorithm 2, we compute the first 200 eigenpairs of (5.2). Here, we choose H = 1/16, β = 2, ̟ = ̟ n = 1 and n = 5. In the first step of one correction step defined by Algorithm 1, 1 multigrid step with 2 Conjugate-Gradient (CG) steps for pre-and post-smoothing is adopted to solve the linear problem (3.6). Since the exact solutions are not known, the adequate accurate approximations are chosen as the exact solutions for our numerical test. Figure 5 shows the corresponding error estimates of |λ i − λ i,hn | for i = 1, · · · , 200 and the CPU time for each eigenpair, respectively. From Figure 5, we can find that Algorithm 2 has the optimal error estimate, and computing different eigenpair needs similar computational work.
In this example, we also test the algebraic error |λ i,hn − λ i,hn | between the numerical approximations by Algorithm 2 and the exact finite element approximations for the first 20 eigenvalues along with the growth of the number of correction steps. The corresponding results are presented in Figure 6, which shows that the algebraic accuracy improves with the growth of ̟ n .
Then, we compute the first 1000 eigenpairs. The corresponding error estimates for the approximate eigenvalues and CPU time for each eigenpair are shown in Figure 7. From Figure 7, we can also find that the eigenvalue-wise parallel method can arrive the theoretical convergence order for the first 1000 eigenpairs, which shows the efficiency of Algorithm 2 and the validity of Theorem 3.3 and Corollary 3.1.
The orthogonality of approximate eigenfunctions by Algorithm 2 is also tested. In Figure 8, we show the biggest values of inner product for the first 100 approximate eigenfunctions according to   different eigenvalues on the finest level of mesh. Form Figure 8, it can be found that Algorithm 2 can keep the orthogonality when the algebraic accuracy is small enough.

Adaptive finite element method
In this example, we consider the following eigenvalue problem (see [12]): Find (λ, u) ∈ R × V such that ∇u 0 = 1 and where Ω = R 3 and |x| = x 2 1 + x 2 2 + x 2 3 . The eigenvalues of (5.3) are λ i,j,k = i + j + k + 3 2 , where i, j, k denote the integral numbers and i, j, k ≥ 0. Since the eigenfunctions are exponential decay, we set Ω = (−4, 4) 3 and the boundary condition u = 0 on ∂Ω in our computation for simplicity. Since the exact eigenfunction with singularities is expected, the adaptive refinement is adopted to couple with Algorithm 2.
In order to check the parallel property of Algorithm 2, we compute the first 200 eigenpairs of (5.3). In this example, we choose H = 1/4, ̟ = ̟ n = 1. In the first step of one correction step defined by Algorithm 1, 1 multigrid step with 2 Conjugate-Gradient (CG) steps for pre-and post-smoothing is adopted to solve the linear problem (3.6). Figure 9 shows the corresponding error estimates of |λ i − λ i,hn | for i = 1, · · · , 200 and the CPU time for each eigenpair. From Figure  9, we can find that Algorithm 2 has the optimal error estimates and different eigenvalue has similar computational work even on the adaptively refined meshes. These results show that Algorithm 2 can also be coupled with the adaptive refinement technique.

Adaptive finite element method for Hydrogen atom
In order to show the potential for electrical structure simulation, in the last example, we consider the following model for Hydrogen atom: Find (λ, u) ∈ R × V such that ∇u 0 = 1 and where Ω = R 3 . The eigenvalues of (5.4) are λ h = − 1 2n 2 with multiplicity n 2 for any positive integer n. Along with the growths of n, it is easy to find that the spectral gap becomes small and the multiplicity large which improve the difficulty of solving the eigenvalue problem. The aim of this example is to show the method in this paper can also compute the cluster eigenvalues and their eigenfunctions. Since the eigenfunction is exponential decay, we set Ω = (−4, 4) 3 and the boundary condition u = 0 on ∂Ω in our computation. Here the adaptive refinement is also adopted to couple with Algorithm 2.
In order to check the parallel property of Algorithm 2, we compute the first 200 eigenpairs of (5.4). In this example, we choose H = 1/4, ̟ = ̟ n = 1. In the first step of one correction step defined by Algorithm 1, 1 multigrid step with 2 Conjugate-Gradient (CG) steps for pre-and post-smoothing is adopted to solve the linear problem (3.6). Figure 10 shows the corresponding error estimates of |λ i − λ i,hn | for i = 1, · · · , 200 and the CPU time for each eigenpair. From Figure  10, we can also find that Algorithm 2 has the optimal error estimate and different eigenvalues has similar computational work. These results also shows that Algorithm 2 can be coupled with the adaptive refinement technique.

Concluding remarks
In this paper, we propose an eigenwise parallel augmented subspace scheme for eigenvalue problems by using coarse space from the multigrid method. In this numerical method, solving the eigenvalue problem in the finest space is decomposed into solving the standard linear boundary value problems and very low dimensional eigenvalue problems. Furthermore, for different eigenvalues, the corresponding boundary value problem and low dimensional eigenvalue problem can be solved in the parallel way since they are independent of each other and there exists no data exchanging. This property means that we do not need to do the orthogonalization in the highest dimensional spaces and the the efficiency and scalability can be improved obviously.