Laguerre spectral method

# A fully diagonalized spectral method using generalized Laguerre functions on the half line

Fu-jun Liu Zhong-qing Wang  and  Hui-yuan Li
###### Abstract.

A fully diagonalized spectral method using generalized Laguerre functions is proposed and analyzed for solving elliptic equations on the half line. We first define the generalized Laguerre functions which are complete and mutually orthogonal with respect to an equivalent Sobolev inner product. Then the Fourier-like Sobolev orthogonal basis functions are constructed for the diagonalized Laguerre spectral method of elliptic equations. Besides, a unified orthogonal Laguerre projection is established for various elliptic equations. On the basis of this orthogonal Laguerre projection, we obtain optimal error estimates of the fully diagonalized Laguerre spectral method for both Dirichlet and Robin boundary value problems. Finally, numerical experiments, which are in agreement with the theoretical analysis, demonstrate the effectiveness and the spectral accuracy of our diagonalized method.

###### Key words and phrases:
Spectral method, Sobolev orthogonal Laguerre functions, elliptic boundary value problems, error estimates.
###### 2000 Mathematics Subject Classification:
76M22, 33C45, 35J25, 65L70
Department of Mathematics, Shanghai Normal University, Shanghai, 200234, China; School of Science, Henan Institute of Engineering, Zhengzhou, 451191, China. Email: liufujun1981@126.com
School of Science, University of Shanghai for Science and Technology, Shanghai, 200093, China. Email: zqwang@usst.edu.cn
State Key Laboratory of Computer Science/Laboratory of Parallel Computing, Institute of Software, Chinese Academy of Sciences, Beijing 100190, China. Email: huiyuan@iscas.ac.cn
The first author was supported by Science and Technology Research Program of Education Department of Henan Province (No. 13A110005); The second author was supported in part by National Natural Science Foundation of China (No. 11571238) and the Research Fund for Doctoral Program of Higher Education of China (No. 20133127110006); The third author was supported by National Natural Science Foundation of China (Nos. 91130014, 11471312 and 91430216).

## 1. Introduction

Spectral methods for solving partial differential equations on unbounded domains have gained a rapid development during the last few decades. An abundance of literature on this research topic has emerged, and their underlying approximation approaches can be essentially classified into three catalogues [4, 27]:

• truncate an unbounded domain to a bounded one and solve the problem on the bounded domain subject to artificial or transparent boundary conditions [22, 26];

• map the original problem on an unbounded domain to one on a bounded domain and use classic spectral methods to solve the new problem [9]; or equivalently, approximate the original problem by some non-classical functions mapped from the classic orthogonal polynomials/functions on a bounded domain [2, 3, 7, 11, 12, 27, 31, 34];

• directly approximate the original problem by genuine orthogonal functions such as Laguerre polynomials or functions on the unbounded domain [6, 10, 13, 14, 15, 16, 17, 18, 19, 20, 24, 30, 32, 33, 35].

The third approach is of particular interest to researchers, and has won an increasing popularity in a broad class of applications, owing to its essential advantages over other two approaches. These direct approximation schemes constitute an initial step towards the efficient spectral methods, which admit fast and stable algorithms for their efficient implementations.

As we know, the Fourier spectral method makes use of the eigenfunctions of the Laplace operator which are orthogonal to each other with respect to the Sobolev inner product involving derivatives, thus the corresponding algebraic system is diagonal [4, 5, 25]. This fact together with the availability of the fast Fourier transform (FFT) makes the Fourier spectral method be an ideal approximation approach for differential equations with periodic boundary conditions. Although the utilization of the genuine orthogonal polynomials/functions in this direct approach usually leads to a highly sparse (e.g., tri-diagonal, penta-diagonal) and well-conditioned algebraic system, however, in many cases, people still want to get a set of Fourier-like basis functions for a fully diagonalized algebraic system [28].

