RICHARDSON EXTRAPOLATION OF ITERATED DISCRETE COLLOCATION METHOD FOR EIGENVALUE PROBLEM OF A TWO DIMENSIONAL COMPACT INTEGRAL OPERATOR

Abstract

In this paper, we consider approximation of eigenelements of a two dimensional compact integral operator with a smooth kernel by discrete collocation and iterated discrete collocation methods. By choosing numerical quadrature appropriately, we obtain convergence rates for gap between the spectral subspaces, and also we obtain superconvergence rates for eigenvalues and iterated eigenvectors. We then apply Richardson extrapolation to obtain further improved error bounds for the eigenvalues. Numerical examples are presented to illustrate theoretical estimates.

1. Introduction

Consider the following integral operator K defined on by

where kernel K(., ., ., .) ∈ C(D) × C(D), D = [a, b] × [c, d] ⊂ ℝ2 is a given rectangular domain. Then K is a compact linear operator on .

We are interested in the following eigenvalue problem: find and λ ∈ ℂ − {0} such that

Many practical problems in science and engineering are formulated as eigenvalue problems (2) of compact linear integral operators K (cf.,[3]). For many years, numerical solution of eigenvalue problems have attracted much attention. During the last some years, significant work has been done in the numerical analysis of the one-dimensional eigenvalue problem for the compact integral operator K. The Galerkin, petrove-Galerkin, collocation, Nyström and degenerate kernel methods are the commonly used methods for the approximation of eigenelements of the compact integral operator K. The analysis for the convergence of Galerkin, petrove-Galerkin, collocation, Nystr¨om and degenerate kernel methods for the one dimensional eigenvalue problems are well documented in ([1], [3], [12], [13], [14]). By replacing the various integrals appearing in these methods by numerical quadrature leads to discrete methods. In ([9]) discrete and iterated discrete Galerkin methods and in ([4]) discrete and iterated discrete collocation methods were discussed for the one dimensional eigenvalue problem (2) with a smooth kernel.

In [15], we were interested to solve the eigenvalue problem of a two dimensional compact integral operator with smooth kernel taking the help of discrete Galerkin and iterated discrete Galerkin methods and obtained the error bounds for approximated eigenelements. Further, to improve the convergence rates for the eigenvalues, we derived an asymptotic expansion for the iterated discrete Galerkin operator and then using Richardson extrapolation we improved the convergence rates for the eigenvalues. Meanwhile, to do so, we replace the various integrals arise in L2 inner product, when the projection is an orthogonal projection and the two-dimensional integral operator K by using Gauss quadrature rule, which improves the computational cost to generate the matrix eigenvalue problem. To avoid such, discrete collocation method receive favorable attention due to lower computational cost in generating matrix eigenvalue problem. In fact, in comparison to discrete Galerkin and discrete petrove Galerkin methods, the discrete collocation method requires much less computational effort in evaluation of its entries defined by integrals. This motivates us to do this work.

In section-2, we develop discrete collocation, iterated discrete collocation methods and theoretical frame work for the eigenvalue problem using interpolatory projections. In section-3, we discuss the convergence rates for the approximated eigenfunctions to the exact eigenfunctions. In section-4, we discuss Richardson extrapolation for eigenvalue problem to improve convergence rates. In section-5, we present numerical results, which agree with the theoretical results. Throughout the paper, we assume c as the generic constant.

 

2. Discrete and Iterated discrete collocation methods

Consider the following compact integral operator K defined on by

where the kernel K(., ., ., .) ∈ C(D) × C(D), D = [a, b] × [c, d] ⊂ ℝ2.

We are interested in the eigenvalue problem (2). Assume λ be the eigenvalue of K with algebraic multiplicity m and ascent ℓ. Let Γ ⊂ ρ(K) be a simple closed rectifiable curve such that σ(K) ∩ intΓ = {λ}, 0 ≠ intΓ, where intΓ denotes the interior of Γ. Now we describe the collocation method for the eigenvalue problem (2).

