c5-7 (779491)
Текст из файла
186Chapter 5.Evaluation of Functions5.7 Numerical Derivativesf 0 (x) ≈f(x + h) − f(x)h(5.7.1)practically suggests the program: Pick a small value h; evaluate f(x + h); youprobably have f(x) already evaluated, but if not, do it too; finally apply equation(5.7.1). What more needs to be said?Quite a lot, actually.
Applied uncritically, the above procedure is almostguaranteed to produce inaccurate results. Applied properly, it can be the right wayto compute a derivative only when the function f is fiercely expensive to compute,when you already have invested in computing f(x), and when, therefore, you wantto get the derivative in no more than a single additional function evaluation.
In sucha situation, the remaining issue is to choose h properly, an issue we now discuss:There are two sources of error in equation (5.7.1), truncation error and roundofferror. The truncation error comes from higher terms in the Taylor series expansion,11f(x + h) = f(x) + hf 0 (x) + h2 f 00 (x) + h3 f 000 (x) + · · ·26(5.7.2)whence1f(x + h) − f(x)= f 0 + hf 00 + · · ·h2(5.7.3)The roundoff error has various contributions. First there is roundoff error in h:Suppose, by way of an example, that you are at a point x = 10.3 and you blindlychoose h = 0.0001.
Neither x = 10.3 nor x + h = 10.30001 is a number withan exact representation in binary; each is therefore represented with some fractionalerror characteristic of the machine’s floating-point format, m , whose value in singleprecision may be ∼ 10−7 . The error in the effective value of h, namely the differencebetween x + h and x as represented in the machine, is therefore on the order of m x,which implies a fractional error in h of order ∼ m x/h ∼ 10−2 ! By equation (5.7.1)this immediately implies at least the same large fractional error in the derivative.We arrive at Lesson 1: Always choose h so that x + h and x differ by an exactlyrepresentable number. This can usually be accomplished by the program stepstemp = x + hh = temp − x(5.7.4)Some optimizing compilers, and some computers whose floating-point chips havehigher internal accuracy than is stored externally, can foil this trick; if so, it isusually enough to declare temp as volatile, or else to call a dummy functiondonothing(temp) between the two equations (5.7.4).
This forces temp into andout of addressable memory.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).Imagine that you have a procedure which computes a function f(x), and nowyou want to compute its derivative f 0 (x).
Easy, right? The definition of thederivative, the limit as h → 0 of1875.7 Numerical Derivativeswhere xc ≡ (f/f 00 )1/2 is the “curvature scale” of the function f, or “characteristicscale” over which it changes. In the absence of any other information, one oftenassumes xc = x (except near x = 0 where some other estimate of the typical xscale should be used).With the choice of equation (5.7.5), the fractional accuracy of the computedderivative is(er + et )/|f 0 | ∼√f (ff 00 /f 0 )1/2 ∼2√f(5.7.6)Here the last order-of-magnitude equality assumes that f, f 0 , and f 00 all sharethe same characteristic length scale, usually the case. One sees that the simplefinite-difference equation (5.7.1) gives at best only the square root of the machineaccuracy m .If you can afford two function evaluations for each derivative calculation, thenit is significantly better to use the symmetrized formf 0 (x) ≈f(x + h) − f(x − h)2h(5.7.7)In this case, by equation (5.7.2), the truncation error is et ∼ h2 f 000 .
The roundofferror er is about the same as before. The optimal choice of h, by a short calculationanalogous to the one above, is nowh∼f ff 0001/3∼ (f )1/3 xc(5.7.8)and the fractional error is(er + et )/|f 0 | ∼ (f )2/3 f 2/3 (f 000 )1/3 /f 0 ∼ (f )2/3(5.7.9)which will typically be an order of magnitude (single precision) or two orders ofmagnitude (double precision) better than equation (5.7.6).
We have arrived at Lesson2: Choose h to be the correct power of f or m times a characteristic scale xc.You can easily derive the correct powers for other cases [1]. For a function oftwo dimensions, for example, and the mixed derivative formula[f(x + h, y + h) − f(x + h, y − h)] − [f(x − h, y + h) − f(x − h, y − h)]∂2f=∂x∂y4h2(5.7.10)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).With h an “exact” number, the roundoff error in equation (5.7.1) is er ∼f |f(x)/h|. Here f is the fractional accuracy with which f is computed; for asimple function this may be comparable to the machine accuracy, f ≈ m , but for acomplicated calculation with additional sources of inaccuracy it may be larger.
Thetruncation error in equation (5.7.3) is on the order of et ∼ |hf 00 (x)|. Varying h tominimize the sum er + et gives the optimal choice of h,sf f√h∼≈ f x c(5.7.5)00f188Chapter 5.Evaluation of Functions1/3#include <math.h>#include "nrutil.h"#define CON 1.4#define CON2 (CON*CON)#define BIG 1.0e30#define NTAB 10#define SAFE 2.0Stepsize is decreased by CON at each iteration.Sets maximum size of tableau.Return when error is SAFE worse than the best sofar.float dfridr(float (*func)(float), float x, float h, float *err)Returns the derivative of a function func at a point x by Ridders’ method of polynomialextrapolation. The value h is input as an estimated initial stepsize; it need not be small, butrather should be an increment in x over which func changes substantially.
An estimate of theerror in the derivative is returned as err.{int i,j;float errt,fac,hh,**a,ans;if (h == 0.0) nrerror("h must be nonzero in dfridr.");a=matrix(1,NTAB,1,NTAB);hh=h;a[1][1]=((*func)(x+hh)-(*func)(x-hh))/(2.0*hh);*err=BIG;for (i=2;i<=NTAB;i++) {Successive columns in the Neville tableau will go to smaller stepsizes and higher orders ofextrapolation.hh /= CON;a[1][i]=((*func)(x+hh)-(*func)(x-hh))/(2.0*hh);Try new, smaller stepfac=CON2;size.for (j=2;j<=i;j++) {Compute extrapolations of various orders, requiringa[j][i]=(a[j-1][i]*fac-a[j-1][i-1])/(fac-1.0);no new function evalfac=CON2*fac;uations.errt=FMAX(fabs(a[j][i]-a[j-1][i]),fabs(a[j][i]-a[j-1][i-1]));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).the correct scaling is h ∼ f xc .It is disappointing, certainly, that no simple finite-difference formula likeequation (5.7.1) or (5.7.7) gives an accuracy comparable to the machine accuracym , or even the lower accuracy to which f is evaluated, f .
Are there no bettermethods?Yes, there are. All, however, involve exploration of the function’s behaviorover scales comparable to xc , plus some assumption of smoothness, or analyticity,so that the high-order terms in a Taylor expansion like equation (5.7.2) have somemeaning. Such methods also involve multiple evaluations of the function f, so theirincreased accuracy must be weighed against increased cost.The general idea of “Richardson’s deferred approach to the limit” is particularlyattractive. For numerical integrals, that idea leads to so-called Romberg integration(for review, see §4.3).
For derivatives, one seeks to extrapolate, to h → 0, the resultof finite-difference calculations with smaller and smaller finite values of h. By theuse of Neville’s algorithm (§3.1), one uses each new finite-difference calculation toproduce both an extrapolation of higher order, and also extrapolations of previous,lower, orders but with smaller scales h. Ridders [2] has given a nice implementationof this idea; the following program, dfridr, is based on his algorithm, modified byan improved termination criterion. Input to the routine is a function f (called func),a position x, and a largest stepsize h (more analogous to what we have called xcabove than to what we have called h). Output is the returned value of the derivative,and an estimate of its error, err.5.7 Numerical Derivatives189The error strategy is to compare each new extrapolation to one order lower, bothat the present stepsize and the previous one.if (errt <= *err) {If error is decreased, save the improved answer.*err=errt;ans=a[j][i];}}free_matrix(a,1,NTAB,1,NTAB);return ans;}In dfridr, the number of evaluations of func is typically 6 to 12, but is allowedto be as great as 2×NTAB.
As a function of input h, it is typical for the accuracyto get better as h is made larger, until a sudden point is reached where nonsensicalextrapolation produces early return with a large error. You should therefore choosea fairly large value for h, but monitor the returned value err, decreasing h if it isnot small. For functions whose characteristic x scale is of order unity, we typicallytake h to be a few tenths.Besides Ridders’ method, there are other possible techniques.
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.















