NASA Contractor Report 178041

lCASE REPORT NO. 86-2

NASA-CR-178041

19860011765

leASE

STABILITY ANALYSIS OF SPECTRAL METHODS

FOR HYPERBOLIC INITIAL-BOUNDARY VALUE SYSTEMS

David Gottlieb

Liviu Lustman

Eitan Tadmor

Contract Nos. NASl-17070 and NASl-18107

January 1986

INSTITUTE FOR CO~WUTER APPLICATIONS IN SCIENCE AND ENGINEERING

NASA Langley Research Center, Hampton, Virginia 23665

Operated by the Universities Space Research Association

NI\SI\

J .~.~' ,~" ~~ ~.~~ f j :-.':

~"'" ..: \ r

L/\hCJ:.-L 1 '~ •• J .0-'

National AeronautiCS and

Space AdminIstration

Langley Research center

Hampton,VirginIa 23665

STABILITY ANALYSIS OF SPECTRALMETHODS

FOR HYPERBOLICINITIAL-BOUNDARY VALUE SYSTEMS

David Gottlieb

Tel-Aviv University

Institute for Computer Applications in Science and Engineering

Brown University

Liviu Lustman

Institute for Computer Applications in Science and Engineering

Eitan Tadmor

Tel-Aviv University

Institute for Computer Applications in Science and Engineering

Bat-Sheva Foundation Fellow

Abstract

We treat a constant coefficient hyperbolic system in one space variable,

with zero initial data. Dissipative boundary conditions are imposed at the

two points x = ±i. This problem is discretized either by a spectral or

pseudospectral approximation in space. We demonstrate sufficient conditions

under which the spectral numerical solution is stable; moreover, these

conditions have to be checked only for scalar equations. The stability

theorems take the form of explicit bounds for the norm of the solution in

terms of the boundary data. The dependence of these bounds on N, the number

of points in the domain (or equivalently the degree of the polynomials

involved), is investigated for a class of standard spectral methods, including

Chebyshev and Legendre collocations.

Research was supported in part by the National Aeronautics and Space

Administration under NASA Contract Nos. NASI-17070 and NASI-18107 while the

authors were in residence at ICASE, NASA Langley Research Center, Hampton, VA

23665-5225. Additional support was also provided in part by AFOSR-83-0089 for

the first author, and by NSF Grant No. DMS85-03294 and ARO Grant No. DAAG29-

85-K-0190 for the third author while in residence at UCLA.

i

INTRODUCTION

In this paper we study the stability of spectral polynomial methods for

the approximation of initlal-boundary value hyperbolic systems with constant

coefficients and dissipative boundary conditions (Assumption I). We show that

any particular spectral method is stable when applied to a system, if it

satisfies certain conditions (Assumptions II and III) for the corresponding

scalar problem.

Our treatment follows closely the approach of Kreiss [4], leading to an

algebraic condition. However, in the spectral method the requirement that the

approximate solution be a polynomial produces a different algebraic problem

from the standard finite difference discretlzatlons.

We should also review briefly the stability results for spectral methods

and hyperbolic systems. The stability of Chebyshev approximation of the

scalar problem ut = u

/u(x=l,t)x= 0

was proved in [2] for the Galerkln and tau methods, and in [7] for the

collocation method. The same proofs hold for the Legendre method; they are

based on energy estimates, showing that the norm

I

f (l+x) w(x) u2(x,t)dx

-i

is bounded by the data. Here the weight w(x) depends on the discretizatlon

used and equals (I - x2)-I/2 for Chebyshev, w(x) = I for Legendre,

Reyna [8] suggested a stable numerical method for strictly hyperbolic

systems with maximal dissipative boundary conditions, based on high mode
-2-

smoothing. His result applies to problemswith variable coefficients. While

his procedureensures £2 stability,our constantcoefficientsproof requires

no smoothing;instead stabilityis obtainedhere in certain weighted

£2 norms.

The first results about hyperbolic systems with general dissipative

boundaryconditionswere reportedin [5],where it is shown that the numerical

solutionsdecay in time, thus mimickingthe differentialsolutions.

In the present paper we combinethe above results into a general stability

theory for the stabilityof spectralmethods. A future paper will discuss the

convergenceof the spectralapproximationsto the exact solution.

Definitions

Consider the first order hyperbolic system of partial differential

equations

_ _u

(1.1a) ~

2--={ A_x , -I < x _ i, t > O.

Here,

= _(x,t) = (u(1)(x,t),...,u(n)(x,t))"

is the vector of unknowns and A is a fixed nxn coefficient matrix. Since

by hyperbolicity A is similar to a real diagonal matrix, we may assume

without loss of generality that it is diagonal:
-3-

(I. Ib)

A I

AI = • • <0, -- • >0.

a% an

The solution of this system is uniquely determined if we prescribe initial

conditions

(l.lc) _(x,0) = 0, -i < x _ 1

and boundary conditions

_l(-l,t) = L _ll(-l,t) + _l(t)

(l.ld)

_II(l,t)= R BI(1,t)+ _II(t), t > 0.

In these formulas,

= _(t)= (_I(t),_II(t))"

is a given n-vector, and

(l.le ) Nu I = Iu (I) ,..o,u (%)). , u II = Iu (4+I) ,..o,u (n)).