Let be the uniform partitions of finite intervals [a, b] and [c, d], respectively, defined by and for m = 0, 1, 2, . . . ,M − 1 and n = 0, 1, 2, . . . ,N − 1. These partitions define a grid for D, : 0 ≤ m ≤ M − 1, 0 ≤ n ≤ N − 1}. Set and m = 0, 1, . . . ,M − 1 and n = 0, 1, . . . ,N − 1. For any given positive integer p and q, let Pp−1,q−1 denotes the space of polynomials of degree p − 1 in x and q − 1 in y, then for 0 ≤ m ≤ M − 1, 0 ≤ n ≤ N − 1,

is the finite element space of dimension MNpq, which is the tensor product space of univariate spline spaces on [c, d]. The use of superscript (-1) in the notation for the above finite element space is to emphasize that it is not a subspace of C(D).

Let sm,i = xm + τih and tn,j = yn + θjk, be the collocation points on [xm, xm+1], m = 0, 1, . . . ,M −1, and [yn, yn+1], n = 0, 1, . . . ,N −1, respectively, where τ0, τ1, . . . , τp−1 and θ0, θ1, . . . , θq−1 are zeros of Legendre polynomials of degree p and q, respectively on [0, 1].

Let and be the interpolatory projection with respect to the nodes {sm,i} and {tn,j}, respectively, that is, for u ∈ L∞([a, b]),

then there holds([3]), ∥Ph∥∞ ≤ c1 < ∞, ∥Pk∥∞ ≤ c2 < ∞ and for any u ∈ Cp[a, b] and u ∈ Cq[c, d], there holds,

Note that

As a consequence, from (6), it follows that

Let ϕim, ψjn denote the Lagrange polynomials of degree p − 1 and q − 1 on the subintervals [xm, xm+1], m = 0, 1, . . . ,M − 1 and [yn, yn+1], n = 0, 1, . . . ,N − 1 respectively, where, j = 0, 1, . . . q − 1, i = 0, 1, . . . , p − 1. Then it follows that and Then we have

Now the collocation method for solving the eigenvalue problem (2) is defined as follows: find ∥uhk∥ = 1 and λhk ∈ ℂ − {0} such that

for i′ = 0, 1, . . . p−1, m′ = 0, 1, . . .M−1, j′ = 0, 1, . . . q−1, n′ = 0, 1, . . .N −1. Using the equation (9) can be converted to the matrix eigenvalue problem. The iterated eigenvector is defined by

To solve the matrix eigenvalue problem and the iterated eigenvector, we need to evaluate various integrals arising from the integral operator K. In practice, numerical quadrature has to be used to compute these integrals. This leads to discrete methods. To do this, let for ƒ, g ∈ C[0, 1],

be the numerical quadrature with weights and quadrature points ci, i = 0, 1, . . . , k′ − 1, dj, j = 0, 1, . . . , l′ − 1, chosen as Gauss points in [0, 1] which satisfy 0 < c0 < c1 < · · · < ck′−1 < 1 and 0 < d0 < d1, . . . ,< dl′−1 < 1 having degree of precision 2k′ − 1, 2l′ − 1, respectively. Then the composite Gauss quadrature rule for any ƒ ∈ C([a, b]), g ∈ C([c, d]) is given by

where xm,i = xm+cih, i = 0, 1, . . . , k′−1, and yn,j = yn+djk, j = 0, 1 . . . , l′−1, be the quadrature points on the subintervals [xm, xm+1], m = 0, 1, . . . ,M − 1 of [a, b] and [yn, yn+1], n = 0, 1, . . . ,N −1 of [c, d], respectively. Now using (12) and (13), we define the composite quadrature rule for g ∈ C(D) by

Let be the Nyström operator defined for

Now replacing the integral operator K by the Nyström operator (14), the matrix eigenvalue problem leads to the discrete collocation method,

By solving this discrete matrix eigenvalue problem (15), we find the eigenvalue and β = [βij, i = 0, 1 . . . , p−1, j = 0, 1, . . . , q −1]. Then the discrete collocation eigenvector is defined by, The discrete matrix eigenvalue problem (15) can be written in operator form as

