﻿ The Asymptotic Eigenvalues of First-Order Spectral Differentiation Matrices Journal of Applied Mathematics and Physics, 2014, 2, 176-188 Published Online April 2014 in SciRes. http://www.scirp.org/journal/jamp http://dx.doi.org/10.4236/jamp.2014.25022 How to cite this paper: Wang, J. and Waleffe, F. (2014) The Asymptotic Eigenvalues of First-Order Spectral Differentiation Matrices. Journal of Applied Mathematics and Physics, 2, 176-188. http://dx.doi.org/10. 4236/ja mp.2014.25022 The Asymptotic Eigenvalues of First-Order Spectral Differentiation Matrices Jue Wang1, Fabian Waleffe2 1Department of Mathematics, Union College, Schenectady, USA 2Department of Mathematics, University of Wisconsin-Madison, Madison, USA Email: wangj@union.edu Received January 2014 Abstract We complete and extend the asymptotic analysis of the spectrum of Jacobi Tau approximations that were first considered by Dubiner. The asymptotic formulas for Jacobi polynomials NP(,),, 1>−αβαβ are derived and confirmed by numerical approxi- mations. More accurate results for the slowest decaying mode are obtained. We explain where the large negative eigenvalues come from. Furthermore, we show that a large negative eigenvalue of order  appears for10−< <α; there are no large negative eigenvalues for collocations at Gauss-Lobatto points. The asymptotic results indicate unstable eigenvalues for 1>α. The eigenvalues for Legendre polynomials are directly related to the roots of the spherical Bessel and Hankel functions that are involved in solving Helmholtz equation inspherical coordinates. Keywords Asymptotic Analysis, Spectral Approximations, Jacobi Polynomials, Collocation, Eigenvalue s 1. Introduction Pseudospectral methods were developed to solve differential equations, where derivatives are computed numer- ically by multiplying a spectral differentiation matrix . Compared to finite difference methods that use local information, pseudospectral methods are global, and have exponential rate of convergence and low dissipation and dispersion errors. However in boundary value problems, they are often subject to stability restrictions . If the grid is not periodic, the spectral differentiation matrices are typically nonnormal , and the nonnormality may have a big effect on the numerical stability and behavior of the methods. On a grid of size N, pseudos- pec-tral methods require a time step restriction of 2()ON−  for hyperbolic and 4()ON−  for parabolic prob-lems. Dubiner  carried an asymptotic analysis for theone-dimensional wave equation. He pointed out the ( )ONboundedness of the eigenvalues of the spectral differentiation matrix with collocation at Legendre points. It was also proposed in  that the 2()ON eigenvalues could be shrunk to ( )ON byreplacing Chebyshev Tau method with J. Wang, F. Waleffe 177 Legendre Tau method. That would mean a time step size increase from 2()ON− to 1()ON−. However, the eigen- values of Chebyshev and Legend respectral differentiation matrices are extremely sensitive to rounding errors and other perturbations . On agrid of size N, machine rounding could lead to errors of size 2()ON. Thus the Le- gendre Tau method which has an 1()ON− time step restriction in theory, is subject to an ( )2ON− restriction in practice. The slowest decaying mode and largest wave numbers are often of interest. They affect stabilities and limit time step sizes of pseudospectral approximations. It is reported in   that there exists a large negative eigenvalue of order 2N for Chebyshev spectral differentiation matrix. Where does that eigenvalue come from? When does it ap- pear? And how does it affect the time step size? In this paper, we will consider the first-order spectral differentiation matrix and examine the behavior of its eigenvalues asymptotically and numerically. Let’s consider the first-order eigenvalue problem ( ), 11,10du u xudxλ=−<<= (1) The spectrum of the differentiation operator is empty. However, the spectrum of the Tau approximations to the eigenvalue problem affects the stability of the associated wave equation txuu= on ( 1,1)− with boundary condi- tion ( )1, 0ut=, whose solution is a leftward translation at speed 1 . Let ()Nux be a polynomial approxima- tion of degree N. ()Nuxsatisfies ()( ), 10N NNNuDuPx uλ−= = where /Dddx≡,()NPx is a polynomial of degreeN. We can write ()( )xxsNNauxeePsdsλλ−= −∫ for some arbi- trary constant a. Integration by parts gives ( )( )()00()1kkxaNNNNNkkkkDP xDP aeuxλλλλλ−= == −∑∑. (2) The non-polynomial contribution (the second term) is eliminated by picking a=−∞ if 0λR. We then obtain ( )( )01kNNNkkDP xuxλλ==∑ (3) The boundary condition ( )10Nu= implies ( )( )110, 0sNePs dsλλ−−∞−=<∫R (4) ( )( )110, 0sNePs dsλλ∞−= >∫R (5) In the Tau method , the polynomial approximation ()Nux is determined from the boundary condition ( )10Nu= and the requirement that ()NPx is orthogonal to all polynomials of degree 1n− with respect to the weight function ( )0wx≥ in the interval ( 1,1)−. For Jacobi weight function ( )()(,)1(1 )wxx xαβαβ=−+,( )( )(,)NNPx Pxαβ=, the Jacobi polynomial of degree N with ,1αβ>−. Assume 0λ≠, otherwise 0λ= corresponds to the trivial solution ( )0Nux=. From (3), the boundary condition ( )10Nu= leads to the characteristic polynomial ( )0Φ(1)NkkNNkDPνν==∑, with 1/νλ=. It is proved  that the eigenvalues lie in the left half-plane for Jacobi polynomial ( )(,)NPxαβ if 11α−< ≤ and 1β>−. The eigenvalues are computed numerically using the three term recurrence relation for corresponding Jacobi polynomials  (see Figure 1). We will show that the theorem is sharp by obtaining asymptotic results indicating unstable eigenvalues for 1α>. In order to obtain the asymptotic behavior for large values of λ, we use the method of steepest descents to deform the integration paths to obtain the dominant contribution from saddle points. In general there are two saddle and two boundary points. The balance between the dominant saddle and boundary contributions leads to an asymptotic equation for λ. The two saddle points collide to form a third-order saddle when /Niλ= ±; nearly merge when 0α= and /1Niλ±; and are too close to the boundary points when ()ONλ>. These cases complicate the analysis. J. Wang, F. Waleffe 178 Figure 1. The limit curve and eigenvalue approximations for Cheby- shev polynomials (12αβ−= =). Thick solid: ( ), 10Fµσ=−=, dash: ( )3 ln1,12NFONNµσ=−=− −, dash-dot: ( )3 ln,1 2NFNµσ== −1ON−, where ()2211,1 lnFσµµσµ σµµ−+=+ ++R. Thin solid: contours of ()Bµ. Take dash-dot branch for ( )0Bµ>, dash for ( )0Bµ<. The paper is organized as follows. In Section 2, we present the asymptotic analysis and numerical results for Chebyshev polynomials. We generalize the results to Jacobi polynomials in Section 3, and derive the approxi- mations of the slowest decaying mode and largest wave numbers. In Section 4, we show that the eigenvalues for Legendre polynomials are directly related to the roots of spherical Bessel and Hankel functions. The analysis and numerical results for collocation methods are explained in Section 5. Finally we conclude in Section 6. 2. Chebyshev Polynomials We use Chebyshev polynomials ( )( )NNPs Ts≡ to illustrate the approach and derive the asymptotic formulas. They correspond to the class of Jacobi polynomials with 1/2αβ== − and are especially relevant in practice. With a change of variable cossφ=, (cos )cosNTNφφ=, and (4) becomes 0coscos( sin)0ieNdλφπφ φφ−+∞−=∫ (6) It suggests that Nλµ=. Define ( )cos iρφµφ φ=−+, and I twice the above integral, () 1202siniiNIed IIππρφπφφ+∞ −∞=+=+∫∫ (7) We will apply the method of steepest descent to deform the integration path without changing the value of the integral, so that it goes through the critical point (saddle point) *φ in such a way that *()ρφR is maximum along the path, and *()ρφR decreases along either direction away from *φ as rapidly as possible. As N→∞, the dominant contributions come from the saddle points and boundary points 0 and 2π. 2.1. Saddle Contributions The saddle points of ( )ρφ satisfy the equation '( )sinsin0iρφ µφ= +=, which implies *sin /iφµ= −. The steepest-descent curves of ( )Neρφ vary as μ varies. We consider two cases: 0µ≠R and 0µ=R. 2.1.1. 0≠µR There are two saddle points ()*σφ, ()*sign 0σσφ= ≠I, and ()*"()0σρφ≠. They satisfy the relationship J. Wang, F. Waleffe 179 ( 1)(1)**φ πφ−= −, ()( )11**( )()ρφ ρφ−= −RR, ()( )11**( )()ρφπ ρφ−= −II, and ( )( )()()22* **11cos1 lniσ σσσµρφµ φφσµµ−+=−+=++ (8) The integration paths follow the constant-phase contours( )*()( )σρφ ρφ=II, i.e . the steepest-descent curves, from 0 to iπ+∞ for I1 and from 2π to iπ−∞ for I2, passing through the saddle points ()*σφ. The steepest curves corresponding to the opposite signs of µI are mirror images about φπ=R axis. Define ( )( )( )( )( )211 2** 2112 1ln11Bµµ ρφρφµµ−−+≡−= ++++RR R (9) and redefine σ to be sign ()Bσµ=. Then *()sIσφ dominates when () 0Bµ≠. Following the standard saddle point approximation ( §7.3), we evaluate the integrals near the saddle points and obtain the dominant contribution ( )( )( )( )*()1/2**()*2~sin "( )Nss iI INeeσσσρφσθσπφφρφ−= ( 10) As N→∞ when () 0Bµ≠. If ( )0Bµ=, ( )()(1)(1)**ss sII Iφφ−= +. 2.1.2. 0Rµ= Let ibµ=. The two saddles points collide and form athird-order saddle point at iµ= ±. 1) 01b<<. In this case ()0Bµ<. ( )11*/ 2cosh(1/)ibφπ±−=±− when 10b−<<, 13/2cosh(1/)ibπ−± when 01b<<. ( )( 1)*ssIIφ−=. 2) 1b>. We have ( )0Bµ=. ( )11*sin(1/)bφ−= −, ( )11*sin(1/)bφπ−−=−−when 1b<−. ( )11*2sin(1/)bφπ−=+−, ( )11*sin(1/)bφπ−−=−− when 1b>. ( )(1)*ssIIφ=. 3) 1b= ±. */2bφππ= + is a third-order saddle point. A similar saddle point approximation leads to( )*1123332*~()6/ 3Γ(1 /3)biiiNsINe eeππ πφφ−−−−+. 2.2. Boundary Contributions We approximate the contour at 0φ= by the straight line ivφ= with [0, ]vε∈ and obtain the dominant contribution at 0φ=, ( )()()()( )'00022'00~~ 0NNvNbeIeevdv Nρερρρ−−−−∫ ( 11) as N→∞. The approximation is the same near 2φπ=. Therefore, ( )( )( )022'0(2) ~20Nbb beII INρπρ−=+− (12) 2.3. Balance between Saddle and Boundary Contributions There are four possible balances from Section 2.1: 1) 0µ≠R and sign ()0Bσµ= ≠, 2) 0µ=R and ( )0Bµ=, 3) 0µR, and the limit curve is a reflection by the imaginary axis. 2.4. The Large Negative Eigenvalue The largest wave number limits the time step size of the numerical approximation. When ()ONλ>, the sad- dle points ()11*sinh Niφλ−=− and ( )21*2sinh Niφπ λ−= +− are too close to the boundaries 00,2φπ=. Thus we have to take the integration path of I1 from 0 up to ( )1*φ, then to iπ+∞, and that of I2 from 2π up to ( )2*φ, then to iπ+∞, going down passing through ( )1*φ− to iπ−∞. Both ( )1*φ and ( )2*φ are dominant sad- dles and ( )1*φ− is subdominant. The integral then becomes J. Wang, F. Waleffe 181 || cos|| 1202sinNiiiIed IIλφφππ λπφφ++∞ −∞= ++=+∫∫ (17) The balance between saddle and boundary contributions implies a polynomial equation for 2/| |xNλ=, 100!2 (21)!2nnnnnxxnn n+∞=−=+∑ (18) The real root of (18) approaches a constant as the degree of the equation increases. Solving the equation of degree 5 gives 2| |/~Nλ0.5855. This large negative eigenvalue is plotted in Figure 2. The approximation agrees well with numerical results. 3. Jacobi Polynomials We now generalize the asymptotic analysis in Section 2 to Jacobipolynomials ( )(,)cosNPαβφwith weight func- tion ()()cos1cos(1cos) ,,1wαβφφφ αβ= −+>−. 3.1. Approximation of Jacobi Polynomials We approximate Jacobi polynomials in two regions, i. e. near 0 and away from 0. Figure 2. The large negative eigenvalue for Chebyshev polynomial, ~λ0.58552N, confirmed by numerical approximations atN=[10:10:80], indicated by *. Figure 3. ○: the slowest decaying eigenmodes for Chebyshevpolynomial giv- en by (16) (left) and Jacobi polynomialwith 0.3, 0.5αβ=−= (right). ◊: Du- biner’s approximation ( Equation (8.6)). *: numerical approximations. The limit curves move to the rightas N increases. J. Wang, F. Waleffe 182 1) If ||nφ πε−≥, ( )( ),1/2 11cos~ (2)()()iNiiN iNPe DeeDeαβ φ φφφφπ−−− −−+ (19) where ( )1112 222(1 )(1)Dzz zαβ αβ++− ++= −+, defined for 1,1z zR<<≤ <∞, || 1l(m)i zDz±→ exists. We de- rived it from Szego  Theorem 12.1.2. 2) If |0| 1φ−, ( )( )( )11,222sincos~ (2)2iN iiNiNPeee eαβααβ φγφγφ φπ−−+−−+ (20) where (1/2) /2γα π=−+. We derived it from Szego  Theorem 12.1.4. 3.2. Limit Curve Using approximation (19), the integral in (4) becomes () 102() siniiNiIeDedππρφ φπφφ+∞ −∞−−= +∫∫, (21) with the same ( )cos iρφµφ φ=−+. A similar analysis gives the saddle contribution ( )()( )( )()**()1/2 *()*sin 2~"( )NsiiI NeeDeσσσσρφ θσφφπρφ−− (22) This is in agreement with the results for Chebyshev polynomials in Section 2, where () ( )1/2siniiDei eαφφφ+−−= for αβ=. Using the asymptotic formula for Bessel function of the first kind , ( )()13222() cos()J zzOzzαγπ−= ++ as ||z→∞ with ( )arg,0zπ δδ≤− >, we obtain ( )()1120002() cos()J NNONNαφ φγπφ−= ++ As N→∞ for 1/30|0|Nφ−−. From (20), we have ()( )( )1,22cos~ 2()NPN NJNαβ ααβ ααφ φφ−+− ( 23) Plug (23) to the integral in (4) and we obtain the boundary contribution 3112 222(0) ~2(2)Γ()bNI Neαβ ααµπα−−−− (24) (0) ~ 0bI when 0α=. 3.2.1. 0≠α The dominant balance is between saddle and boundary contributions. This gives the same limit curve as (14). The slowest decaying eigenmodes are plotted in Figure 3 for 0.3, 0.5αβ=−=. The theorem in  proves that the eigenvalues lie in the left half-plane for Jacobi polynomial ( ),NPαβ if 11α−< ≤. We have derived asymptotically that the eigenvalue becomes unstable for 1α>. 3.2.2. 0=α The dominant balance is between two saddle contributions from ( )1*φ±, and this leads to a new limit curve 222112 1ln011µµµ−+++ =++R (25) J. Wang, F. Waleffe 183 The intervals [, )ii±∞ are excluded because only one saddle point contributes. The limit curve and eigenvalue approximations for Chebyshevpolynomials of the 2nd kind (1/2αβ= =) and Legendre polynomials (0αβ= =) are plotted in Figures 4 and 5. The eigenvalues huddle around the limit curve. Note that even at a small N = 10, the eigenvalue approximations for Legendre polynomials lie exactly on the dashed curve, which is true for all the other figures. Figure 4. The limit curve and eigenvalue approximations for Cheby- shev polynomials of the 2nd kind (1/2αβ= =). Thick solid: ( ), 10Fµσ=−=, dash: ( )1 ln1,12NFONNµσ=−=− −, dash-dot: ( )1 ln1,1 2NFONNµσ==−−. Figure 5. The limit curve and eigenvalue approximations for Legendre poly- nomials (0αβ= =). Thick solid: 222112 1ln011µµµ−+++ =++R, dash: 2222211 1112 1lnln211 11Nµµµµµ−+ ++++ =++ −+R. J. Wang, F. Waleffe 184 3.3. Olver Paths When 0α= and | |1iµ±, the two saddle points ( )11*sin iiφµ−= − and ( )11*sin iiφπ µ−−=−− nearly merge. They coincide at */2φπ= when iµ= −, and 3 /2π when iµ=. The expansions given by the or- dinary method of steepest descents are not uniformly valid for | |1iµ±. We need to construct uniform ex- pansions near i± and deform the contours into Olver paths . Apply the cubic transformation, ()()31, ()3Aρφµζη µζµ=−+, to map the saddle points of ( )ρφ in the φ-plane to the saddle points η± in the ζ-plane. It is shown  that the mapping is uniformly regular and one-to-one for all µ in a neighborhood of i±. The integral (21) becomes 31() 13() sinsNAidI eDeddζ ηζφφφζζ±−+ −−=∫ ( 26) where ± is a contour running from −∞ to /3ieπ±∞ when 1||iNµ−±, through the saddle point ζη=. The first order approximation gives 514 21263333 6~2Ai2()iN iisIN eNeiππ πππµ−±± (27) for 1||iNµ−±, where Ai is the Airy function. The eigenvalues are related to the zeros of the Bessel and Hankel functions of half-integer order 1/2N+ . Further discussions are presented in Section 4. The slowest decaying eigenmode’s wave number satisfies ( )21 1336Ai2() 0iNei ONπµ±−±+ =. (28) Ai (z) has a maximal zero 00,0cc−> . Then Equation (28) gives 11 233 630~2() [(1)]iNNecONiNOπλµ−−=−+ +. (29) This is also the largest wave number (see Section 3.4). Hence when 0α=, the largest wave number is of or-der N instead of 2N. 3.4. The Large Negative Eigenvalue The formulas for saddle points contributions are not valid when the saddle points are too close to the boundaries at ()ONλ>. The approximation of Jacobi polynomials has to be replaced by Equation (20). If λ is not real, we can take the same integration paths as in Section 3.2 and obtain the saddle contribution. The boundary con- tribution remains the same with λ used first and then replaced by Nµ. The balance between saddle and boundary contributions when 0λ≠, or between two saddle contributions when 0λ= gives the same limit curve as before. If λ is real negative, we take the integration paths as in Section 2.4. The balance leads to a polynomial equ- ation for /| |xNλ=, 2 2210 0011cos1sin 0222 !2!2 (21)nnnnnnn nnn xxcx n nnnαα πααπ+∞ ∞∞−== = −  ++−−+=  +  ∑ ∑∑ (30) where ( )2/221!!nckπ= − when n is even 2nk=; 2knck= when n is odd 21nk= +. The lowest order approximation by solving Equation (30) of degree 2 gives 2/ 2cot(/2)cot2(3/ 2)(1/ 2)~3/2xπαπ παπααα−+++ −+ (31) For 10α−< <. Figure 6 shows the ratio 2/Nλ obtained from (31) and by solving Equation (30) of degree 4, respectively, confirmed by the numerical eigenvalues at 40N= for different values of α. The real root of the polynomial of degree > 2 is negative for 0 1/2α<< and 1/x approaches −∞ as α approaches 1/2 from the left. Thus, there is a large negative eigenvalue of 2()ON for 10α−< < and no large negative ei- J. Wang, F. Waleffe 185 Figure 6. The large negative eigenvalue for 10α−<< scaled by 2N, compared with the numerical approximations for Jacobipoly- nomials at 40N= and α= [−0.001, −0. 1: −0.1: −0.9], indicated as *. Solid: from solving Equation (30) of degree 4, dash: from Equa- tion (31). genvalue for 0 1/2α<≤. The integral diverges for 1/2α>. Therefor there is a large negative eigenvalue of 2()ON for 10α−< <. When |0| 1α−, the assumption ()ONλ> is no longer valid. In this case, the negative eigenvalue is simply described by the zero of the limit curve given in (14). When 0α=, the balance is between two saddle points. It gives same limit curve as (20). This is not consistent with the assumption ()ONλ>. So (29) also gives the largest magnitude eigenvalue which is ()ON. 4. Legendre Polynomials and Spherical Bessel Functions Denote Legendre polynomials by ( ),NNPPαβ=, 0αβ= =. As indicated in Section 3.3, the eigenvalues are di- rectly related to the roots of spherical Bessel and Hankel functions. The emission and scattering of electromagnetic radiation involves solving the vector wave equation . The solutions are then expressed as an expansion in the orthogonal spherical waves, known as the multipole expan- sion. Consider the scalar wave equation 210UUtc∂∆−=∂, where (,)UU t=x is the velocity potential, (,,)xyz=x and c is the speed of sound. For time-harmonic acoustic waves of the form ( ){ },()itU tueω−=xxR with frequency 0ω>, each Fourier harmonic ()ux satisfies the Helmholtz equation 20u ku∆+ = where /kcω= is the wave number. It can be rewritten in terms of spherical coordinates(, ,)rθϕ. With the expansion ( )(,)nnuf rYθϕ=∑, the radial function ( )nfr satisfies the equation 2 22"2' [(1)]0nn nrfrfkrn nf+ +−+= (32) known as the spherical Bessel differential equation  (10.1.1). With a change of variable ( )( )nng rrf r=, ng satisfies the Bessel differential equation of half integer order 1/2n+. Now let /riz k= and replace nf by w, n by N, then Equation (32) becomes the modified spherical Bessel equation  (10.2.1) 22"2' [(1)]0z wzwzNNw+ −++= It has a solution, the modified spherical Bessel function of the third kind  (10.2.4), ( )( )122NNkzK zzπ+=, where ( )12NKz+ is the modified Bessel function of the third kind of order 1/2N+. Using  (9.6.23) and integration by parts N times, we obtain J. Wang, F. Waleffe 186 ()( )()101zszztNN NkzePs dseePtdt∞∞− −−== +∫∫ Similar to deriving Equation (3), ( )( )01kzNNNkkDPekz zz−==∑ Thus the eigenvalues λ in (3) are the roots of ( )Nkz, which are also the roots of Hankel function of order 1/2N+, spherical and modified spherical Hankel functions using the relations  (10.1.1),(10.2.15). 5. Collocation The Chebyshev collation method and stability analysis was discussed in , and then extended to general Gauss-Radau collocation methods  for the one-dimensional wave equation. The Gauss-Radau and Gauss- Lobatto Jacobi methods were compared and it is shown that the latter is asymptotically stable for 1αβ≤≤ . In particular the only Gauss-Lobattoultraspherical method which is marginally stable corresponds to 1αβ= =. We take the discretization 0111Nxx x−=< <<= of the interval [−1,1], and approximate ()ux by a po- lynomial ()Nux of degree N that satisfies ( )11,for,,1 0N NNNDuuxxxuλ−= == (33) ()Nuxinterpolates the values 01 1, ,,0Naa a− at 0111Nxx x−=< <<= and can be written as ( )10()NN jjJuxal x−==∑, where { }()jlx are the cardinal functions. Substituting Nu into (33) gives the matrix system with eigenvalue λ. Collocation at zeros of orthogonal polynomials (Gauss points) is identical to the Tau method with Gauss quadrature. We consider three sets of Gauss-Lobattopoints 11[, ]TNxx−=x: 1) Chebyshev extreme points:cosjjxNπ= −; 2) Chebyshev extreme points of the 2nd kind: jxj= th zero of 'NU; 3) Legendre extreme points: jxj= th zero of 'NP. They can all be described as zeros of '( )NQx, where NQ is Chebyshev, Chebyshev of the 2nd kind and Le- gendre polynomial of degree N, respectively. Defineresidual ()() ()NNRxDu xλ= −. It is a polynomial of degree N that vanishes at 1x=− and zeros of '( )NQx. Thus NR takes the form ()(1)'( )NNRxxQx= +. In general, for Jacobi polynomial ( ),NPαβ, from , ( )()( )()( )( ),1, 1111 112NN NdRxxPNxPdxαβα βαβ++−= +=++++ Using (19) and (20), we obtain approximations of ( )cosNRφ for φ away from 0 and near 0. The balance remains valid with ()()( )colcol col22, 1,11iiiD eDeeφφφα αββ= =+=++ However, for Legendre extreme points, the balance is now between saddle and boundary contributions sincecol 1α=, and the limit curve is the same as the other two. Figure 7 shows the eigenvalues obtained using the Tau and collocation methods, together with their limit curves for Chebyshev (1/2αβ== −), Chebyshev of the 2nd kind (1/2αβ= =) and Legendre (0αβ= =) polynomials. Collocations at Chebyshev and Legendre ex-treme points are stable. Collocation at Chebyshev extreme points of the 2nd kind becomes unstable. In general, collocation at Gauss-Lobatto points is stable when 11α+≤, i.e. 0α≤, unstable when 0α>. There is no large negative eigenvalue with Chebyshev extreme points used. This can be derived asymptotically using the same method as before, or simply follows from col 1αα= +. 6. Conclusions We have presented pseudospectral approximations of the first-order spectral differentiation matrices with zero J. Wang, F. Waleffe 187 Figure 7. From left to right: Tau method at 1/2αβ==−, collocation at Chebyshev extreme points, Tau method at 1/2αβ= =, collocation at Chebyshev extreme pointsof the 2nd kind, Tau method at 0αβ= =, collocation at Legendre extreme points. boundary condition. We use the method of steepest descent to deform the integration path without changing the value of the integral, in order to obtain the dominant contributions from saddle and boundary points. This ap- proach leads to the asymptotic formulas for the eigenvalues of Jacobi Tau method. Numerical examples with Chebyshev of the 1st and 2nd kind and Legendre polynomials are presented. They agree well with the asymptotic analysis, even at small sizes. The approximations for the slowest decaying modes give more accuracy than those obtained in . The larg- est wave numbers have raised interest in stabilities in pseudospectral approximations. We show that a large neg- ative eigenvalue of order 2N appears for 10α−< <. The collocation methods are also examined. There are no large negative eigenvalues for collocations at Gauss-Lobatto points. For Jacobi polynomials, the eigenvalues lie in the left half-plane if 11α−< ≤. We show that the theorem is sharp by obtaining asymptotic results that indicate unstable eigenvalues for 1α>. The eigenvalues for Legen- dre polynomials are related to the roots of Bessel and Hankel functions of half-integer order, spherical and mod- ified spherical Bessel and Hankel functions. These results complete Dubiner’s earlier analysis . References  Canuto, C., Hussaini, M.Y., Quarteroni, A. and Zang, T.A. (1988) Spectral Methods in Fluid Dynamics. Springer, New York. http://dx.doi.org/10.1007/978-3-642-84108-8  Gottlieb, D. (1981) The Stability of Pseudospectral Chebyshev Methods. Mathematics of Computation, 36, 107-118. http://dx.doi.org/10.1090/S0025-5718-1981-0595045-1  Trefethen, L.N. and Embree, M. (2005) Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Opera-tors. Princeton University Press.  Trefethen, L.N. and Trummer, M.R. (1987) An Instability Phenomenon in Spectral Methods. SIAM Journal on Nu-merical Analysis, 24, 1008-1023. http://dx.doi.org/10.1137/0724066  Weideman, J.A.C. and Trefethen, L.N. (1988) The Eigenvalues of Second-Order Spectral Differentiation Matrices. SIAM Journal on Numerical Analysis, 25, 1279-1298. http://dx.doi.org/10.1137/0725072  Dubiner, M. (1987) Asymptotic Analysis of Spectral Methods. Journal of Scientific Computing, 2, 3-31. http://dx.doi.org/10.1007/BF01061510  Tal-Ezer, H. (1986) Spectral Methods in Time for Hyperbolic Equations. SIAM Journal on Numerical Analysis, 23, 11-26. http://dx.doi.org/10.1137/0723002  Jackiewicz, Z. and Welfert, B.D. (2003) Stability of Gauss-RadauPseudospectral Approximations of the One-Dimensional Wave Equation. Journal of Scientific Computing, 18, 287-313. http://dx.doi.org/10.1023/A:1021121008091  Csordas, G., Charalambides, M. and Waleffe, F. (2005) A New Property of a Class of Jacobi Polynomials. Proceedings of the AMS, 133, 3351-3560.  Weideman, J.A.C. and Reddy, S.C. (2000) A Matlab Differentiation MatrixSuite. ACM Transactions on Mathematical Software, 26, 465-519. http://dx.doi.org/10.1145/365723.365727  Arfken, G.B. and Weber, H.J. (1995) Mathematical Methods for Physicists. Academic Press. J. Wang, F. Waleffe 188  Szego, G. (1939) Orthogonal Polynomials. AMS Colloquium Publication, 23.  Waston, G.N. (1995) A Treatise on the Theory of Bessel Functions. 2nd Edition, Cambridge University Press.  Olver, F.W.J. (1970) Why Steepest Descents? SIAM Review, 12, 228-247. http://dx.doi.org/10.1137/1012044  Chester, C., Friedman, B. and Urse l l, F. (1957) An Extension of The Method of Steepest Descents. Proc. Cambridge Philos. Soc., 53, 599-611. http://dx.doi.org/10.1017/S0305004100032655  Driver, K.A. and Temme, N.M. (1999) Zero and Pole Distribution of Diagonal Padé Approximants to the Exponential Function. Questiones Mathematicae, 22, 7-17. http://dx.doi.org/10.1080/16073606.1999.9632055  Abramowitz, M. and Stegun, I.A. (Eds.) (1965) Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables. Dover Publications.  Colton, D. and Kress, R. (1998) Inverse Acoustic and Electromagnetic Scattering Theory. 2nd Edition, Springer. http://dx.doi.org/10.1007/978-3-662-03537-5  Doha, E.H. (2002) On the Coefficients of Eifferentiated Expansions and Derivatives of Jacobi Polynomials. J. Phys. A: Math. Gen., 35, 3467-3478. http://dx.doi.org/10.1088/0305-4470/35/15/308