is a partition of _ into its inflow and outflow components - corresponding
-4-

to the partition of A in (1.1b) - while L and R are constant

reflection matrices of order k x (n-k) and (n-k) x k, respectively.

The system (l.la) - (1.1e) is a well-posed problem in the sense that it

satisfies an a priori energy estimate which we now describe. Define the

spatial norm associated with a pair of positive weight functions

=

I i

(l.2a) IIu(x)II2- IIu(x)_I2= S lul(x)l2 ml(x)dx + S lull(x)l2 mll(x)dx. -i -I

Here and elsewhere in the paper we denote by Ivl the Euclidean norms of a

vector v; similarly IAI = suplAvl/Iv I.

Then, for well-posedness one must obtain the following inequality, with some

no _ 0:

(I .2b)

(n- n o ) f e-2(n-n0 )t Ilu(x,t)112dt < const f e-2(n-n0)tls(t)I2 dt

t=0 t=0

for all n > n O •

The parameter n O measures the exponential time growth of the solution. We

are interested in the case where there is no growth in time, and therefore

postulate the following dissipativity requirement:
-5-

Assumption I. There exists a constant _ > 0 such that

(1.3) IRI.ILI _ 1- _ < I.

The inequality (1.3) guarantees that waves originating at one of the

boundaries are not amplified when reflected at the other one. Consequently,

there is no time growth; in Appendix A we show that there is a weight _ such

that for any _ > 0

(1.4) B f e-2_t llu(x,t)tl2 dt < const f e-2nt l_(t)l 2 dt,

t=O t=O

i.e., n 0 in formula (1.2b) may be taken as zero in this dissipative case.

We study polynomial spectral and pseudospectral discretlzatlons of

(i.I). In any such approximation, one seeks a vector of N-degree polynomials

_ _N(X,t) = (v l)(x,t),.--,VN (n)(x.t) • )_

such that

(l.5a) _-_= A-ff_+ Q(x) _.

Here, Q(x) is a diagonal matrix of the form

qll(x) l(n_Z)×(n_A )
-6-

where ql(x), qll(x) are the N-degree polynomials which characterize the

specific (pseudo)spectral method employed, and

is an n-vector to be determined by the set of boundary conditions (l.ld):

I vI(-1,t)= L + gI(t)

(1.5c)

_ll(l,t) R vl(l,t) + _ll(t).

We shall call ql(x), qll(x) the forcing polynomials, since they appear as

inhomogeneous terms in the discretization of the originally homogeneous system

(I.i). Some examples are in order.

In the spectral Galerkin Chebyshev method, one has:

N

I 1

q (x) = _ T0(x) + _ (-I) k rk(x )

k=1

(1.6)

N

1

qll(x) = _ T0(x) + Z Tk(X)

k=l

with Tk(X) denoting the k-th Chebyshev polynomial.

In the pseudospectral Chebyshev method, one may collocate only at the

interior extrema of TN+I(X), yielding:

(1.7) ql(x) = qII(x) = T_+l(X)
-7-

or use as collocation points these extrema together with the downstream

boundaries, in which case:

i I

q (x) = (I-x) T_(x)

(1.s)

qll(x) = (l+x) T_(x).

Similarly, in the pseudospectral Legendre case one has

(1.9) ql(x) = qll(x) = P_+l(X)

or alternatively

(1.10) l q1(x)= (l-x)P (x)

(qll(x) = (l+x) P_(x)

with Pk(X) denoting the k-th Legendre polynomial. These and other examples

are outlined in [2,3].

In this paper we provide a stability study for the (pseudo)spectral method

(1.5). We define stability in terms of an a priori energy estimate analogous

to the differential one:

Definition (I.I). (Stability). The approximation (l.Sa) - (1.5c) is

stable if there exist a weighting pair _(x) and constants = and no > 0

such that for all q > n o we have

(I.11)

_ -2(q-n0)t _ -2(n-n0)t 12 (q-n O ) f e llvN(x,t)U2 dt _ const N 2_ f e l_(t) dt.

t=0 t=0
-8-

We note that the exponential time growth factor here, nO , need not be the

same as the one in the differential estimate (1.2b). The other constant, _,

accounts for a possible algebraically increasing dependence on the

discretization parameter N.

To present our results, let us introduce into (1.5) a new variable

e-_t _N(X,t) and, following the procedure of [4], Fourier-transform the

resulting equations with respect to time. Denote by _ the real dual

variable corresponding to t, and by

i - VN(X,S) =_( e-nt VN(X,t))

(1.12)

^

- g(s) =_Ie -nt _(t)) , s = n + i_