Next we define the iterated discrete collocation eigenvector by Clearly we see that and

This is the iterated discrete collocation method.

Next we discuss the convergence of approximated eigenvalues and eigenvectors to the exact eigenvalues and eigenvectors of the integral operator K. To do this, first we set the following notations: Set K(s, t, x, y) = Ks,t(x, y). For K(., ., ., .) ∈ C(i,j)(D) × C(i′,j′)(D), denote

For any α, β, γ, δ ∈ ℕ, we set

Then we have

and

In the following theorem we give the error bounds for the Nystr¨om operator.

Theorem 2.1 ([15]). Let K be an integral operator with a kernel K(., ., ., .) ∈ C(2k′,2l′)(D) × C(2k′,2l′)(D) and Khk be the Nyström operator defined by (14), then for any u ∈ C(2k′,2l′)(D), the following holds

where c is independent of h and k.

Definition 2.2 ([1]). Let be a Banach space and, T and Tn are bounded linear operators from into . Then {Tn} is said to be ν-convergent to T, if

We quote the following lemma which is useful in proving the existence of eigenvalue and eigenvectors in discrete and iterated discrete collocation methods.

Lemma 2.3 ([2]). Let be a Banach space and S ⊂ is a relatively compact set. Assume that T and Tn are bounded linear operators from into satisfying ∥Tn∥ ≤ c for all n ∈ ℕ, and for each x ∈ S, ∥Tnx−Tx∥ → 0 as n → ∞, where c is a constant independent of n. Then ∥Tnx−Tx∥ → 0 uniformly for all x ∈ S.

Theorem 2.4. PhPkKhk and KhkPhPk are ν-convergent to K.

Proof. Since Khk, Ph and Pk are uniformly bounded, it follows that ∥PhPkKhk∥∞ ≤ ∥Ph∥∞∥Pk∥∞∥Khk∥∞ < ∞ and ∥KhkPhPk∥∞ ≤ ∥Khk∥∞∥Ph∥∞∥Pk∥∞ < ∞. Now using (8) and Theorem 2.1, we see that

This shows that PhPkKhk point wise converges to K.

Let be a closed unit ball in . Since K is a compact operator, the set S = {Kx : x ∈ B} is a relatively compact set in X. By Lemma 2.3, we have

Since PhPk is bounded and Khk compact, S′ = {PhPkKhkx : x ∈ B} is a relatively compact set. Thus

as h, k → 0. Combining all these results leads to the first result that PhPkKhkx is ν-convergent to K. The proof of KhkPhPk is ν-convergent to K follows by similar steps as in above.

Since PhPkKhk and KhkPhPk are ν-convergent to K, the spectrum of both PhPkKhk and KhkPhPk inside Γ consists of m eigenvalues say counted accordingly to their algebraic multiplicities inside Γ with ascent ℓ (cf., [3,14]). Let

denote their arithmetic mean and we approximate λ by Let

be the spectral projections of K associated with their corresponding spectra inside Γ. Similarly, be the spectral projections of PhPkKhk and KhkPhPk, respectively. Let be the ranges of the spectral projections respectively.

To discuss the closeness of eigenfunctions of the integral operator K and those of the approximate operators, we recall (cf., [3]) the concept of gap between the spectral subspaces. For nonzero closed subspaces let let

then

is known as the gap between

We quote the following three Lemmas, which give the error bounds for the eigenelements.

Theorem 2.5 ([1], [13]). Let PhPkKhk be ν-convergent to K. Then for sufficiently large M,N, there exists a constant c independent of M,N, we have

Theorem 2.6 ([13], [16]). Let KhkPhPk be ν-convergent to K. Then for sufficiently large M,N, there exists a constant c independent of M,N, we have

Theorem 2.7 ([1], [13]). If KhkPhPk is ν-convergent to K then for sufficiently large M,N, there exists a constant c independent of M,N such that for j = 1, 2...., m,

 

