Asymptotic Analysis and the Numerical Solution of Partial Differential Equations

Asymptotic Analysis and the Numerical Solution of Partial Differential Equations

ByHans G. Kaper, Marc Garbey

Edition 1st Edition

First Published 1991

eBook Published 25 February 1991

Pub. location Boca Raton

Imprint CRC Press

Pages 286 pages

eBook ISBN 9780429081996

SubjectsMathematics & Statistics

Export Citation

#### Get Citation

Kaper, H., Garbey, M. (1991). Asymptotic Analysis and the Numerical Solution of Partial Differential Equations. Boca Raton: CRC Press, https://doi.org/10.1201/b16933

Integrates two fields generally held to be incompatible, if not downright antithetical, in 16 lectures from a February 1990 workshop at the Argonne National Laboratory, Illinois. The topics, of interest to industrial and applied mathematicians, analysts, and computer scientists, include singular per

TABLE OF CONTENTS

Robert E. O'Malley, Jr.

Bythe Connection Formulas for the Transition

Uo(t) (e;) -+ g( u)

By< e< < < B, the monotonic solution U(t) will < t < h < t < + + = that all the time will be spent at e; if j the other to the literature (e.g., Fife [6]) to find applications to reaction- and Vasil'eva's study [2] of The conservation equation

(t) according to the implicit function theorem, pro-

By(50) the singular solution (51) and the implicit solutions (52) be uniquely solved for U = r.

of contrasting spatial of Internal Layers and Diffusive Interfaces, SIAM, Philadelphia,

ByJ. Math. F. Butuzov and A. B. Vasil'eva, The asymptotic theory 28, no. 2 (1988), 26-36. and C. E. Pearson, Ordinary Differential Equations, Blaisdell, Waltham,

ByR. C. Y. Chin, G. W. Hedstrom Garry Rodrigue, Edna Reiter