the transforms of e-nt v N and e-nt _, respectively (these functions are set

to zero for t < 0). The resulting equation for v^= C I(x's), v II(x,s))"

is then:

I ^I d$1

s_ = A I aT+_ qI(x) I

(1.13a)

s$11 All aS11~ II

~ dx + qll(x) %

where the n-vector _ _ _(s) = (TI, TII )" is to be determined from the

boundary conditions:

£1(-l,s) = L £11(-l,s)+ _l(s)

(1.13b)
-9-

Using Parseval's relation we conclude:

LEMMA 1.2. The approximation is stable if there exist a weighting pair

m(x) and constants _ and no > 0 such that for all s with

Re s = n > _0 the following holds

(1.14) (n -_0)II£N(X,S)II2 _ const N 2_ l (s)l2

In sections 2 and 3 we discuss necessary conditions for the stability

estimate (1.14). One such condition is the obvious requirement that the

(pseudo)spectral method must be stable for a scalar problem. Another one is

that the homogeneous two-polnt boundary value problem (1.13) should have no

eigenvalues with Re s > _0" In section 4 we show that a strengthened version

of these necessary conditions guarantees stability. The novelty of our

sufficient stability criterion is that it deals exclusively with the

properties of the scalar model and there is no need to consider the

complicated coupling of such scalar equations through the boundary conditions

(1.5c).

In the last sections we demonstrate various important cases in which our

sufficient conditions for stability are satisfied and, thus, deduce stability

theorems for pseudospectral discretizations of the initial-boundary value

system (I.I).
-I0-

2. THE SCALARPROBLEM

The system (l.13a) decouples into n scalar equations of the form

^ dv

(2.1a) sv = a _-_ + q(x)T, -I < x < 1.

^ _N(X '

Here v _ s) and a, stand for any of the corresponding components of

_N(X,S) and the diagonal of A, q(x) equals either ql(x) for a < 0, or

qll(x) for a > 0, and T = T(s) is to be determined by boundary conditions.

In the scalar case, these boundary conditions amount to prescribing the

upstream values:

IVN(-l,s) = g(s), a < 0,

(2.1b) )

VN(l,s) g(s), a > 0.

If formula (1.14) is valid, there must be a similar stability estimate

holds for (2.1). We make the necessary

Assumption II. (Scalar stability). There exist a weighting pair ml(x),

mll(x) and constants = and q0 _ 0, such that for all s with

Re s = n > n o , we have

(2.2a) (_ -q0 ) llVN(X,s)ll21_ const N2_IVN(-l,s)l 2 a < 0

(2.2b) (_ -q0 ) llVN(X,s)ll211 < const N 2= JVN(l,s)J2 , a > 0.
3. THE EIGENVALUE PROBLEM

Connected with (1.13) is the following eigenvalue problem:

Let _ be the space of all polynomials of degree < N. We say that

E _ is an eigenfunction of (1.13) corresponding to an eigenvalue s, if

- I II-.

_ _N(X,S) = i_ ,_ ) is a nontrivial solution of

s_I A I d#I

I = _ + qI(x)!I

(3.1a)

d_II II II

s_II = A ll __dx + q (x)!

with boundary conditions

_I(-l,s)= L_II(-l,s)

(3.1b)

 II(1,sR I(1 ) ,s).

In order to determine whether s # 0 is an eigenvalue, we proceed as

follows. Equations (3.1a) form a system of ordinary differential equations

in x, depending on a parameter s

(3.2) _Slnx n - A d_)_ = Q(x)%.

This can be solved by a formal series

oo

