c19-0 (779615), страница 2
Текст из файла (страница 2)
The matrix A is called “tridiagonalwith fringes.” A general linear second-order elliptic equationa(x, y)∂2u∂u∂u∂2 u+c(x,y)+b(x,y)+ d(x, y)∂x2∂x∂y2∂y2∂ u+ f(x, y)u = g(x, y)+ e(x, y)∂x∂y(19.0.11)will lead to a matrix of similar structure except that the nonzero entries will notbe constants.As a rough classification, there are three different approaches to the solutionof equation (19.0.10), not all applicable in all cases: relaxation methods, “rapid”methods (e.g., Fourier methods), and direct matrix methods.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).∆832Chapter 19.Partial Differential Equationsincreasing iincreasing j−4111 −41•••eachblock(L + 1) ×(L + 1)1••••1 −4•1•1 −411111 −41••••••••••••••1•••••••••••••••••••••••••••••J+1blocks••••••••••••••••••••••1 −41•111 −411−41••11 −41•••••••1 −4•111 −4J + 1 blocksFigure 19.0.3.
Matrix structure derived from a second-order elliptic equation (here equation 19.0.6). Allelements not shown are zero. The matrix has diagonal blocks that are themselves tridiagonal, and suband super-diagonal blocks that are diagonal. This form is called “tridiagonal with fringes.” A matrix thissparse would never be stored in its full form as shown here.Relaxation methods make immediate use of the structure of the sparse matrixA.
The matrix is split into two partsA=E−F(19.0.12)where E is easily invertible and F is the remainder. Then (19.0.10) becomesE·u= F·u+b(19.0.13)The relaxation method involves choosing an initial guess u(0) and then solvingsuccessively for iterates u(r) fromE · u(r) = F · u(r−1) + b(19.0.14)Since E is chosen to be easily invertible, each iteration is fast. We will discussrelaxation methods in some detail in §19.5 and §19.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).−4119.0 Introduction833A·x=b(19.0.15)directly. The degree to which this is practical depends very strongly on the exactstructure of the matrix A for the problem at hand, so our discussion can go no fartherthan a few remarks and references at this point.Sparseness of the matrix must be the guiding force. Otherwise the matrixproblem is prohibitively large. For example, the simplest problem on a 100 × 100spatial grid would involve 10000 unknown uj,l ’s, implying a 10000 × 10000 matrixA, containing 108 elements!As we discussed at the end of §2.7, if A is symmetric and positive definite(as it usually is in elliptic problems), the conjugate-gradient algorithm can beused.
In practice, rounding error often spoils the effectiveness of the conjugategradient algorithm for solving finite-difference equations. However, it is usefulwhen incorporated in methods that first rewrite the equations so that A is transformedto a matrix A0 that is close to the identity matrix. The quadratic surface defined by theequations then has almost spherical contours, and the conjugate gradient algorithmworks very well. In §2.7, in the routine linbcg, an analogous preconditionerwas exploited for non-positive definite problems with the more general biconjugategradient method. For the positive definite case that arises in PDEs, an example ofa successful implementation is the incomplete Cholesky conjugate gradient method(ICCG) (see [6-8]).Another method that relies on a transformation approach is the strongly implicitprocedure of Stone [9]. A program called SIPSOL that implements this routine hasbeen published [10].A third class of matrix methods is the Analyze-Factorize-Operate approach asdescribed in §2.7.Generally speaking, when you have the storage available to implement thesemethods — not nearly as much as the 108 above, but usually much more than isrequired by relaxation methods — then you should consider doing so.
Only multigridrelaxation methods (§19.6) are competitive with the best matrix methods. For gridslarger than, say, 300 × 300, however, it is generally found that only relaxationmethods, or “rapid” methods when they are applicable, are possible.There Is More to Life than Finite DifferencingBesides finite differencing, there are other methods for solving PDEs. Mostimportant are finite element, Monte Carlo, spectral, and variational methods. Unfortunately, we shall barely be able to do justice to finite differencing in this chapter,and so shall not be able to discuss these other methods in this book.
Finite elementmethods [11-12] are often preferred by practitioners in solid mechanics and structuralSample 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).So-called rapid methods [5] apply for only a rather special class of equations:those with constant coefficients, or, more generally, those that are separable in thechosen coordinates.
In addition, the boundaries must coincide with coordinate lines.This special class of equations is met quite often in practice. We defer detaileddiscussion to §19.4. Note, however, that the multigrid relaxation methods discussedin §19.6 can be faster than “rapid” methods.Matrix methods attempt to solve the equation834Chapter 19.Partial Differential Equationsengineering; these methods allow considerable freedom in putting computationalelements where you want them, important when dealing with highly irregular geometries. Spectral methods [13-15] are preferred for very regular geometries and smoothfunctions; they converge more rapidly than finite-difference methods (cf. §19.4), butthey do not work well for problems with discontinuities.Ames, W.F.
1977, Numerical Methods for Partial Differential Equations, 2nd ed. (New York:Academic Press). [1]Richtmyer, R.D., and Morton, K.W. 1967, Difference Methods for Initial Value Problems, 2nd ed.(New York: Wiley-Interscience). [2]Roache, P.J. 1976, Computational Fluid Dynamics (Albuquerque: Hermosa). [3]Mitchell, A.R., and Griffiths, D.F. 1980, The Finite Difference Method in Partial Differential Equations (New York: Wiley) [includes discussion of finite element methods]. [4]Dorr, F.W.
1970, SIAM Review, vol. 12, pp. 248–263. [5]Meijerink, J.A., and van der Vorst, H.A. 1977, Mathematics of Computation, vol. 31, pp. 148–162. [6]van der Vorst, H.A. 1981, Journal of Computational Physics, vol. 44, pp. 1–19 [review of sparseiterative methods]. [7]Kershaw, D.S. 1970, Journal of Computational Physics, vol.
26, pp. 43–65. [8]Stone, H.J. 1968, SIAM Journal on Numerical Analysis, vol. 5, pp. 530–558. [9]Jesshope, C.R. 1979, Computer Physics Communications, vol. 17, pp. 383–391. [10]Strang, G., and Fix, G. 1973, An Analysis of the Finite Element Method (Englewood Cliffs, NJ:Prentice-Hall). [11]Burnett, D.S.
1987, Finite Element Analysis: From Concepts to Applications (Reading, MA:Addison-Wesley). [12]Gottlieb, D. and Orszag, S.A. 1977, Numerical Analysis of Spectral Methods: Theory and Applications (Philadelphia: S.I.A.M.). [13]Canuto, C., Hussaini, M.Y., Quarteroni, A., and Zang, T.A. 1988, Spectral Methods in FluidDynamics (New York: Springer-Verlag). [14]Boyd, J.P.
1989, Chebyshev and Fourier Spectral Methods (New York: Springer-Verlag). [15]19.1 Flux-Conservative Initial Value ProblemsA large class of initial value (time-evolution) PDEs in one space dimension canbe cast into the form of a flux-conservative equation,∂F(u)∂u=−∂t∂x(19.1.1)where u and F are vectors, and where (in some cases) F may depend not only on ubut also on spatial derivatives of u.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