Xn(x; x )/f.} )(1 x 1. Thus, up to terms of order m Xn(O; fn, Xn(1;

Byand sn = that in the subcase (15) we have + 0( the coefficient Xn is the solution of the problem the purely convected solution

Byu(1) (32) the coefficients a and b and the source term and solves the resulting equation exactly. Here, we adapt this method to the adaptations to be made when turning points are

Bythe recurrence coefficients an and bn for n 0. We note that the integrands of ( 45-46) are essentially polynomials of order 2n and 2n+ 1 and, therefore, their numerical 44) and ( 46) with (44) into (43) gives q (45) to that this is a forward marching process. It is well known that the calculation of numerical the delicate nature of the evaluation of highly oscillatory the process of numerically generating orthogonal polynomials-given just (48) ( 49) parts the integrals

ByP. A. Lagerstrom, L. N. Howard and C.-S. Liu, eds., Academic Press, New York, 1967. J. Kevorkian and Verlag, New York, 1981. P. A. Lagerstrom and R. G. Casten, Basic concepts underlying singular perturbation

g(y), 0 u, g(y),

Bydata g, Jo, of the type given by (1) frequently arise in fluid dynamics and heat transfer the boundary functions are not continuous at the corners of D. These of the differential equation which, if absent, would change the character < y < That is, on most of the outflow boundary x the boundaries y thin layers of fluid near solid boundaries, i.e., boundary layers. For a fundamental discussion of the the classical

fo(x) g(y)

Bywith£= 5 x 10-was solved by using the algorithms given by (29)and (31). The boundary conditions for the problem were taken to be ft(x) the values = 15 x 10-4, we computed = 10-to be Then, taking £= X£ and £= of k, we carry out algorithm (31) using

lo [0(x, h(r)dr lo 8y 8rk rx !!_0(x -r,y) fo(r)dr. lo 8y 8rk lo By[0(x,

By10-5 10- the collaboration of Ray Chin and Gerry Hedstrom during this research Proof of Theorem 2. We first show that - (!1), n;::: 0, and is given by 0(x, +2 rx !!_0(x-r,y-c-lf 0(x, y

{1/,.fo s(x,y)= lo {1/,.fo [)2 k C (n) so that fJxk C (n). More-

By(A.6) parts twice, and using (A.3), the same arguments as before for computing ::,we that the proof of (A.4) is complete. fork- that &x

txF(U) (x,t) E !1

ByIRx)O,T[, that V is piecewise smooth. We also assume that F and P are smooth that Pis suitable viscosity matrix the the

S(t))/E t S1t)jE to)/f f locally linear S1t)jE to)/f3f4 lxf(u) lx ( G (!, u, D..t,D..x) u(x, tk) given.

Byof Zone Order of Residual Local Coordinates 0( f1f2) S(t))/E1 0(1) of shock 0( f 1/4) the following Cauchy problem: (12)

ByFigure 7: Viscous numeric compared with inviscid exact, t .00---

Thomas Hagstrom

Byin Numerical Computations to Their Numerical Simulation Partial Differential Equations

In C is positive definite. Then there exists 8 of (5) are real when s is real and s > 8. iv and, by an integration by parts, obtain sIn Mv vHTv B(v,A)= vHDoov)+vHcv)

ByStokes equations linearized about a parallel flow. In fact, s* 0 may be important in our asymptotic analysis of linearized viscous flow both for very low and very high Reynolds numbers now turn to the analysis of (5). In what follows we shall drop the tildes and assume that (1-4) have been given in the partially symmetrized form above. Proof. Following [3] we write A (29) } , (30) (31)

of the ex- tinction limit of curved flames, Combust. Sci. and Tech. 64 (1989), 187-198. of solu- tions of semilinear elliptic equations, preprint. of solutions of semilinear

ByF. Benkhaldoun, B. Larrouturou, and B. Denet, Numerical investigation H. Berestycki and L. Nirenberg, Monotonicity, symmetry and antisymmetry Berestycki and L. Nirenberg, Some qualitative properties

{(x, y) : 0

ByTo illustrate the basic features of the hybrid method, we consider the following simple two-dimensional example. define u(x,y,£) as the solution to the problem £Sin(x) cos(y) 2 Esin(x) sin(y) 0, (x, y) E V, (5) with u 0 for (x, y) E 8V. Here V < x, < and 8V denotes the boundary of V. In (5), the symbol Vdenotes the usual two-dimensional Laplacian operator, and the subscript denotes differentiation with respect to x. The solution seems to be positive over the whole domain V for positive values of f. It is invariant under the 180-degree rotation

Byof the simple PDE problem of Section 3 as computed numerically (dots) and by the hy- brid method (solid lines) using from N = toN = terms. = the radius of convergence, R, of the series expansion. that the hybrid method might be a useful tool for the of determining an appropriate decomposition of the domain for a purely

of the limit cycle of the van der Pol equation, SIAM J. Appl. Math. 42 (1982),

Bysolution. Obviously, this difficult problem needs much more investigation, e.g., either development of a numerical method that will yield the steady solution for in the range C. M. Andersen and J. F. Geer, Power series expansions for the frequency and period 678-693. C. M. Andersen and J. F. Geer, Investigating a hybrid perturbation Galerkin technique NASA Langley, !CASE Report No. 88-65, 1988.

A. Louise Perkins Pamela Cook

ByOn the Equations of Physical Oceanography Transonics and Asymptotics Evolution to Detonation in a Nonuniformly Heated Reactive Medium

dtu+ aV'P+ F P1Y' · u a is the specific volume, P is the pressure, and 1 is the ratio of specific u · V' is the total time

ByCoupling hydrodynamics and thermodynamics, consider an adiabatic, inviscid fluid. It can be described by conservation of momentum, a continuity equation, and an energy equation coupled with an equation of state. That is, where u E R heats. Here F is the Coriolis and gravity forces, and derivative.

ByThe initial conditions for a physical model must be realistic and must not violate any assumptions that the model equations encompass. For primitive equation ocean model- ing, the initial conditions must be near geostrophic balance to avoid unrealistically large gravity waves. This has been achieved by spinning up a model to statistical equilibrium, and then applying forcing until the desired initial conditions are achieved. Alternatively, a

S(x,y,z) Sx(K</Jx- -<Px]s

Bybut for the flow is basically subsonic (elliptic) with an embedded supersonic (hyperbolic) the hyperbolic region is not known a priori and must part of the solution. In most cases the flow cannot compress smoothly from to subsonic; thus, a shock (discontinuity in velocity) occurs. the shock surface, with normal S}, then across shocks we have 1'+12

Bythe nonlinear sonic equation (26) with K the variational equation. Thus -(1' ( -(1' in order to satisfy the boundary conditions on the wing < x < 1,

involving lift,

ByThe author acknowledges the contributions of J. D. Cole to this work. H. K. Cheng and M. M. Hafez, Transonic equivalence rule: A nonlinear problem J. Fluid. Mech. 72, no. 1 {1975), 161-187. J.D. Cole, Modern developments in transonic flow, SIAM J. Appl. Math. 29, no. 4

=Yo, u -Ao. t) = u(Lo, t) = 0.

Byto be 70(0) =To, 4{(0) = The rigid-wall boundary conditions are the reference state. Time unit is chosen to

1P"'

ByBefore leaving this section, we point out that the induction equations have been the subject of a number of studies, mostly numerical, for a variety of initial and boundary data. Initiation of chemical reaction by a shock, and the resulting thermal runaway, have been examined by Clarke and Cant [4] and Jackson and Kapila [8] for a piston-driven shock and by Cant [2] and Short and Dold [12] for a shock driven by a column of compressed and as a result, the induction phase is dominated by constant-volume heating. The wave processes play a secondary role, managing to influence the solution only at order A. The asymptotic result is -ln(1-1t) O(A), U'"" O(A). (35)

of the interaction between a shock wave

ByP. A. Blythe and D. G. Crighton, Shock-generated ignition: the induction zone, preprint, 1989. R. S. Cant, Ph.D. thesis, Cranfield Institute of Technology, 1984. J. F. Clarke, On the theoretical modelling and an explosive gas mixture, College of Aeronautics Report 7801, Cranfield Institute

Byat different times on the track of the triple points found on the to the radius of curvature of the shock in that region. Thus, if the characteristic 1-D reaction zone length and is the local shock curvature, weak If the of and a quasi-steady detonation flow, then it is possible to derive an equation(s) the evolution of the smooth part of the shock surface in space and time. The asymp- to the local shock The functional relationship between Dn and depends on the detailed

A plane, starts

Bylayers, where the asymptotic expansions for u in terms of"' have different form, hence the singular character of the perturbation problem. The largest layer in the u near the shock and ends near the equilibrium point. Klein and Stewart [5] call this the "main reaction layer"; it is a regular perturbation of the 1-D solution. The second layer is called the "transonic layer"; it lies close to the equilibrium point where perturbations Here we cite some research issues that we know are being pursued at the time of this writing. First, there is a great deal of interest in broadening the scope of these calculations to more realistic models for the equations of state and different reaction rate laws that might be imposed by more complex chemistry or heat release mechanisms. Some of the general features of the analysis, represented in this account, may not survive such a

interaction of edge rarefactions with finite reaction zones, J. Fluid Mechanics 171 of two-dimensional detonation with detona- tion shock dynamics,

ByJ. B. Bdzil and D. S. Stewart, Two dimensional, time dependent detonation: The (1986) 1-26. J. B. Bdzil and D. S. Stewart, Modeling Physics of Fluids A, 1, no. 7 (1988) 1261.

t) Oo[Vb]wo(p, v, t)dp. I/:). Vb(p, TJ)}( t) exp( iTJv)dTJdp t) denotes the Fourier transform of the function f( v, t) in (28) in the v direction. p -TJ and 0 I v,

Bythe definition of -oo, v, t) where}( and 0 that the result.

m ).