(3.3) _(x,) = Isl - A d)-I Q(x)_ - [ s-k-I Ak Q(k)(x)_. k=O
-12-

We note that since Q(x) consists of N-degree polynomials, the sum on the

right include only the first N+I terms. Rewriting this in the usual

partitioned form, we have

I_ 3.1. Let pl(x,s) _ pl(x,s;Al) and pll(x,s) _ pll(x,s;All)

denote the diagonal matrices

N -k-I

el(x,s) = _ s (ql(xl)(k) (All k

k=0

(3.4a)

N

ell(x, s) = I s-k-I (qll(x)) (k) (All)k"

k--0

Then the polynomial solution of (2.2a) is given by

_l(x,s) = pl(x,s)Tl

(3.4b) {

II

(x,s)= P_(x,s)TI_.

Substituting this expression into the boundary conditions (3.1b), we end

up with the homogeneous linear system:

-ReI(1,s) PII(1,s)/

Denoting the coefficient matrix by D(s)
-13-

(3.6a) D(s) = _ -RPI(I's) l(-l,s) -LpII(-I,s)) pII(I, s)

we arrive at

I_MMA 3.2. A complex number s _ 0 is an elgenvalue of (3.1), if and

only if,

Det[D(s_ = 0.

Suppose now that an elgenvalue s exists, with Re s = n > n0 _ 0. Then

the corresponding eigenfunctlon is a nonvanishlng solution to the homogeneous

problem (1.13) with g = 0, in contradiction to (1.14). For future reference

we state this as

LEMMA 3.3. A necessary condition for the stability estimate (1.14) to

hold is that the elgenvalue problem (3.1) has no eigenvalues s with

Re s > n0 > 0, i.e.,

(3.6b) Det[D(s)] # 0, Re s = n > no > 0.
-14-

4. THE STABILITY ESTIMATE

Here we present a stability theorem for hyperbolic systems with initial

and boundary value conditions.

In the previous two sections we exhibited two necessary conditions for the

stability of (1.5). The first of these, (2.2), is a property of the

individual scalar equations which constitute system (l.5a). The second one,

(3.6), involves the coupling of these equations through the boundary

conditions. Thus, in principle, in order to test this second condition, one

should check the boundary detrminant of D(s) for each pair of the reflection

matrices L and R. However, (3.6) is implied by a much simpler condition,

which deals with the scalar equations separately and does not depend upon L

and R. Moreover, this strengthened condition turns out to guarantee

stability, as we shall presently prove.

We make the following

Assumption III. The polynomial solution of the scalar inflow problem:

dvI

(4.1a) s _I = a d__x__+ql(x), a < 0

satisfies

(4.1b) Ivl(l,s)l _ Ivl(-l,s)l, Re s = n > n o > 0.

Similarly, in the outflow scalar case,

^II dill

(4.2a) s v = a d-_+ qll(x), a > 0,
-15-

we have

(4.2b) Ivll(-l,s)l < Ivll(l,s)l, Re s = n > n0 _ 0.

The estimates (4.1b) and (4.2b) are discrete analogues of differential

estimates, as shown in Appendix A, cf. a.ll.

We can now prove

LEMMA 4.1. Let Assumptions I, II, and III hold. Then, there exists a

constant n o _ 0 such that for all s with Re s = n > n o we have

(4.3) Det[D(s)] # 0.

Proof. Consider the diagonal matrix pl(x,s) of (3.4a). Its j-th

diagonal element p_J)(x,s) is given by

That is, we have

(4.4a) Pl(x,s) = "-o

PN^(%)(x's)

where p(J) - p(NJ)(x,s) satisfies the inflow scalar equation
-16-

dP(j) ql(x)

(4.5a) sp(j) = aj d--x----+ ' aj < 0, i < j < £.

Similarly, the diagonal matrix

(4.4b) pll(x,s) = "..

_(n),

_N _x,s)

consists of solutions to the outflow scalar equations

(4.5b) sp(j) = a. dp(j)-+ qll(x ) a. > 0, £ < j < n. 3 dx 3

Fix s with Re s = n > _0 > 0. We first claim that both matrices P(-l,s)

^(j)

and P(l,s) are nonsingular. Indeed, if for some index j, I _ j < I, PN

vanishes on the left boundary

3 ,s) = O,

then by the scalar stability Assumption I, see (2.2a), the solution of (4.5a)

vanishes everywhere,

li"N("p3)(x,s)ilI = 0.

L0

This, however, contradicts formula (4.5a). In the same way, one shows that

P(l,s) is regular. Hence, the matrix D(s) can be factored (here the ratio

of two diagonal matrices is understood as the diagonal of quotients)
-17-

(4.6)

-RPI(I,s) PII(1,s

-R pl(l's) I pll(l,s)

el(-I ,s)

Next we prove that the matrix E given by

pII(-I, s) 1 l_ I L If(

P 1,s)

(4.7a) E(s) =

R el(lz_ I

Pl(-l,s)

is regular, with a uniformly bounded inverse:

(4.7b) IE-l(s)l < K --1 + l,RI + lIlL, Re s = _ > n0 > 0.

It is here that we make use of Assumptions I and III. Denoting

= R pl(l's)

Pl(-l,s)

(4.8)

i = n pII(-l's)

pII(l,s)
-18-

we have, by (1.3), (4.1b) and (4.2b):

1 1

'",,T- RL)-II' < 1 - IRI ILI <

(4.9)

1 1

-RL)-II'< i- IRIILl<_

l(i

Hence, E can be inverted

(4.10a) E-l(s) = (ii If(ILl 0 (T-0LR)-I)

and the estimate (4.7b) follows:

(4.10b) IE-1(s)< (2 + IRI2 + ILl2)I/2"-_I < I + IRI_ + ILl•

Thus, the matrix D in (4.6) is regular, being the product of regular

matrices, and the lemma is proved.

The estimate (4.10b) is, in fact, the last ingredient needed for a

stability proof for (1.5). We can state

_O_ 4.3. (Stability). Consider the hyperbolic system (l.la) -

(l.le), satisfying the dlssipatlvity Assumption I. Then its (pseudo)spectral

approximation (l.5a) - (1.5c) is stable, provided that Assumptions II and III

are fulfilled.
-19-

Proof. System (l.5a) consists of _ inflow equations

dv (j) I )

(4.11a) sv (j) = a. --+ q (x)T (j 1 _ j _

j dx

and n-I outflow ones

d$ (j) (j)

(4.11b) s$(j) = a. --+ qll(x)T £ < j _ n. j dx

By the scalar stability Assumption II, each of these scalar equations is

stable. That is, there exist a weighting pair m(x) = (ml(x),mII(x))" and

constants _ and no > 0 such that for all s with Re s = n > nO we have

for arbitrary vectors r

(4.12a) (_ -_0)IIvN(J)(x,s)1121 _ const N 2_ IvN(J)(-l,s)I2 1 < j < 4,

(4.12b) (n -_0).IVN(J)(x,s)H211 _ const N2=Iv(NJ)(l,s)12 _ < j < n.

_0

Using the spatial norm of (l.2a)

"$N(X,S)'I2 - HSN(X,S)" 2m = II^IVN(X'S)H21 + ,v i(x,s) 211

these inequalities can be added together and rewritten in concise form

Re s = n > no > 0.
-20-

It remains to estimate the two squared boundary terms inside the brackets on

the right. To this end we rewrite our solution in the form, cf. (4.5) and

(4.11)

I ^I s) P!(x,s)iI

ZN(X, =

(4.14)

^II _N (x,s) pll(x,s)iII

Inserting this in the boundary conditions (1.5c) we find, as in (3.5), that

= _(s) is determined uniquely by

(4.15) II : D-I (s) _II (s)]

Using (4.14), the quantities we want to estimate can be expressed in terms of

the vector TN

From formulas (4.6) and (4.16) we deduce
-21-

and, since by (4.7b) E(s) has a uniformly bounded inverse, we end up with:

(4.18) lSNI(-1,s) 2 + l~N$11(l,s)l 2 _ K] (s)2!l

Re s = _ > _0 _ 0.

Thus, (4.13) and (4.18) add up to (1.14), which is the definition of

stability.

We note that once factorized as in (4.6), the boundary matrix D(s) is

invertible, if and only if the matrix E(s) is. Both Assumptions I and III

were introduced in order to further guarantee that E(s) is uniformly

bounded; in fact, one could assume, instead, that an s-uniform estimate like

(4.7b) holds. This corresponds to the Uniform Kreiss Condition which

characterizes the stability of difference approximations to initial-boundary

value systems. The merit of Assumption III, however, is that it deals with

the scalar problem only, rather than with the intricate coupled boundary

matrix E(s).

5. THE OIEBYS_V, LEGENDRE, AND GEGENBAUER PSEUDOSPECTRAL METHODS

In this section we address ourselves to the question of the stability of

the most commonly used pseudospectral methods, namely Legendre and Chebyshev

collocation. We employ the stability criterion stated in Theorem 4.3, i.e.,

we shall verify that the scalar stability Assumptions II and III hold.
-22-

The scalar model is, of course:

A

^ dvN

(5.1) svN = a d--_ + q(x) T(s)

^

with v given at x = i(-i) if a > 0 (a < 0), respectively.

Our results may be stated in a more general context, involving Gegenbauer

polynomials, of which Chebyshev and Legendre polynomials are particular

a+I/2 are (suitably normalized) cases. The Gegenbauer polynomials CN

orthogonal polynomials with respect to the weight

(5.2) _(x) = (I - x2) _

in the interval -I _ x _ I. Thus, the Legendre polynomial PN is a multiple

of C# , and the Chebyshev polynomial TN is a multiple of CoN.

We shall make extensive use of the Lobatto quadrature rule, which is valid

whenever p(x) is a polynomial of degree < 2N+2:

1 N+I

(5.3) f m(x) p(x) dx = [ mj p(xj) _. > 0.

-I j=O ' 3

Here, x0 = 1 and XN+ 1 = -i denote the endpoints of the interval of

integration and the interior points

(5.4) i > x I > x2 > ... XN_ 1 > xN > -i

3

are the zeros of the Gegenbauer polynomial CN . In particular, for the
-23-

Chebyshev method, these points are the zeros of T_+I, while for the Legendre

method they are the zeros of PN+I - this follows from the identity:

(5.5) d a 2a C_+I d--_CN= N-I "

We consider collocating the equations at the interior points

3

a+ -f

xj, 1 < j < N. Then the forcing polynomial q(x) is proportional to CN

and may be normalized as

_+3

cN 2(x)

(5.6) q(x) - 3 '

C N

in order to have q(1) = I.

THEOREM 5.1. (Scalar stability).

Let -I < a < 0; then we have for all s with Re s = _ > 0:

1 A

(5.7a) n f (l+x) _(x) ]v(x,s)] 2 dx < const A(N) [v(l,s)] 2, a > 0,

-1

1

(5.7b) _ f (l-x) ,n(x) [v(x,s)] 2^ dx < const A(N) ]v(-l,s)] 2,^ a < 0.

-I

Here, A(N) is a power growth bound:

N-2_ -I < a < 0,

(5.7c) A(N) = I

(N 2 = 0 (Legendre),
-24-

and the corresponding, weighting pair is given by

I

(5.7d) I mII = (l-x) (l+x) _(l-x2)__.

Proof. Let us study first the case a > 0 - without loss of generality

we may take a = I. Consider the Nth degree polynomial

(5.8) = - q(x)$(I,s).

We have:

A

(5.9a) z(xj,s) = v(xj,s), 1 < j < N,

(5.9b) z(l,s)= 0,

and, in view of (5.1):

(5.9c) sz(xj,s) =_xd z(xj,s) + v(l,s) q'(xj) for i < j < N.

Multiply (5.9c) by (l+xj)z (xj,s)_j and sum (terms with j = 0 and j = N + 1

can be included, as their contribution vanishes)
-25-

N d ,s) (5.10) l t+xj)Iz(xj,s)l 2 =N+I _j(t+xj)(_- x Izl2)(xj j=1 s_j( j--J0

N+I ,

+ v(l,s) _ mj(l+xj) z (xj,s) qN(xj). j=0

On the left of this equality we have z(xj,s) = v(xj,s), while the sums on the

right are in fact exact with the corresponding integrals in (5.3) and we

obtain:

(5.11)

N I 1

d 2 v(l,s) f (l+x)_z*

Sj=ly" (l+xj)_j [v(xj,s[2 = -If (l+s)_ _xx [zl dx + ._I q_ dx.

Integrating by parts the first term on the rlght-hand side and using the

Cauchy-Schwartz inequality on the second, we conclude:

(5.12)

N

B J=l[ (l+xj)mj Iv(xj,s)l 2 =

= (l+x)m Izl2 1 - f 1((l+x)_)" Iz12 dx + Re(v(l,s) f (l+x)mz q_ dx)

-I -I -I

I 1 12fl((l+x)_q.)2 --if ((1+x)_)"Izl2 +_if C(1+x)_)I"z2l+ Iv(1,s) -1 [(l+x)_)"dx

= Iv(1,s)l 2 A(N),

where
-26-

(5.13) A(N) = f [(l+x)m(x) ). dx. -I

Note that by adding the boundary terms to the sum on the left, it can be

replaced by the corresponding integral in (5.7), and that all the integrals

involved are proper integrals, because of the limitations on _. It remains

to show that A(N) is bounded by some power of N. A coarse bound can be

obtained immediately, by remarking that

Ixl l

and therefore [q(x) l is bounded by 1 in -I € x _ I. By a well known

extremal property of Chebyshev polynomial, [6 ,Theorem 2.24], the derivatives

cannot grow faster than the corresponding derivatives of TN, and, thus,

lq'l _ N2- Hence,

I

(5.15) A(N) < N 2 F(_), F(a)= f (l+x)a+2 (l-x)a+l

-I 1 + (2_-l)x dx.

The more delicate estimates of (5.7c) are relegated to Appendix B. The proof

of (5.7b) follows along similar lines.

We now show that the collocation methods we have described satisfy

Assumption III.

A

T_EOREM 5.2. Let v be a polynomial in x which satisfies

^

= d v(x,s) + (x) (5.16) sv(x,s) _-_ q
-27-

A

where q(x) is defined by (5.6). If v(x=l,s 0) = O, then Re so < O. l__n_n

^

fact, v satisfies (4.1b), i.e., Assumption III holds.

Proof. In a manner similar to (5.10) - (5.12) we can deduce

^

(5.17) Re •so(l+xj)_ j [v(xj,So)

^

l2 = (l+xj)mj _-_ d Iv12(xj,s0) , 0 < j < N+I,

(5.18)

1 I 1

Re sO f (l+x)m Ivl2 dx = f (l+x)_d iv21dx =- f Ivl2((l+x)m)" dx.

-i -I -I

From the last formula, it is obvious that Re sO < 0. Because of the symmetry

of q(x),

(5.19) q(x) m qN(x) = (-I) N qN(-X),

A

it is clear that $(-l,s) = (-I) N v(l,-s). Consider now the rational

function of s, with real coefficients:

^

(5.20) f(s) - _(l,-s) •

v(l,s)

The limit of If(s) I exists uniformly as s . = and equals one, while on the

imaginary axis

(5.21) If(s)[ = I,
-28-

being the quotient of a number and its conjugate. Therefore, since

f(s) = (-I) N _(-l,s) is regular for Re s > 0, it is bounded by one in

v(1,s)

magnitude - which is Assumption III.

With minor modifications the same proof covers the case a < 0.

Finally, we can summarize this section by applying Theorem 4.3 to the

scalar results we have obtained that far.

THEOREM 5.3. (Stability of Gegenbauer Collocation for Systems).

Consider the hyperbolic system (I.i), satisfying Assumption I. Discretize

2

it by collocating at the zeros of CN , where -I < _ < 0. The resulting

spectral method is stable, in the sense of Lemma 1.2, with parameters defined

by (5.7c) and (5.7d). In particular, Chebyshev and Legendre collocation

methods are is stable.
-29-

APPENDIX A

In this section we prove several a priori bounds for the solution of the

differential system (I.I) and its scalar counterparts.

Introducing into (1.1) the new variable e-qt _(x,t) and taking the

Fourier transform with respect to time, we obtain the equations

du I

(A. la) s ^Iu = AI ~

N dx

dill

^II All ~

(A. Ib) s _ = dx ' -I < x < i,

with boundary conditions

_l(-l,s) = L _ll(-1,s) + _l(s)

(A.2)

^I

_ll(l,s) = R _ (l,s) + _ll(s).

Set ml(x) = 1 + cx with 0 < g < I yet to be determined, and multiply

I ^I*

(A.la) by m _ . After integration by parts, we find:

. : 1 r1 ul, A I uI

x=-I -i ~ ~

where q = Re s. Since sA I is negative, the second term on the right-hand

side does not exceed 1 - s IAI llul(x,s)ll21 ~ and hence

^I 21 ~ I1

(A.3) (q - qO ) I1_ (x,s)ll < (l+sx)ul*(x,s) AI _l(x,s)

x=-i
-30-

with n O =_ _ IAI.

We can draw two conclusions from (A.3). First, by applying it to each of the

scalar inflow equations separately, we have for all s with Re s = n > _, a

positive left-hand side

^(") 2

0 < aj[(l+_) lu(J)(1,s)l 2 - (1-c) lu 3 (_1,s) I ], 1 < j ¢ £,

that is, as aj < 0

u(1)(x's) )]

(A.4) U I (l,s) _ 1 - _ , U I = • .

uI(_1,s)I+c •u(O(x,s.

Secondly, since (i + sx)A I is negative, the inequality (A.3) yields:

(A.5) (_ -_0 ) ,;_l(x,s),j21< (I - _) IAII lul(-l,s)l 2.

We can treat similarly the outflow part of the system. Choosing

II

= 1 - _x, we find that

[A 1

u(_+1)(x,s)

ulI(-l.s) _ i - g UII • •

_il . ,

(A.6) ills) I + _ ' = u(n)(x,s

and
-31-

(A.7) (B-n0),'un(x,s),2_i < (I-OIA_Il_II(1,s)l 2 to

Finally, we rewrite the boundary conditions in the form

( )c ic

u (-1,s) g (s)

I L ull(l,s) ~

= , Re s = n > _.

R ' I

UI(-I,s)

As in Lemma 4, the matrix on the left is invertible; in view of (a.4) and

(A.6) the inverse is bounded in the norm

I I -L UII(-l's)i I

, I - c

UII(I s) 1 + ([R I + ILl) 1 +

< (I-_2 '

-R UI(I's) 1- IRI ILl _1+€_

UI(_I,s) "

provided _ is chosen so that the denominator is positive. Then we have:

(A.8) I_I(-1,s)12+ I_I_(1,s)l 2 € eonst[l_I(s)12+

and together with (A.5) and (A.7) we reach

(A.9) (n - _0 ) II (x,s)ll2 < const l_(s)l A 2, Re s = n > no •

Hence, the energy estimate referred to in (l.2b) follows from Parseval's

relation, with any _0 satisfying
-32-

(A.10) nO - 1 s- _ JA[ _ (]R] ]L]) 1/2 - 1 [A["

2(]R[ JL[) I/2

We conclude by noting that in the dissipative case, where ]R] [L[< I,

one may take q0 = _ = 0 and the formulas (A.4),(A.6), then lead to

[u(J)(1,s)] _ ]u(J)(-1,s)], i < j <

(A.II)

[u(J)(-l,s)[ ( [u(J)(l,s)], £ < j < n,

for Re s = n > 0.
-33-

APPENDIX B

Here we evaluate (5.13), in order to obtain the asymptotic behavior of our

estimates for large N. We shall compute

11(l+x)_q_) 2 1 [(l+x)(l-x2) = _C_N+I+I/2(x))"]2

(B.I) A(N) =-I[ C(1+x)_)" dx =-I[ [(l+x)(l-x2)_]" ['-N+I (fi_+i/2),,(1)]2 dx,

separately for = = 0, the Legendre spectral method, and for -I < _ < 0. In

the first case, we integrate-by-parts

(B.2)

" - [ PN+IIP_+I(I+x)2)" dx. -i[ (l+x) 2 (e_+l)2 dx = eN+l(X) eN+l(X)(l+x2) I -I

Once we obtain an integral involving PN+I' we can use Lobatto quadrature to

exploit the fact that PN+I vanishes at the interior nodes. Then use the

explicit estimates

dm

C_xm PN)(1) = 0(N 2m)

and

_0 = 0(N-2)

to deduce that A(N) = 0(N2).

In the second case we proceed by first estimating
-34-

11 i,I i:i 2

(B.3) A(N) = (l-x)_+2 (l-x) =+I -N+I )"(x

-i I + (2c_-l)x fi_+ I/2 -- dx

-N+I 1 )J

)(x)1

4 I/2 /i (l+x) a+2 (l-x) a r a+ i/2 2 dx.

min[l_,, 1- ,=1)-1 [ICN_: 1/2 )_(1)J

We then substitute from the Gegenbauer differential equation

(B.4) (1-x2)(_N+ Ca+I1/2 )" = (2a+l) x r-N+I Ca+I/2 )" - (N+I)(N+2a+2)CN: 1/2

to reach

1 [ /11 (l_x2)_(l+x)(2_+l , r^_+1/2 . a+I/2 ),,

)Xk_N+ 1 )"[CN+1 dx

(B.5) A(N) ~ i(c_+ 1/2 12 _-N+I )'(i) -

I j

- (N+I)(N+2_+2) / (l-x2)_(l+x) C_+:_ ICa+l_ )" dx -N+I _-N+I -I

In this expression, the second integral vanishes by orthogonality, and the

first may be evaluated by Lobatto quadrature

(B.6) A(N) ~ 1 " 2m0" "(C_+1/2N+I )"(1)(C_+I N+I/2"")'(I).

3

As (l+x)C N vanishes at the nodes xj, I < j < N+I, we can compute m0:
-35-

3

_+

As (l+x)C N vanishes at the nodes xj, 1 _ j < N+I, we can compute m_:u

3

_+

I (l-x2)a (l+x) CN (x)dx

(B.7) m0 = f 3 '

-I _+ --

2CN 2(I)

which we estimate by the means of (5.14). Finally, we obtain

(B.8) A(N) < 0(N-2e)

in accordance to formula (5.7c).
-36-

References

[I] P. J. Davis and P. Rabinowitz, Methods of Numerical Integration, Academic

Press, New York, 1975.

[2] D. Gottlieb and S. A. Orszag, Numerical Analysis of Spectral Methods:

Theory and Applications, SLAM, Philadelphia, 1977.

[3] R. G. Voigt, D. Gottlieb, and M. Y. Hussaini (eds.), Spectral Methods for

Partial Differential Equations, SLAM, Philadelphia, 1984.

[4] H. O. Kreiss, "Initial boundary value problems for hyperbolic systems,"

Comm. Pure Appl. Math., Vol. 23 (1970), pp. 277-298.

[5] L. Lustman, "The time evolution of spectral discretizations of hyperbolic

systems," SINUM, to appear.

[6] T. J. Rivlin, The Chebyshev Polynomials, Wiley-lnterscience, New York,

1974.

[7] D. Gottlieb and E. Turkel, "Spectral methods for time-dependent partial

differential equations," NASA Contractor Report 172241, October 1983.

[8] L.G.M. Reyna, Ph.D. Thesis, Part I, California Institute of Technology,

1982.

Standard Bibliographic Page

1. Report No. NASA CR-l 7804 J. 12. Government Accession No. 3. Recipient's Catalog No.

Tr.ASP. 'Rpnnrt Nn Rn-?

4. Title and Subtitle 5. Report Date

STABILITY ANALYSIS OF SPECTRAL METHODS January 1986

FOR HYPERBOLIC INITIAL-BOUNDARY VALUE SYSTEMS 6. Performing Organization Code

7. Author(s) 8'l6;!~rmingOrganization Report No. David Gottlieb, Liviu Lustman, and Eitan Tadmor

9. Performing Organization Narne and Address 10. Work Unit No.

Institute for Computer Applications in Science

and Engineering 11. Contract or Grant No.

Mail Stop 132C, NASA Langley Research Center NASl-17070; NASl-18l07 u~~~.. ~~ UA ?-:te:.e:.,,_,,??c;

12. Spoifsoring'Agency Name and Address 13. Type of Report and Period Covered

National Aeronautics and Space Administration

Washington, D.C. 20546 14. ::Sponsormg ;:tgency-voae

'lO'l-':n -R1-01

15. Supplementary Notes

Langley Technical Monitor: Submitted to SIAM J. Numer. Anal.

J. C. South

Final Report

16. Abstract

We treat a constant coefficient hyperbolic system in one space variable,

with zero initial data. Dissipative boundary conditions are imposed at the

two points x = ±l. This problem .is discretized by a spectral approximation

in space. We demonstrate sufficient conditions under which the spectral

numerical solution is stable - moreover, these conditions have to be checked

only for scalar equations. The stability theorems take the form of explicit

bounds for the norm of the solution in terms of the boundary data. The

dependence of these bounds on N, the number of points in the domain (or

equivalently the degree of the polynomials involved), is investigated for a

class of standard spectral methods, including Chebyshev and Legendre

collocations.

17. Key Words (Suggested by Authors(s)) 18. Distribution Statement

hyperbolic systems, 64 - Numerical Analysis

boundary conditions,

spectral methods,

stability Unclassified - unlimited

19. Security Classif.(of this report) 120. Security Classif.(of this page) 21. No. of Pages 122. Price

Unclassified Unclassified 38 A03

For sale by the National Technical Information Service, Springfield, Virginia 22161

NASA-Langley, 1986