3. Convergence Rates

In this section we discuss the convergence rates for the approximated eigenvalues and eigenvectors to the exact eigenvalues and exact eigenvectors of the integral operator K. To do this, first we prove the following Lemma.

Lemma 3.1. Let Khk be the Nyström operator defined by (14) with a kernel K(., ., ., .) ∈ C(2k′,2l′)(D)×C(2k′,2l′)(D), k′ ≥ p, l′ ≥ q. Then for u ∈ C(2p,2q)(D), the following hold

Proof. Using the estimates (6), we have

and

Combining these estimates with the identity (7), proof of (22) follows.

To prove the estimate (23), let us denote Since H(t) and are orthogonal polynomials of degree p and q, respectively, and the numerical quadratures defined by (10) and (11) have degree of precision 2k′ and 2l′, respectively, it follows that, for k′ ≥ p and l′ ≥ q,

Since Ph is the interpolatory projection interpolating at u(s, t) in the first variable s at the points sm,0, sm,1, . . . , sm,p−1 in the subintervals [xm, xm+1], m = 0, 1, . . . ,M − 1, we have for s ∈ [xm, xm+1], t ∈ [c, d],

where δ(p,0)u(s, t) = [sm,0, sm,1, . . . , sm,p−1, s; t]u be the Newton divided difference of u in first variable. Similarly, since Pk is the interpolatory projection interpolating at u(s, t) in the second variable t at the points tn,0, tn,1, . . . , tn,q−1 in the subintervals [yn, yn+1], n = 0, 1, . . . ,N − 1, we have for t ∈ [yn, yn+1], s ∈ [a, b],

where δ(0,q)u(s, t) = [s; tn,0, tn,1, . . . , tn, q−1, t]u be the Newton divided difference of u in second variable. Now using the identity (7), we have

For the first term in the above, for any (s, t) ∈ D, using (26), we obtain

where gi(xm,i, yn,j) = Ks,t(xm,i, yn,j)δ(p,0)u(xm,i, yn,j). The Taylor’s expansion of gi(xm,i, yn,j) = gi(xm + cih, yn,j) at the point xm is given by

where ξm ∈ [xm, xm+1]. Using (30) in the estimate (29), we obtain

Using the estimate (24) in (31), it follows that

On the similar mechanism, for the second term in (28), we can prove that

Let δ(p,q)u(s, t) = [sm,0, sm,1, . . . , sm,p−1, s; tn,0, tn,1, . . . , tn,q−1, t]u be p and qth Newton divided difference of u in first and second variables, respectively. Then we have

Using this in the third term of (28), we have

where gi,j(xm,i, yn,j) = Ks,t(xm,i, yn,j)δ(p,q)u(xm,i, yn,j ). The Taylor’s series expansion for gi,j(xm,i, yn,j) = gi,j(xm + cih, yn + djk) at the point xm and yn is given by

where r1 = p+q. Using (34) in the above equation and by adjustment with the the estimates (24) and (25), we obtain

where is a constant independent of h and k. Combining the estimates (28), (32) and (35), the result (23) follows. This completes the proof.

Theorem 3.2. Let K be a compact integral operator with a kernel K(., ., ., .) ∈ C(2k′,2l′)(D) × C(2k′,2l′)(D), k′ ≥ p, l′ ≥ q and Khk be the Nyström operator defined by (14). Then the following hold

Proof. Replacing u by Ku in (20), and using the estimate (19), we obtain

where c is a constant independent of h and k. This completes the proof of (36). Now replacing u by Ku in (22) and using the estimate (19), we see that

this proves the estimate (37). Again replacing u by Ku in (23), then we obtain

this proves the estimate (38).

Theorem 3.3. Assume that all the conditions of theorem 3.2 hold. Then the following hold

Proof. Since

the proof follows from the above Theorem 3.2.

In the following Theorem we give the superconvergence results for the eigenvalues and eigenvectors