Byon a two-dimensional toroidal lattice, in which each node has four directions. The the LB1 are described as follows: if a state has exactly two particles, the new state (prior to advection) is determined by sending the particles into the the kth element of the collision operator, C(), for LB )(b) the sequel), the k indices are evaluated modulo four. The kth element of the Analysis The Lattice Boltzmann Methods ·The lattice Boltzmann methods that we shall discuss are related to the dynamics of the mean occupation numbers in the associated lattice gas methods. Let ( ·) denote the

u:n.+l vfit [(D(Ul+ (Ul+t,j- UlJ) (D(UT".) [(D(UlJ+ D(U[j)) (U;j+ UlJ) (D(U;j) D(U;j_ (U;j-

Byto compare the lattice Boltzmann-computed solutions, in which D(u) is defined in Eq. (9) and Eq. (10) for LB1 and LB2, respectively: 2(fix)2 + (urn.- 2(fiy) U;'J-t)].

u1(x,y)

ByAsin(2?rX)sin(27ry)+B, where 0.45 and B 1/2; (b) u(x,y;t) at t 1/32, using LB1.

ByThe first case regards Tables 1 and 2. Here, the initial condition satisfies the conditions of monotonici ty for LB 1· We look at the ratio between the norm of the difference between the finite difference- and lattice Boltzmann-computed solutions for grid sizes a factor of two apart, expecting the ratio to be near four. The L-norm results support the theoretical second-order convergence in that norm ofLB1. and even suggest a slight superconvergence. entries in Tables 7 and 8. It should be noted that 128 entries do uphold the ratio of sixteen for errors on grid sizes a factor of four apart, e.g., at t 1/32 the ratio, IIE< 16.042 based on Table 7. The apparent discrepancy may be restricted to smaller grid sizes.

Avner Friedman Gunilla Kreiss

ByPartial Differential Equations Blow-up of Solutions of Nonlinear Heat and Wave Equations Convergence to Steady State of Solutions of Viscous Conservation Laws

T the maximal number such that the solution exists for 0 u(x, t)--+ oo if t--+ T. tn to t T) such that (xm, tm)· Such points (x T) are called blow-up points. (x) has at most two local maxima, see [14, 8]). f 4) that C > 0,

By> 0 independent of n. It follows that the solution exists for 0 < t < T (6) it follows that there must exist points (x 1 the number of blow-up points is finite (Chen and Matano [9]; for earlier that requires that u > 1, does the blow-up set have Lebesgue measure zero? and Kohn [17-19] have studied the behavior of u near a blow-up point ( x + yJT-

0(8) and R(x) R(O) 0(8) uv(x),Rv(x) with an +0(e -ixifv), Rv R+O( e-ixifv) R(x) p'(x, t), u U(x) u'(x, t), _xs=O.

BySince the nozzle is almost straight, we have U(x) on each side of the shock. In the viscous case we assume that there is a steady solution interiorlayer at x 0 instead of a shock and that on each side of the shock layer. in Section 2 we shall linearize around the steady solutions and derive eigenvalue problems. We start with the inviscid case and introduce Again, as in Section 2, we must allow perturbations of the shock position and consider all [ ( p' ) )+ _+ RU2 ) )+ (21) Here are smooth matrices on each side of the shock. Let where is a constant. Then = -1 ¢2(1) (22)

2RU)¢1x=xa jxa ax (U 2RU)¢dx v¢2xlx=xa· (U, R)¢1x=xo = 0((1>.1 8)v), and to first approximation we have (U, R)¢1x=xo J', 2RU)¢1x=xo -lxo

ByL:o dx

of the shock layer there is exactly one solution of (46} lsi+ h:Ao- O(v).

Bytrue for and Do evaluated on either side of the shock layer. The lemma follows. Lemma 2 On the right side > 0. Then 0. Proof. Note that the corresponding constant coefficient half-plane problem, -oo < x' 0, t' 0, v2(0) Dol= 0, Re(s) > 0, > 0, the solutions of (46) must also J. Goodman, and A. Majda, Multiple steady states for 1-D transonic flow, J. Sci. Stat Comp. 5, no. 1 (1984), 21-41.

equation, Appl. Num. Math. 2 (1986), 161-179. Equation, Academic Press, New York, 1989.

ByH.-0. Kreiss and J. Lorenz, Initial-Boundary Value Problems and the Navier-Stokes T. P. Liu, Nonlinear stability and instability Comm. Math. Phys. 83 (1982), 243-260.

f( u) uP+u, showed that the solution p is supercritical). In fact, for balls of certain sizes they were able to show that at least

Bypis greater than or equal top*, where is the celebrated Sobolev (critical) exponent, > 2. the case of the Dirichlet problem for a ball of finite radius is in general not unique if p that, on the other if < p < pis subcritical). the same problem with

J, such that there is at most one positive solution of the differential 8uf 8a satisfies the [q(r)f'(u)]w w'(O)

By(0, b) satisfying the conditions u'(O) 0 and u(b) of u(b) with the initial height a. The partial derivative w (10) that w(b) < 0. Intuitively, w(b) < 0 means that, when the initial height a increases, the solution curve moves down at b, so that the to the left. It follows that only one initial height can produce a zero at

Jt f( (p, q) satisfy these require-

Byand elegant. Their to concrete examples is another story. The power of a criterion is measured of problems it can cover, but its usefulness depends on how easily it can be applied. In practice, one keeps a repertoire of test problems, usually not yet to benchmark any new criterion. In many instances MAPLE the equation into an equivalent but much that, if the coefficient of the linear term of the transformed equation is nondecreasing in [0, c] and nonincreasing in the Dirichlet problem for (13) has at most one solution. the number given by The exact domain of admissible pairs has no simple explicit characterization; it is the conditions to some nonlinear inequalities, which can then be evaluated numerically. to find out whether uniqueness is still valid for the values of p and q not covered

ByBackground Our second case study comes from the area of low-Mach number combustion. of chemical reactions transforming fuel the reactions require a large activation energy, so the the temperature of the fuel is sufficiently high. That

AxRR2 R2-d

Bythe case where is a simple eigenvalue, strictly less than all other eigenvalues. that this can be accomplished by a judicious choice of p.) that, in general, the dynamical system lives on the -scale. i¢jh .... ) 1 2 (43) the evolution of the amplitudes is governed by the system of equations (44) + + + the phase functions satisfy (46) (47)

-Ax. The system goes

Bythat interesting phenomena occur as Ac or Ac to unstable subcritical (Fig. 2(e)), or a standing mode from unstable subcritical to stable supercritical (Fig. 2( c)). These transitions are observable, as they the analysis was done with

ABOUT THIS BOOK

CONTENTS