c6-11 (779512), страница 3
Текст из файла (страница 3)
Theargument ranges are 0 ≤ φ ≤ π/2, 0 ≤ k sin φ ≤ 1.{float rf(float x, float y, float z);float s;s=sin(phi);return s*rf(SQR(cos(phi)),(1.0-s*ak)*(1.0+s*ak),1.0);}#include <math.h>#include "nrutil.h"float elle(float phi, float ak)Legendre elliptic integral of the 2nd kind E(φ, k), evaluated using Carlson’s functions RD andRF . The argument ranges are 0 ≤ φ ≤ π/2, 0 ≤ k sin φ ≤ 1.{float rd(float x, float y, float z);float rf(float x, float y, float z);float cc,q,s;s=sin(phi);cc=SQR(cos(phi));q=(1.0-s*ak)*(1.0+s*ak);return s*(rf(cc,q,1.0)-(SQR(s*ak))*rd(cc,q,1.0)/3.0);}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).= sin φRF (cos2 φ, 1 − k 2 sin2 φ, 1)6.11 Elliptic Integrals and Jacobian Elliptic Functions269#include <math.h>#include "nrutil.h"s=sin(phi);enss=en*s*s;cc=SQR(cos(phi));q=(1.0-s*ak)*(1.0+s*ak);return s*(rf(cc,q,1.0)-enss*rj(cc,q,1.0,1.0+enss)/3.0);}Carlson’s functions are homogeneous of degree − 12 and − 32 , soRF (λx, λy, λz) = λ−1/2 RF (x, y, z)RJ (λx, λy, λz, λp) = λ−3/2 RJ (x, y, z, p)(6.11.22)Thus to express a Carlson function in Legendre’s notation, permute the first threearguments into ascending order, use homogeneity to scale the third argument to be1, and then use equations (6.11.19)–(6.11.21).Jacobian Elliptic FunctionsThe Jacobian elliptic function sn is defined as follows: instead of consideringthe elliptic integralu(y, k) ≡ u = F (φ, k)(6.11.23)consider the inverse functiony = sin φ = sn(u, k)(6.11.24)Equivalently,Z snu=0dyp2(1 − y )(1 − k 2 y2 )(6.11.25)When k = 0, sn is just sin.
The functions cn and dn are defined by the relationssn2 + cn2 = 1,k 2 sn2 + dn2 = 1(6.11.26)The routine given below actually takes mc ≡ kc2 = 1 − k 2 as an input parameter.It also computes all three functions sn, cn, and dn since computing all three is noharder than computing any one of them. For a description of the method, see [8].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).float ellpi(float phi, float en, float ak)Legendre elliptic integral of the 3rd kind Π(φ, n, k), evaluated using Carlson’s functions RJ andRF .
(Note that the sign convention on n is opposite that of Abramowitz and Stegun.) Theranges of φ and k are 0 ≤ φ ≤ π/2, 0 ≤ k sin φ ≤ 1.{float rf(float x, float y, float z);float rj(float x, float y, float z, float p);float cc,enss,q,s;270Chapter 6.#include <math.h>#define CA 0.0003Special FunctionsThe accuracy is the square of CA.emc=emmc;u=uu;if (emc) {bo=(emc < 0.0);if (bo) {d=1.0-emc;emc /= -1.0/d;u *= (d=sqrt(d));}a=1.0;*dn=1.0;for (i=1;i<=13;i++) {l=i;em[i]=a;en[i]=(emc=sqrt(emc));c=0.5*(a+emc);if (fabs(a-emc) <= CA*a) break;emc *= a;a=c;}u *= c;*sn=sin(u);*cn=cos(u);if (*sn) {a=(*cn)/(*sn);c *= a;for (ii=l;ii>=1;ii--) {b=em[ii];a *= c;c *= (*dn);*dn=(en[ii]+a)/(b+a);a=c/b;}a=1.0/sqrt(c*c+1.0);*sn=(*sn >= 0.0 ? a : -a);*cn=c*(*sn);}if (bo) {a=(*dn);*dn=(*cn);*cn=a;*sn /= d;}} else {*cn=1.0/cosh(u);*dn=(*cn);*sn=tanh(u);}}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).void sncndn(float uu, float emmc, float *sn, float *cn, float *dn)Returns the Jacobian elliptic functions sn(u, kc ), cn(u, kc ), and dn(u, kc ).
Here uu = u, whileemmc = kc2 .{float a,b,c,d,emc,u;float em[14],en[14];int i,ii,l,bo;6.12 Hypergeometric Functions271CITED REFERENCES AND FURTHER READING:Erdélyi, A., Magnus, W., Oberhettinger, F., and Tricomi, F.G. 1953, Higher TranscendentalFunctions, Vol. II, (New York: McGraw-Hill). [1]Gradshteyn, I.S., and Ryzhik, I.W. 1980, Table of Integrals, Series, and Products (New York:Academic Press).
[2]Carlson, B.C. 1977, SIAM Journal on Mathematical Analysis, vol. 8, pp. 231–242. [3]Bulirsch, R. 1965, Numerische Mathematik, vol. 7, pp. 78–90; 1965, op. cit., vol. 7, pp. 353–354;1969, op. cit., vol. 13, pp. 305–315. [8]Carlson, B.C. 1979, Numerische Mathematik, vol. 33, pp.
1–16. [9]Carlson, B.C., and Notis, E.M. 1981, ACM Transactions on Mathematical Software, vol. 7,pp. 398–403. [10]Carlson, B.C. 1978, SIAM Journal on Mathematical Analysis, vol. 9, p. 524–528. [11]Abramowitz, M., and Stegun, I.A. 1964, Handbook of Mathematical Functions, Applied Mathematics Series, Volume 55 (Washington: National Bureau of Standards; reprinted 1968 byDover Publications, New York), Chapter 17. [12]Mathews, J., and Walker, R.L. 1970, Mathematical Methods of Physics, 2nd ed. (Reading, MA:W.A.
Benjamin/Addison-Wesley), pp. 78–79.6.12 Hypergeometric FunctionsAs was discussed in §5.14, a fast, general routine for the the complex hypergeometric function 2 F1 (a, b, c; z), is difficult or impossible. The function is defined asthe analytic continuation of the hypergeometric series,a(a + 1)b(b + 1) z 2ab z++···c 1!c(c + 1)2!a(a + 1) . . . (a + j − 1)b(b + 1) .
. . (b + j − 1) z j+···+c(c + 1) . . . (c + j − 1)j!(6.12.1)This series converges only within the unit circle |z| < 1 (see [1]), but one’s interestin the function is not confined to this region.Section 5.14 discussed the method of evaluating this function by direct pathintegration in the complex plane. We here merely list the routines that result.Implementation of the function hypgeo is straightforward, and is described bycomments in the program. The machinery associated with Chapter 16’s routine forintegrating differential equations, odeint, is only minimally intrusive, and need noteven be completely understood: use of odeint requires one zeroed global variable,one function call, and a prescribed format for the derivative routine hypdrv.The function hypgeo will fail, of course, for values of z too close to thesingularity at 1.
(If you need to approach this singularity, or the one at ∞, usethe “linear transformation formulas” in §15.3 of [1].) Away from z = 1, and formoderate values of a, b, c, it is often remarkable how few steps are required tointegrate the equations. A half-dozen is typical.2 F1 (a, b, c; z)=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).Carlson, B.C. 1987, Mathematics of Computation, vol. 49, pp. 595–606 [4]; 1988, op.
cit., vol. 51,pp. 267–280 [5]; 1989, op. cit., vol. 53, pp. 327–333 [6]; 1991, op. cit., vol. 56, pp. 267–280.[7].















