c5-12 (779496)
Текст из файла
200Chapter 5.Evaluation of Functions#define NFEW ..#define NMANY ..float *c,*d,*e,a,b;Economize NMANY power series coefficients e[0..NMANY-1] in the range (a, b) into NFEWcoefficients d[0..NFEW-1].In our example, by the way, the 8th through 10th Chebyshev coefficients turn out tobe on the order of −7 × 10−6 , 3 × 10−7 , and −9 × 10−9 , so reasonable truncations (forsingle precision calculations) are somewhere in this range, yielding a polynomial with 8 –10 terms instead of the original 13.Replacing a 13-term polynomial with a (say) 10-term polynomial without any loss ofaccuracy — that does seem to be getting something for nothing.
Is there some magic inthis technique? Not really. The 13-term polynomial defined a function f (x). Equivalent toeconomizing the series, we could instead have evaluated f (x) at enough points to constructits Chebyshev approximation in the interval of interest, by the methods of §5.8. We wouldhave obtained just the same lower-order polynomial. The principal lesson is that the rateof convergence of Chebyshev coefficients has nothing to do with the rate of convergence ofpower series coefficients; and it is the former that dictates the number of terms needed in apolynomial approximation. A function might have a divergent power series in some regionof interest, but if the function itself is well-behaved, it will have perfectly good polynomialapproximations. These can be found by the methods of §5.8, but not by economization ofseries.
There is slightly less to economization of series than meets the eye.CITED REFERENCES AND FURTHER READING:Acton, F.S. 1970, Numerical Methods That Work; 1990, corrected edition (Washington: Mathematical Association of America), Chapter 12.Arfken, G. 1970, Mathematical Methods for Physicists, 2nd ed. (New York: Academic Press),p. 631. [1]5.12 Padé ApproximantsA Padé approximant, so called, is that rational function (of a specified order) whosepower series expansion agrees with a given power series to the highest possible order.
Ifthe rational function isMXR(x) ≡ak xkk=01+NXk=1(5.12.1)bk xkSample 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).c=vector(0,NMANY-1);d=vector(0,NFEW-1);e=vector(0,NMANY-1);pcshft((-2.0-b-a)/(b-a),(2.0-b-a)/(b-a),e,NMANY);pccheb(e,c,NMANY);...Here one would normally examine the Chebyshev coefficients c[0..NMANY-1] to decidehow small NFEW can be.chebpc(c,d,NFEW);pcshft(a,b,d,NFEW);5.12 Padé Approximants201then R(x) is said to be a Padé approximant to the seriesf (x) ≡∞Xc k xk(5.12.2)k=0ifR(0) = f (0)dkdkR(x)=f (x),kdxkdxx=0x=0(5.12.3)k = 1, 2, .
. . , M + N(5.12.4)Equations (5.12.3) and (5.12.4) furnish M + N + 1 equations for the unknowns a0 , . . . , aMand b1 , . . . , bN . The easiest way to see what these equations are is to equate (5.12.1) and(5.12.2), multiply both by the denominator of equation (5.12.1), and equate all powers ofx that have either a’s or b’s in their coefficients. If we consider only the special case ofa diagonal rational approximation, M = N (cf.
§3.2), then we have a0 = c0 , with theremaining a’s and b’s satisfyingNXbm cN −m+k = −cN +k ,k = 1, . . . , N(5.12.5)k = 1, . . . , N(5.12.6)m=1kXbm ck−m = ak ,m=0(note, in equation 5.12.1, that b0 = 1). To solve these, start with equations (5.12.5), whichare a set of linear equations for all the unknown b’s. Although the set is in the form of aToeplitz matrix (compare equation 2.8.8), experience shows that the equations are frequentlyclose to singular, so that one should not solve them by the methods of §2.8, but rather byfull LU decomposition. Additionally, it is a good idea to refine the solution by iterativeimprovement (routine mprove in §2.5) [1].Once the b’s are known, then equation (5.12.6) gives an explicit formula for the unknowna’s, completing the solution.Padé approximants are typically used when there is some unknown underlying functionf (x). We suppose that you are able somehow to compute, perhaps by laborious analyticexpansions, the values of f (x) and a few of its derivatives at x = 0: f (0), f 0 (0), f 00 (0),and so on.
These are of course the first few coefficients in the power series expansion off (x); but they are not necessarily getting small, and you have no idea where (or whether)the power series is convergent.By contrast with techniques like Chebyshev approximation (§5.8) or economizationof power series (§5.11) that only condense the information that you already know about afunction, Padé approximants can give you genuinely new information about your function’svalues. It is sometimes quite mysterious how well this can work. (Like other mysteries inmathematics, it relates to analyticity.) An example will illustrate.Imagine that, by extraordinary labors, you have ground out the first five terms in thepower series expansion of an unknown function f (x),f (x) ≈ 2 +11 249 3175 4x+x −x +x +···981874878732(5.12.7)(It is not really necessary that you know the coefficients in exact rational form — numericalvalues are just as good.
We here write them as rationals to give you the impression thatthey derive from some side analytic calculation.) Equation (5.12.7) is plotted as the curvelabeled “power series” in Figure 5.12.1. One sees that for x >∼ 4 it is dominated by itslargest, quartic, term.We now take the five coefficients in equation (5.12.7) and run them through the routinepade listed below. It returns five rational coefficients, three a’s and two b’s, for use in equation(5.12.1) with M = N = 2.
The curve in the figure labeled “Padé” plots the resulting rationalfunction. Note that both solid curves derive from the same five original coefficient values.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).and also202Chapter 5.Evaluation of Functions10f(x) = [7 + (1 + x)4/3]1/38f (x)power series (5 terms)4Padé (5 coefficients)exact200246810xFigure 5.12.1. The five-term power series expansion and the derived five-coefficient Padé approximantfor a sample function f (x).
The full power series converges only for x < 1. Note that the Padéapproximant maintains accuracy far outside the radius of convergence of the series.To evaluate the results, we need Deus ex machina (a useful fellow, when he is available)to tell us that equation (5.12.7) is in fact the power series expansion of the functionf (x) = [7 + (1 + x)4/3 ]1/3(5.12.8)which is plotted as the dotted curve in the figure. This function has a branch point at x = −1,so its power series is convergent only in the range −1 < x < 1.
In most of the rangeshown in the figure, the series is divergent, and the value of its truncation to five terms israther meaningless. Nevertheless, those five terms, converted to a Padé approximant, give aremarkably good representation of the function up to at least x ∼ 10.Why does this work? Are there not other functions with the same first five terms intheir power series, but completely different behavior in the range (say) 2 < x < 10? Indeedthere are.
Padé approximation has the uncanny knack of picking the function you had inmind from among all the possibilities. Except when it doesn’t! That is the downside ofPadé approximation: it is uncontrolled. There is, in general, no way to tell how accurateit is, or how far out in x it can usefully be extended. It is a powerful, but in the end stillmysterious, technique.Here is the routine that gets a’s and b’s from your c’s.
Note that the routine is specializedto the case M = N , and also that, on output, the rational coefficients are arranged in a formatfor use with the evaluation routine ratval (§5.3). (Also for consistency with that routine,the array of c’s is passed in double precision.)#include <math.h>#include "nrutil.h"#define BIG 1.0e30void pade(double cof[], int n, float *resid)Given cof[0..2*n], the leading terms in the power series expansion of a function, solve thelinear Padé equations to return the coefficients of a diagonal rational function approximation tothe same function, namely (cof[0] + cof[1]x + · · · + cof[n]xN )/(1 + cof[n+1]x + · · · +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).65.12 Padé Approximants203cof[2*n]xN ). The value resid is the norm of the residual vector; a small value indicates awell-converged solution. Note that cof is double precision for consistency with ratval.{indx=ivector(1,n);q=matrix(1,n,1,n);qlu=matrix(1,n,1,n);x=vector(1,n);y=vector(1,n);z=vector(1,n);for (j=1;j<=n;j++) {y[j]=x[j]=cof[n+j];for (k=1;k<=n;k++) {q[j][k]=cof[j-k+n];qlu[j][k]=q[j][k];}}ludcmp(qlu,n,indx,&d);lubksb(qlu,n,indx,x);rr=BIG;do {rrold=rr;for (j=1;j<=n;j++) z[j]=x[j];mprove(q,qlu,n,indx,y,x);for (rr=0.0,j=1;j<=n;j++)rr += SQR(z[j]-x[j]);} while (rr < rrold);*resid=sqrt(rr);for (k=1;k<=n;k++) {for (sum=cof[k],j=1;j<=k;j++)y[k]=sum;}for (j=1;j<=n;j++) {cof[j]=y[j];cof[j+n] = -x[j];}free_vector(z,1,n);free_vector(y,1,n);free_vector(x,1,n);free_matrix(qlu,1,n,1,n);free_matrix(q,1,n,1,n);free_ivector(indx,1,n);Set up matrix for solving.Solve by LU decomposition and backsubstitution.Important to use iterative improvement, sincethe Padé equations tend to be ill-conditioned.Calculate residual.If it is no longer improving, call it quits.Calculate the remaining coefficients.sum -= x[j]*cof[k-j];Copy answers to output.}CITED REFERENCES AND FURTHER READING:Ralston, A.
Характеристики
Тип файла PDF
PDF-формат наиболее широко используется для просмотра любого типа файлов на любом устройстве. В него можно сохранить документ, таблицы, презентацию, текст, чертежи, вычисления, графики и всё остальное, что можно показать на экране любого устройства. Именно его лучше всего использовать для печати.
Например, если Вам нужно распечатать чертёж из автокада, Вы сохраните чертёж на флешку, но будет ли автокад в пункте печати? А если будет, то нужная версия с нужными библиотеками? Именно для этого и нужен формат PDF - в нём точно будет показано верно вне зависимости от того, в какой программе создали PDF-файл и есть ли нужная программа для его просмотра.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