The main purpose of this paper is to construct the Fourier-like Sobolev orthogonal basis functions [8, 21] for elliptic boundary value problems on the half line . For this purpose, we shall first extend the definition of Laguerre polynomials and Laguerre functions for to allow being any real number. The resulting generalized Laguerre functions are proven to be the eigenfunctions of certain high order Sturm-Liouville differential operators (see Lemma 2.6 of this paper). Moreover, they are complete and mutually orthogonal in for any nonnegative integer with respect to an equivalent Sobolev inner product (see (2.23) of this paper).

Since the problem is dependent on the inner product originated from the coercive bilinear form of the elliptic equation, it does not necessarily coincide with the equivalent Sobolev inner product, further efforts should be paid to obtain the Fourier-like basis functions for a fully diagonalized spectral approximation, in spite of the Sobolev orthogonality of . Starting with , stable and efficient algorithms are then proposed to construct the Fourier-like basis functions for the non-homogeneous Dirichlet and Robin boundary value problems of the second order elliptic equations. In the sequel, both the exact solution and the approximate solution can be represented as infinite and truncated Fourier series in , respectively. Although the fully diagonalized spectral methods are studied for second order equations, they can be readily generalized to solve -th order equations by starting with .

An ideal spectral approximation to differential equations may guarantee an optimal error estimate in its convergence analysis. To match this requirement, various orthogonal projections involving different orders of derivatives and boundary conditions have been designed and studied case by case, which frequently make the numerical analysis in spectral method a tedious task. Moreover, the traditional routine to measure the approximation error is first to establish the norm defined by a second-order self-adjoint differential operator, and then estimate the upper bound of the approximation error with the induced norms. However, this practical approach usually fails to characterize the function space in which the orthogonal projection has an optimal error estimate.

To conquer these difficulties, we need a unified definition of the orthogonal spectral projections with a systematic numerical analysis. Fortunately, the Sobolev orthogonality of the generalized Laguerre functions with a negative integer enables us to define the unified orthogonal projection from to the finite approximation space for all nonnegative integer , ignoring the specific value of . More importantly, such an orthogonal projection interpolates the endpoint function values up to the -th derivative, i.e, for any and . This endpoint interpolation property ensures , thus makes applicable to both the Dirichlet and Robin boundary value problems, and available to multi-domain spectral methods. Besides, owing to the clarity of the orthogonality structure of the generalised Laguerre functions, one can not only derive an optimal order of the convergence for the approximated function, but also get a generic characterization of the function space where the orthogonal projection has an optimal error estimate.

Therefore, the second purpose of this paper is to establish such a unified orthogonal Laguerre projection, and apply it to the convergence analysis on the fully diagonalized Laguerre spectral method for both the Dirichlet and Robin boundary value problems of second order elliptic equations.

The remainder of the paper is organized as follows. In Section 2, we first make conventions on the frequently used notations, and then introduce generalized Laguerre polynomials and functions with arbitrary index . The fully diagonalized Laguerre spectral methods and the implementation of algorithms are proposed in Section 3 for the Dirichlet and Robin boundary value problems of second order elliptic equations. Section 4 is then devoted to the convergence analysis of the unified orthogonal projection together with our Laguerre spectral methods. Finally, numerical results are presented in Section 5 to demonstrate the effectiveness and accuracy of the proposed diagonalized Laguerre spectral methods, which are in agreement with our theoretical predictions.

## 2. Generalized Laguerre polynomials and functions

### 2.1. Notations and preliminaries

Let and be a weight function which is not necessary in . We define

 L2ϖ(Λ)={v | v is measurable on Λ and~{}∥v∥ϖ<∞},

with the following inner product and norm,

 (u,v)ϖ=∫Λu(x)v(x)ϖ(x)dx,∥v∥ϖ=(v,v)12ϖ,∀u,v∈L2ϖ(Λ).

For simplicity, we denote and For any integer , we define

 Hmϖ(Λ)={v | ∂kxv∈L2ϖ(Λ), 0≤k≤m},

with the following semi-norm and norm,