Theorem 3.4. Suppose K is a compact integral operator with a kernel function K(., ., ., .) ∈ C(2k′,2l′)(D) × C(2k′,2l′)(D), k′ ≥ p, l′ ≥ q, and suppose that λ be the eigenvalue of K with algebraic multiplicity m and ascent ℓ. Let {KhkPhPk} and {PhPkKhk} be a sequence of bounded operators on ,which converges to K in ν−convergence. Then

In particular, for any we have

For j = 1, 2, . . . , m,

Proof. The proof follows directly using the Theorems 2.5, 2.6, 2.7 and 3.3.

Remark: From Theorem 3.4, we observe that discrete collocation eigenvectors converges with the order of convergence O(max{hmin{p,2k′}, kmin{q,2l′}}) where as iterated discrete collocation eigenvectors and eigenvalues converges with the order of convergence O(max{h{min{2p,2k′}, k{min{2q,2l′}}). This shows that iterated discrete eigenvectors gives superconvergence results over the discrete collocation eigenvectors. By choosing the degree of precisions of the numerical quadrature rules sufficiently large, i.e., 2k′ ≥ 2p and 2l′ ≥ 2q on [a, b] and [c, d], respectively, we obtain the superconvergence results for the eigenvalues and eigenvectors in the discrete collocation and iterated discrete collocation methods.

 

4. Richardson Extrapolation

In this section, we derive an asymptotic error expansions (cf., [10], [11]) for the iterated discrete collocation operator KhkPhPk and an asymptotic error expansion of arithmetic mean of approximate eigenvalues. We then apply Richardson extrapolation to obtain improved error bounds for the eigenvalues.

Lemma 4.1. (Euler-MacLaurin summation formulae)([7]).

Let f(x, y) ∈ Cr+1(D), 0 ≤ τ ≤ 1, 0 ≤ θ ≤ 1. Then

where Bi(t) are Bernoulli polynomials of degree i.

Theorem 4.2 ([15]). Let K be a compact integral operator with a kernel K(., ., ., .) ∈ C(D) × Cr+1(D) and Khk be the Nyström operator defined by (14), then there holds

where D2i, ε2j, and F2i,2j are bounded linear operators independent of h and k.

Lemma 4.3 ([7]). Assume that u(x, y) ∈ Cr+1(D). Let Ph and Pk be the interpolatory projections defined by (4) and (5), respectively. Then for any (x, y) ∈ Imn, the following holds

where

Proposition 4.4. Assume that the kernel K(., ., ., .) ∈ Cr+1(D) × Cr+1(D) and u ∈ Cr+1(D), then the following holds

where R2i,2j, S2i,2j and T2i,2j are bounded linear operators on Cr+1(D).

Proof. Using the definition of Khk defined by (14) and the Lemma 4.3, we have

Now we consider I1,

Using Euler-Maclaurin summation formula in I1, we obtain

where and S(Bb) is the numerical quadrature of Ba defined as in (11). Since τ0, τ1, . . . , τp−1 are symmetric points in the interval [0, 1], i.e., τj = 1 − τp−1−j for 0 ≤ j ≤ p − 1, we have

and

Using these estimates we obtain Φμ(1−ci) = (−1)μΦμ(ci), i = 0, 1, 2, . . . , k′−1. Also note that Ba(ci) = (−1)aBa(1 − ci), ci = 1 − ck′−1−i and wk′−i−1 = wi, for 0 ≤ i ≤ k′ − 1. Hence we have

From this, it follows that Aμ,a = 0, when μ + a = odd. Since the quadrature rule (10) is exact for polynomial of degree less than 2k′ and is orthogonal to all polynomial of degree less than p and is a polynomial of degree μ − p, then for μ + a < 2p, k′ ≥ p, we have

and Thus we obtain

Since dj, j = 0, 1, . . . , l′ − 1 are the Gauss points in the interval (0, 1) and the quadrature rule (3) is an symmetric quadrature rule, we have dl′−j−1 = 1 − dj and for j = 0, 1, . . . , l′ − 1. Noting that the Bernoulli polynomials have the property that Bb(1 − d) = (−1)bBb(d), we have

From this, it follows that S(Bb) is zero when b is odd. Also we have S(B0) = 1. Now since the degree of precision of this quadrature rule (10) is 2l′−1, it follows that for 1 ≤ b ≤ 2l′ − 1,

Thus we have

Combining the estimates (46) and (45) with (44) we obtain

where,

Similarly we can prove that

where,

Similarly for I3, we can prove that

where,

Now combining the estimates for I1, I2 and I3 with (43), we obtain the following asymptotic expansion This completes the proof.

Theorem 4.5. Let K be a compact integral operator with the kernel K(., ., ., .) ∈ Cr+1(D) × Cr+1(D) with p = k′ and q = l′. Then the following holds

where U2i, V2j and W2i,2j are bounded linear operators on Cr+1(D).

Proof. Combining Theorems 4.4 and 4.2, we have

where U2i = R2i,0−D2i, V2j = S0,2j −ε2j, and W2i,2j = R2i,2j + S2i,2j − T2i,2j − F2i,2j are bounded linear operators on Cr+1(D). This completes the proof.

In rest of the paper, we choose the domain D = [a, b]×[a, b] and the partition Also we choose the numerical quadratures (10) and (11) to be same, i.e., k′ = l′. Then we have the following corollary.

Corollary 4.6. Let Khk be the Nystr¨om operator defined by (14) with p = q = k′ = l′. Assume that the kernel K(., ., ., .) ∈ C2p+2(D) × C2p+2(D). Then the following holds

where C2p = U2p + V2p is a bounded linear operators on C2p+2(D).

Richardson Extrapolation For Eigenvalue problem: By the similar way as followed in [15], we obtain the following theorem which gives an asymptotic error expansion of the arithmetic mean of eigenvalues by iterated discrete collocation method.

Theorem 4.7 ([15]). Let λ be the eigenvalue of K with algebraic multiplicity m and be the arithmetic mean of the eigenvalues Then the following holds

where Q2p = C2pPS − KU2p is a bounded linear operator independent of h.

According to the asymptotic expansion (47), the Richardson extrapolation for eigenvalue approximation should be the following. We first divide each subinterval with respect to the partitions of into two halves which makes up a new partitions denoted by

Here and D = [a, b] × [a, b]. We then have following asymptotic expansion for eigenvalue approximation with respect to this new partitions,

From the asymptotic expansions (47) and (48), the Richardson extrapolation for the eigenvalue approximation is defined by

In the following Theorem we give the superconvergence rates for the eigenvalue approximation using Richardson extrapolation.

Theorem 4.8. Assume that conditions of Theorem 4.7 hold and the Richardson extrapolation is defined by (49). Then the following error estimate holds

 

5. Numerical Example

Consider the eigenvalue problem (2), for the integral operator K (3) for various smooth kernels K(s, t, x, y).

Let be the space of piecewise constant functions (p=q=1) on [0, 1]×[0, 1] with respect to the initial uniform partitions with We choose numerical quadrature as the one-point composite Gaussian quadrature formula which is exact for all polynomials of degree less than 2, that is k = l = 1.

The quadrature points and weights are given by

and respectively.

For different kernels and for different values of M, we compute the discrete collocation eigen vector iterated eigen vector and eigenvalue and approximated eigenvalue in Richardson extrapolation Denote

where is the step length. For M = 2, 4, 8, 16, we compute α, β, γ and δ which are listed in the following Table.

Since k′ = l′ = 1, we get the theoretical convergence of the order of eigenvector is 1, the orders of the iterated eigenvector and eigenvalues are 2, and the order of eigenvalue in Richardson extrapolation is 4. In the following Table 1 and Table 2, the numerical results agrees with the theoretical results.

Table 1Eigenvector error bounds

Table 2Eigenvalue error bounds

Example. K(s, t, x, y) = s sin(t) + xey, [a, b] = [0, 1], [c, d] = [0, 1].