c16-3 (779595)
Текст из файла
722Chapter 16.Integration of Ordinary Differential Equations}nrerror("Too many steps in routine odeint");}CITED REFERENCES AND FURTHER READING:Gear, C.W. 1971, Numerical Initial Value Problems in Ordinary Differential Equations (EnglewoodCliffs, NJ: Prentice-Hall).
[1]Cash, J.R., and Karp, A.H. 1990, ACM Transactions on Mathematical Software, vol. 16, pp. 201–222. [2]Shampine, L.F., and Watts, H.A. 1977, in Mathematical Software III, J.R. Rice, ed. (New York:Academic Press), pp. 257–275; 1979, Applied Mathematics and Computation, vol. 5,pp. 93–121.Forsythe, G.E., Malcolm, M.A., and Moler, C.B.
1977, Computer Methods for MathematicalComputations (Englewood Cliffs, NJ: Prentice-Hall).16.3 Modified Midpoint MethodThis section discusses the modified midpoint method, which advances a vectorof dependent variables y(x) from a point x to a point x + H by a sequence of nsubsteps each of size h,h = H/n(16.3.1)In principle, one could use the modified midpoint method in its own right as an ODEintegrator. In practice, the method finds its most important application as a part ofthe more powerful Bulirsch-Stoer technique, treated in §16.4.
You can thereforeconsider this section as a preamble to §16.4.The number of right-hand side evaluations required by the modified midpointmethod is n + 1. The formulas for the method arez0 ≡ y(x)z1 = z0 + hf(x, z0 )zm+1 = zm−1 + 2hf(x + mh, zm )y(x + H) ≈ yn ≡1[zn + zn−1 + hf(x + H, zn )]2for m = 1, 2, . . . , n − 1(16.3.2)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).(*rkqs)(y,dydx,nvar,&x,h,eps,yscal,&hdid,&hnext,derivs);if (hdid == h) ++(*nok); else ++(*nbad);if ((x-x2)*(x2-x1) >= 0.0) {Are we done?for (i=1;i<=nvar;i++) ystart[i]=y[i];if (kmax) {xp[++kount]=x;Save final step.for (i=1;i<=nvar;i++) yp[i][kount]=y[i];}free_vector(dydx,1,nvar);free_vector(y,1,nvar);free_vector(yscal,1,nvar);return;Normal exit.}if (fabs(hnext) <= hmin) nrerror("Step size too small in odeint");h=hnext;16.3 Modified Midpoint Method723yn − y(x + H) =∞Xαi h2i(16.3.3)i=1where H is held constant, but h changes by varying n in (16.3.1).
The importanceof this even power series is that, if we play our usual tricks of combining steps toknock out higher-order error terms, we can gain two orders at a time!For example, suppose n is even, and let yn/2 denote the result of applying(16.3.1) and (16.3.2) with half as many steps, n → n/2. Then the estimatey(x + H) ≈4yn − yn/23(16.3.4)is fourth-order accurate, the same as fourth-order Runge-Kutta, but requires onlyabout 1.5 derivative evaluations per step h instead of Runge-Kutta’s 4 evaluations.Don’t be too anxious to implement (16.3.4), since we will soon do even better.Now would be a good time to look back at the routine qsimp in §4.2, andespecially to compare equation (4.2.4) with equation (16.3.4) above.
You will seethat the transition in Chapter 4 to the idea of Richardson extrapolation, as embodiedin Romberg integration of §4.3, is exactly analogous to the transition in going fromthis section to the next one.Here is the routine that implements the modified midpoint method, which willbe used below.#include "nrutil.h"void mmid(float y[], float dydx[], int nvar, float xs, float htot, int nstep,float yout[], void (*derivs)(float, float[], float[]))Modified midpoint step.
At xs, input the dependent variable vector y[1..nvar] and its derivative vector dydx[1..nvar]. Also input is htot, the total step to be made, and nstep, thenumber of substeps to be used. The output is returned as yout[1..nvar], which need notbe a distinct array from y; if it is distinct, however, then y and dydx are returned undamaged.{int n,i;float x,swap,h2,h,*ym,*yn;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).Here the z’s are intermediate approximations which march along in steps of h, whileyn is the final approximation to y(x + H). The method is basically a “centereddifference” or “midpoint” method (compare equation 16.1.2), except at the first andlast points. Those give the qualifier “modified.”The modified midpoint method is a second-order method, like (16.1.2), but withthe advantage of requiring (asymptotically for large n) only one derivative evaluationper step h instead of the two required by second-order Runge-Kutta.
Perhaps thereare applications where the simplicity of (16.3.2), easily coded in-line in some otherprogram, recommends it. In general, however, use of the modified midpoint methodby itself will be dominated by the embedded Runge-Kutta method with adaptivestepsize control, as implemented in the preceding section.The usefulness of the modified midpoint method to the Bulirsch-Stoer technique(§16.4) derives from a “deep” result about equations (16.3.2), due to Gragg. It turnsout that the error of (16.3.2), expressed as a power series in h, the stepsize, containsonly even powers of h,724Chapter 16.Integration of Ordinary Differential Equations}CITED REFERENCES AND FURTHER READING:Gear, C.W. 1971, Numerical Initial Value Problems in Ordinary Differential Equations (EnglewoodCliffs, NJ: Prentice-Hall), §6.1.4.Stoer, J., and Bulirsch, R.
1980, Introduction to Numerical Analysis (New York: Springer-Verlag),§7.2.12.16.4 Richardson Extrapolation and theBulirsch-Stoer MethodThe techniques described in this section are not for differential equationscontaining nonsmooth functions. For example, you might have a differentialequation whose right-hand side involves a function that is evaluated by table look-upand interpolation. If so, go back to Runge-Kutta with adaptive stepsize choice:That method does an excellent job of feeling its way through rocky or discontinuousterrain.
It is also an excellent choice for quick-and-dirty, low-accuracy solutionof a set of equations. A second warning is that the techniques in this section arenot particularly good for differential equations that have singular points inside theinterval of integration. A regular solution must tiptoe very carefully across suchpoints. Runge-Kutta with adaptive stepsize can sometimes effect this; more generally,there are special techniques available for such problems, beyond our scope here.Apart from those two caveats, we believe that the Bulirsch-Stoer method,discussed in this section, is the best known way to obtain high-accuracy solutionsto ordinary differential equations with minimal computational effort.
(A possibleexception, infrequently encountered in practice, is discussed in §16.7.)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).ym=vector(1,nvar);yn=vector(1,nvar);h=htot/nstep;Stepsize this trip.for (i=1;i<=nvar;i++) {ym[i]=y[i];yn[i]=y[i]+h*dydx[i];First step.}x=xs+h;(*derivs)(x,yn,yout);Will use yout for temporary storage of derivah2=2.0*h;tives.for (n=2;n<=nstep;n++) {General step.for (i=1;i<=nvar;i++) {swap=ym[i]+h2*yout[i];ym[i]=yn[i];yn[i]=swap;}x += h;(*derivs)(x,yn,yout);}for (i=1;i<=nvar;i++)Last step.yout[i]=0.5*(ym[i]+yn[i]+h*yout[i]);free_vector(yn,1,nvar);free_vector(ym,1,nvar);.
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