For any real we define the space and its norm by function space interpolation as in [1]. In cases where no confusion arises, may be dropped from the notations whenever Specifically, we shall use the weight functions and in the subsequent sections.

We denote by the collection of real numbers, by and the collections of nonnegative and negative integers, respectively. Further, we let be the space of polynomials of degree .

Let . We also define the characteristic functions for ,

 χn(α)={−α,α+n∈Z−,0,α+n∈(−1,+∞).

For short we write .

### 2.2. Generalized Laguerre polynomials

It is well known that, for the classical Laguerre polynomials admit an explicit representation (see [29]):

 Lαk(x)=k∑ν=0(α+ν+1)k−ν(k−ν)!ν!(−x)ν,x∈Λ,k≥0, (2.1)

where we use the Pochhammer symbol for any and .

The classical Laguerre polynomials can be extended to cases with any and the same representation as (2.1), which are referred to as the generalized Laguerre polynomials (cf. [23]). Obviously, the generalized Laguerre polynomials constitute a complete basis for the linear space of real polynomials as well, since for all .

The generalized Laguerre polynomials fulfill the following recurrence relations.

###### Lemma 2.1.

For any , it holds

 {Lα0(x)=1,Lα1(x)=−x+α+1,(k+1)Lαk+1(x)=(2k+α+1−x)Lαk(x)−(k+α)Lαk−1(x),k≥1. (2.2)
###### Proof.

The recurrence relation (2.2) for can be derived from those of the classic Laguerre polynomials for by the continuation method. Here, we also give a concrete proof by the representation (2.1). Using the expression (2.1), we obtain that for integer

 (2k+α+1)Lαk(x)−(k+α)Lαk−1(x)−(k+1)Lαk+1(x)=(2k+α+1)k∑ν=0(−1)νΓ(k+α+1)ν!(k−ν)!Γ(α+ν+1)xν−(k+α)k−1∑ν=0(−1)νΓ(k+α)ν!(k−ν−1)!Γ(α+ν+1)xν−(k+1)k+1∑ν=0(−1)νΓ(k+α+2)ν!(k−ν+1)!Γ(α+ν+1)xν.

Then a direct computation shows that

 (2k+α +1)Lαk(x)−(k+α)Lαk−1(x)−(k+1)Lαk+1(x) = k+1∑ν=1(−1)ν−1Γ(k+α+1)(ν−1)!(k−ν+1)!Γ(α+ν)xν=xLαk(x).

The desired result is now derived. ∎

###### Lemma 2.2.

For any and , it holds

 Lαk(x)=Lα+1k(x)−Lα+1k−1(x), (2.3) ∂xLαk(x)=−Lα+1k−1(x), (2.4) x∂xLαk(x)=kLαk(x)−(k+α)Lαk−1(x), (2.5) Lαk(x)=∂xLαk(x)−∂xLαk+1(x), (2.6)

where for any

###### Proof.

The recurrence relations (2.3)-(2.5) can be obtained readily by using similar arguments as in Lemma 2.1. Moreover, by (2.3) and (2.4), it is easy to derive (2.6). ∎

###### Lemma 2.3.

For any , the generalized Laguerre polynomials satisfy the Sturm-Liouville equation

 x−αex∂x(xα+1e−x∂xLαk(x))+λkLαk(x)=0,k≥0, (2.7)

or equivalently,

 x∂2xLαk(x)+(α+1−x)∂xLαk(x)+λkLαk(x)=0,k≥0, (2.8)

with the corresponding eigenvalue .

###### Proof.

Lemma 2.3 can be proved by the continuation method from the Sturm-Liouville equation of the classic Laguerre polynomials for . Also one can give a proof by using the representation (2.1). We omit the details here. ∎

We are interested in those generalized Laguerre polynomials with an integer index .

###### Lemma 2.4.

For any , we have

 Lαk(x)=(−x)−α(k+α)!k!L−αk+α(x),k≥χ(α). (2.9)

And for any , the following orthogonality relation holds:

 ∫ΛLαk(x)Lαm(x)xαe−xdx=γαkδk,m,γαk=Γ(k+α+1)k!,k,m≥χ(α), (2.10)

where is the Kronecker symbol.

###### Proof.

The identity (2.9) comes directly from [29]. The orthogonality relation (2.10) is known for classic Laguerre polynomials with ; while for , (2.10) can be obtained immediately from (2.9) together with (2.10) for . ∎

We now conclude this subsection with some generalized Laguerre polynomials for .

 k=0k=1k=2k=3…k≥χ(α)α=−11−x12x(x−2)−16x(x2−6x+6)…−1kxL1k−1(x)α=−21−x−112x2−16x2(x−3)…1k(k−1)x2L2k−2(x)α=−31−x−212x2+x+1−16x3…−1k(k−1)(k−2)x3L3k−3(x)…………………

### 2.3. Generalized Laguerre functions

In this subsection, we shall introduce the generalized Laguerre functions with arbitrary parameters and and present some properties.

The generalized Laguerre functions are defined by

 lα,βk(x)=e−12βxLαk(βx),∀α∈R,β>0, (2.11)

and the multiplication of and the leading term of is simply referred to as the leading term of .

According to (2.9), for any , we have

 lα,βk(x)=(−βx)−α(k+α)!k!l−α,βk+α(x),k≥χ(α), (2.12)

which means that is a zero of with the multiplicity , i.e.,

 ∂νxlα,βk(0)=0,k≥χ(α),ν=0,1,…,−α−1. (2.13)

Due to (2.2)-(2.6), the generalized Laguerre functions satisfy the following recurrence relations:

###### Lemma 2.5.

For any , it holds that

 βxlα,βk(x)=−(k+1)lα,βk+1(x)+(2k+α+1)lα,βk(x)−(k+α)lα,βk−1(x), (2.14) lα,βk(x)=lα+1,βk(x)−lα+1,βk−1(x), (2.15) ∂xlα,βk(x)=−βlα+1,βk−1(x)−β2lα,βk(x)=−β2[lα+1,βk(x)+lα+1,βk−1(x)], (2.16) x∂xlα,βk(x)=k+12lα,βk+1(x)−α+12lα,βk(x)−k+α2lα,βk−1(x), (2.17) ∂xlα,βk(x)−∂xlα,βk+1(x)=12β(lα,βk(x)+lα,βk+1(x)). (2.18)

Hereafter, we use the convention that whenever .

The generalized Laguerre functions are eigenfunctions of certain singular Sturm-Liouville differential operators.

###### Lemma 2.6.

For any , it holds that

 n∑ν=0(−1)ν(nν)β2n−2ν22n−2νx−α∂νx(xα+n∂νxlα,βk)=βn2nλαk,nlα,βk,k≥0, (2.19)

where satisfies the following recurrence relation,

 λαk,0=1,λαk,n=(k+α+1)λα+1k,n−1+kλα+1k−1,n−1,n≥1,k≥0. (2.20)
###### Proof.

We prove (2.19) and (2.20) by induction. It is obvious that (2.19) holds for . Moreover, by virtue of (2.7) and (2.11) we have

 −x−α∂x(xα+1∂xlα,βk(x))+β24xlα,βk(x)=β2(2k+α+1)lα,βk(x)=β2λαk,1lα,βk(x),k≥0,

which gives (2.19) and (2.20) for .

We now assume that (2.19) and (2.20) hold for an integer . Then by the recursive formula of binomial coefficients together with (2.15) and (2.16),

 I:= n+1∑ν=0(−1)ν(n+1ν)β2n+2−2ν22n+2−2νx−α∂νx(xn+1+α∂νxlα,βk) = n+1∑ν=0(−1)ν[(nν)+(nν−1)]β2n+2−2ν22n+2−2νx−α∂νx(xn+1+α∂νxlα,βk) = β24xn∑ν=0(−1)ν(nν)β2n−2ν22n−2νx−(α+1)∂νx[xn+(α+1)∂νx(lα+1,βk−lα+1,βk−1)] +β2x−α∂xn+1∑ν=1(−1)ν−1(nν−1)β2n−2(ν−1)22n−2(ν−1)∂ν−1x[xn+(α+1)∂ν−1x(lα+1,βk+lα+1,βk−1)].

Thus by the induction assumption, (2.14), (2.17) and (2.15), we derive that

 I= βn+22n+2x[λα+1k,nlα+1,βk−λα+1k−1,nlα+1,βk−1]+βn+12n+1x−α∂x[xα+1(λα+1k,nlα+1,βk+λα+1k−1,nlα+1,βk−1)] = βn+12n+1λα+1k,n[(α+1+β2x)lα+1,βk+x∂xlα+1,βk]+βn+12n+1λα+1k−1,n[(α+1−β2x)lα+1,βk−1+x∂xlα+1,βk−1] = βn+12n+1λα+1k,n(k+α+1)[lα+1,βk−lα+1,βk−1]+βn+12n+1λα+1k−1,nk[lα+1,βk−lα+1,βk−1] = βn+12n+1[(k+α+1)λα+1k,n+kλα+1k−1,n]lα,βk,

which is exactly (2.19) and (2.20) with in place of . This ends the proof. ∎

For any , and , define the bilinear form on ,

 aα,βn(u,v)=n∑ν=0(nν)β2n−2ν22n−2ν(∂νxu,∂νxv)wα+n. (2.21)

It is obvious that is an inner product on if .

###### Theorem 2.1.

The generalized Laguerre functions for are mutually orthogonal with respect to the weight function ,

 aα,β0(lα,βk,lα,βm)=(lα,βk,lα,βm)wα=β−α−1γαkδk,m,k,m≥χ(α). (2.22)

More generally, for any and ,

 aα,βn(lα,βk,lα,βm)=βn−α−1γαk,nδk,m,k,m≥χn(α), (2.23)

where the positive numbers satisfy the recurrence relation

 γα+nk,0=γα+nk,γαk,n=12[γα+1k,n−1+γα+1k−1,n−1],n≥1, (2.24)

under the convention that whenever .

###### Proof.

The orthogonality (2.22) is an immediate consequence of (2.10). Meanwhile, the recursive formula (2.21) of binomial coefficients together with (2.15) and (2.16) yields

 aα,βn+1 (lα,βk,lα,βm)=n+1∑ν=0[(nν−1)+(nν)]β2n+2−2ν22n+2−2ν(∂νxlα,βk,∂νxlα,βm)wα+n+1 = +β24n∑ν=0(nν)β2n−2ν22n−2ν(∂νx[lα+1,βk−lα+1,βk−1],∂νx[lα+1,βm−lα+1,βm−1])wn+(α+1) (2.25) = β22n∑ν=0(nν)β2n−2ν22n−2ν[(∂νxlα+1,βk,∂νxlα+1,βm)wn+(α+1)+(∂νxlα+1,βk−1,∂νxlα+1,βm−1)wn+(α+1)] =

To complete the proof of (2.23), we proceed by induction on By (2.3) we get

 aα,β1(lα,βk,lα,βm)=β22[β−(α+1)−1γα+1k,0+β−(α+1)−1γα+1k−1,0]δk,m=β−α2[γα+1k,0+γα+1k−1,0]δk,m, (2.26)

if either (a). and ; or (b). and This exactly gives (2.23) for with .

Assume that the result (2.23) for with holds. We now verify the result with Clearly, by (2.3) we have

 aα,βp+1 (lα,βk,lα,βm)=β22[βp−(α+1)−1γα+1k,p+βp−(α+1)−1γα+1k−1,p]δk,m=βp−α2[γα+1k,p+γα+1k−1,p]δk,m,

if either (a). and ; or (b). and . This statement implies the result (2.23) for with This ends the proof. ∎

The normalization constants and the eigenvalues are closely related. In effect, for any , , we get that

 (2.27)

where the third inequality sign is obtained by integration by parts combined with (2.13).

Moreover, for sufficiently large , an induction procedure starting with (2.20) reveals

 λαk,n=(2k+α+1)n+O(kn−2), (2.28)

which implies

 γαk,nγαk,n+1=22k+α+1+O(k−3). (2.29)

The following eigenvalues and normalization constants are of our particular interest,

 λ0k,0=1,λ−nk,n=2n(k−n+1)n,k,n∈N0, (2.30) γ0k,0=1,γ−nk,n=12nmin(k,n)∑ν=0(nν),k,n∈N0. (2.31)

## 3. Fully diagonalized spectral methods

In this section, we propose the fully diagonalized spectral methods using generalized Laguerre functions for solving differential equations on the half line. The main idea is to find a system of Sobolev orthogonal functions [8, 21] with respect to the coercive bilinear form arising from differential equation, such that both the exact solution and the approximate solution can be explicitly expressed as a Fourier series in the Sobolev orthogonal functions. Although we only consider in this section non-homogenous Robin/Drichlet boundary value problems of a second order equation, one can extend the fully diagonalized spectral methods for solving partial differential equations of an arbitrary high order.

### 3.1. Robin boundary value problems

Consider the second order elliptic boundary value problem:

 {−u′′(x)+γu(x)=f(x),γ≥0,x∈Λ,−u′(0)+μu(0)=η,limx→+∞u(x)=0,μ≥0. (3.1)

A weak formulation of (3.1) is to find such that

 Aγ,μ(u,v):=μu(0)v(0)+(u′,v′)+γ(u,v)=(f,v)+ηv(0),∀v∈H1(Λ). (3.2)

The Lax-Milgram lemma guarantees a unique solution to (3.2) if .

Let

 XβN:={e−12βxp(x):p∈PN}={l−1,βk:0≤k≤N}.

The generalized Laguerre spectral scheme for (3.1) is to find , such that

 Aγ,μ(uN,vN)=(f,vN)+ηvN(0),∀vN∈XβN. (3.3)

For an efficient approximation scheme, one usually chooses the generalized Laguerre functions as the basis functions for problem (3.3). However, we are eager for an ideal approximation scheme whose (total) stiff matrix, in analogue to the Fourier spectral method for periodic problem, is diagonal. Obviously, the utilization of the basis functions leads to a tridiagonal algebraic system. To this end, we shall construct new basis functions which are mutually orthogonal with respect to the Sobolev inner product instead of defined in Theorem 2.1.

###### Lemma 3.1.

Let be the Sobolev orthogonal Laguerre functions such that and

 Aγ,μ(Rβk,Rβm)=ρkδk,m,k,m∈N0. (3.4)

Then satisfy the following recurrence relation,

 Rβ0(x)=l−1,β0(x),Rβk(x)=l−1,βk(x)−dk−1Rβk−1(x),∀ k≥1, (3.5)

where and

 dk−1=β4ρk−1−γβρk−1,ρk=−d2k−1ρk−1+β2+2γβ,k≥1.
###### Proof.

By the orthogonality assumption (3.4) of ,

 l−1,βk(x)=Rβk(x)+k−1∑m=0Aγ,μ(l−1,βk,Rβm)ρmRβm(x).

Meanwhile, by (3.2) and (2.21), for any ,

 Aγ,μ(l−1,βk,Rβm)=a−1,β1(l−1,βk,Rβm)+μl−1,βk(0)Rβm(0)+(γ−β24)(l−1,βk,Rβm).

Both the first and the second terms in the righthand side above are zero due to the orthogonality relation (2.23) of and the homogeneity boundary condition (2.13) for . Further by (2.15) and the orthogonality relation (2.22) for ,

 Aγ,μ(l−1,βk,Rβm)= (γ−β24)(l