c19-1 (779616), страница 4
Текст из файла (страница 4)
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).sv∆tξ = −isin k∆x ±∆x844Chapter 19.Partial Differential Equationsin Figure 19.1.6.
This mesh drifting instability is cured by coupling the two meshesthrough a numerical viscosity term, e.g., adding to the right side of (19.1.31) a smallcoefficient ( 1) times unj+1 − 2unj + unj−1. For more on stabilizing differenceschemes by adding numerical dissipation, see, e.g., [4].The Two-Step Lax-Wendroff scheme is a second-order in time method thatavoids large numerical dissipation and mesh drifting. One defines intermediatevalues uj+1/2 at the half timesteps tn+1/2 and the half mesh points xj+1/2 . Theseare calculated by the Lax scheme:n+1/2uj+1/2 =1 n∆t(u(F n − Fjn )+ unj ) −2 j+12∆x j+1(19.1.37)n+1/2Using these variables, one calculates the fluxes Fj+1/2 .
Then the updated valuesare calculated by the properly centered expressionun+1j= unj −un+1j∆t n+1/2n+1/2Fj+1/2 − Fj−1/2∆x(19.1.38)n+1/2The provisional values uj+1/2 are now discarded. (See Figure 19.1.7.)Let us investigate the stability of this method for our model advective equation,where F = vu. Substitute (19.1.37) in (19.1.38) to get= unj − αun+1j11 n(uj+1 + unj ) − α(unj+1 − unj )221 n1− (uj + unj−1 ) + α(unj − unj−1 )22(19.1.39)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).Figure 19.1.6. Origin of mesh-drift instabilities in a staggered leapfrog scheme. If the mesh pointsare imagined to lie in the squares of a chess board, then white squares couple to themselves, black tothemselves, but there is no coupling between white and black. The fix is to introduce a small diffusivemesh-coupling piece.19.1 Flux-Conservative Initial Value Problems845two-step Lax Wendroffhalfstep pointsx or jFigure 19.1.7.
Representation of the two-step Lax-Wendroff differencing scheme. Two halfstep points(⊗) are calculated by the Lax method. These, plus one of the original points, produce the new point viastaggered leapfrog. Halfstep points are used only temporarily and do not require storage allocation on thegrid. This scheme is second-order accurate in both space and time.whereα≡v∆t∆x(19.1.40)Thenξ = 1 − iα sin k∆x − α2 (1 − cos k∆x)(19.1.41)so|ξ|2 = 1 − α2 (1 − α2 )(1 − cos k∆x)2(19.1.42)The stability criterion |ξ|2 ≤ 1 is therefore α2 ≤ 1, or v∆t ≤ ∆x as usual.Incidentally, you should not think that the Courant condition is the only stabilityrequirement that ever turns up in PDEs.
It keeps doing so in our model examplesjust because those examples are so simple in form. The method of analysis is,however, general.Except when α = 1, |ξ|2 < 1 in (19.1.42), so some amplitude damping doesoccur. The effect is relatively small, however, for wavelengths large compared withthe mesh size ∆x. If we expand (19.1.42) for small k∆x, we find|ξ|2 = 1 − α2 (1 − α2 )(k∆x)4+...4(19.1.43)The departure from unity occurs only at fourth order in k. This should be contrastedwith equation (19.1.16) for the Lax method, which shows that|ξ|2 = 1 − (1 − α2 )(k∆x)2 + .
. .for small k∆x.(19.1.44)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited.
To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).t or n846Chapter 19.Partial Differential EquationsFluid Dynamics with ShocksAs we alluded to earlier, the treatment of fluid dynamics problems with shockshas become a very complicated and very sophisticated subject. All we can attemptto do here is to guide you to some starting points in the literature.There are basically three important general methods for handling shocks.
Theoldest and simplest method, invented by von Neumann and Richtmyer, is to addartificial viscosity to the equations, modeling the way Nature uses real viscosityto smooth discontinuities. A good starting point for trying out this method is thedifferencing scheme in §12.11 of [1]. This scheme is excellent for nearly all problemsin one spatial dimension.The second method combines a high-order differencing scheme that is accuratefor smooth flows with a low order scheme that is very dissipative and can smooththe shocks. Typically, various upwind differencing schemes are combined usingweights chosen to zero the low order scheme unless steep gradients are present, andalso chosen to enforce various “monotonicity” constraints that prevent nonphysicaloscillations from appearing in the numerical solution.
References [2-3,5] are a goodplace to start with these methods.The third, and potentially most powerful method, is Godunov’s approach. Hereone gives up the simple linearization inherent in finite differencing based on Taylorseries and includes the nonlinearity of the equations explicitly. There is an analyticsolution for the evolution of two uniform states of a fluid separated by a discontinuity,the Riemann shock problem. Godunov’s idea was to approximate the fluid by alarge number of cells of uniform states, and piece them together using the Riemannsolution.
There have been many generalizations of Godunov’s approach, of whichthe most powerful is probably the PPM method [6].Readable reviews of all these methods, discussing the difficulties arising whenone-dimensional methods are generalized to multidimensions, are given in [7-9] .CITED REFERENCES AND FURTHER READING:Ames, W.F. 1977, Numerical Methods for Partial Differential Equations, 2nd ed. (New York:Academic Press), Chapter 4.Richtmyer, R.D., and Morton, K.W. 1967, Difference Methods for Initial Value Problems, 2nd ed.(New York: Wiley-Interscience). [1]Centrella, J., and Wilson, J.R.
1984, Astrophysical Journal Supplement, vol. 54, pp. 229–249,Appendix B. [2]Hawley, J.F., Smarr, L.L., and Wilson, J.R. 1984, Astrophysical Journal Supplement, vol. 55,pp. 211–246, §2c. [3]Kreiss, H.-O. 1978, Numerical Methods for Solving Time-Dependent Problems for Partial Differential Equations (Montreal: University of Montreal Press), pp. 66ff. [4]Harten, A., Lax, P.D., and Van Leer, B. 1983, SIAM Review, vol.
25, pp. 36–61. [5]Woodward, P., and Colella, P. 1984, Journal of Computational Physics, vol. 54, pp. 174–201. [6]Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).In summary, our recommendation for initial value problems that can be castin flux-conservative form, and especially problems related to the wave equation,is to use the staggered leapfrog method when possible.
We have personally hadbetter success with it than with the Two-Step Lax-Wendroff method. For problemssensitive to transport errors, upwind differencing or one of its refinements shouldbe considered.19.2 Diffusive Initial Value Problems847Roache, P.J. 1976, Computational Fluid Dynamics (Albuquerque: Hermosa). [7]Woodward, P., and Colella, P. 1984, Journal of Computational Physics, vol. 54, pp. 115–173. [8]Rizzi, A., and Engquist, B.
1987, Journal of Computational Physics, vol. 72, pp. 1–69. [9]Recall the model parabolic equation, the diffusion equation in one spacedimension,∂∂u∂u=D(19.2.1)∂t∂x∂xwhere D is the diffusion coefficient. Actually, this equation is a flux-conservativeequation of the form considered in the previous section, withF = −D∂u∂x(19.2.2)the flux in the x-direction. We will assume D ≥ 0, otherwise equation (19.2.1) hasphysically unstable solutions: A small disturbance evolves to become more and moreconcentrated instead of dispersing. (Don’t make the mistake of trying to find a stabledifferencing scheme for a problem whose underlying PDEs are themselves unstable!)Even though (19.2.1) is of the form already considered, it is useful to considerit as a model in its own right.
The particular form of flux (19.2.2), and its directgeneralizations, occur quite frequently in practice. Moreover, we have already seenthat numerical viscosity and artificial viscosity can introduce diffusive pieces likethe right-hand side of (19.2.1) in many other situations.Consider first the case when D is a constant. Then the equation∂u∂2u=D 2∂t∂x(19.2.3)can be differenced in the obvious way: n− unjun+1uj+1 − 2unj + unj−1j=D∆t(∆x)2(19.2.4)This is the FTCS scheme again, except that it is a second derivative that has beendifferenced on the right-hand side. But this makes a world of difference! TheFTCS scheme was unstable for the hyperbolic equation; however, a quick calculationshows that the amplification factor for equation (19.2.4) isk∆x4D∆t2(19.2.5)sinξ =1−(∆x)22The requirement |ξ| ≤ 1 leads to the stability criterion2D∆t≤1(∆x)2(19.2.6)Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5)Copyright (C) 1988-1992 by Cambridge University Press.Programs Copyright (C) 1988-1992 by Numerical Recipes Software.Permission is granted for internet users to make one paper copy for their own personal use.
Further reproduction, or any copying of machinereadable files (including this one) to any servercomputer, is strictly prohibited. To order Numerical Recipes books,diskettes, or CDROMsvisit website http://www.nr.com or call 1-800-872-7423 (North America only),or send email to trade@cup.cam.ac.uk (outside North America).19.2 Diffusive Initial Value Problems.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
