Press, Teukolsly, Vetterling, Flannery - Numerical Recipes in C

Numerical Recipes in CThe Art of Scientific ComputingSecond EditionWilliam H. PressHarvard-Smithsonian Center for AstrophysicsSaul A. TeukolskyDepartment of Physics, Cornell UniversityWilliam T. VetterlingPolaroid CorporationBrian P. FlanneryEXXON Research and Engineering CompanyCAMBRIDGE UNIVERSITY PRESSCambridge New York Port Chester Melbourne SydneyPublished by the Press Syndicate of the University of CambridgeThe Pitt Building, Trumpington Street, Cambridge CB2 1RP40 West 20th Street, New York, NY 10011-4211, USA10 Stamford Road, Oakleigh, Melbourne 3166, Australiac Cambridge University Press 1988, 1992Copyright except for §13.10 and Appendix B, which are placed into the public domain,and except for all other computer programs and procedures, which arec Numerical Recipes Software 1987, 1988, 1992, 1997Copyright All Rights Reserved.Some sections of this book were originally published, in different form, in Computersc American Institute of Physics, 1988–1992.in Physics magazine, Copyright First Edition originally published 1988; Second Edition originally published 1992.Reprinted with corrections, 1993, 1994, 1995, 1997.This reprinting is corrected to software version 2.08Printed in the United States of AmericaTypeset in TEXWithout an additional license to use the contained software, this book is intended asa text and reference book, for reading purposes only.

A free license for limited use of thesoftware by the individual owner of a copy of this book who personally types one or moreroutines into a single computer is granted under terms described on p. xvii. See the section“License Information” (pp. xvi–xviii) for information on obtaining more general licensesat low cost.Machine-readable media containing the software in this book, with included licensesfor use on a single screen, are available from Cambridge University Press. See theorder form at the back of the book, email to “orders@cup.org” (North America) or“trade@cup.cam.ac.uk” (rest of world), or write to Cambridge University Press, 110Midland Avenue, Port Chester, NY 10573 (USA), for further information.The software may also be downloaded, with immediate purchase of a licensealso possible, from the Numerical Recipes Software Web Site (http://www.nr.com).Unlicensed transfer of Numerical Recipes programs to any other format, or to any computerexcept one that is specifically licensed, is strictly prohibited.

Technical questions,corrections, and requests for information should be addressed to Numerical RecipesSoftware, P.O. Box 243, Cambridge, MA 02238 (USA), email “info@nr.com”, or fax781 863-1739.Library of Congress Cataloging in Publication DataNumerical recipes in C : the art of scientific computing / William H. Press. . . [et al.]. – 2nd ed.Includes bibliographical references (p.) and index.ISBN 0-521-43108-51.

Numerical analysis–Computer programs. 2. Science–Mathematics–Computer programs.3. C (Computer program language) I. Press, William H.QA297.N866 1992519.4028553–dc2092-8876A catalog record for this book is available from the British Library.ISBNISBNISBNISBNISBN00000521 43108521 43720521 43724521 57608521 5760752535BookExample book in CC diskette (IBM 3.5 , 1.44M)CDROM (IBM PC/Macintosh)CDROM (UNIX)ContentsPreface to the Second Edition1Preface to the First EditionxivLicense InformationxviComputer Programs by Chapter and SectionxixPreliminaries11.0 Introduction1.1 Program Organization and Control Structures1.2 Some C Conventions for Scientific Computing1.3 Error, Accuracy, and Stability23xi151528Solution of Linear Algebraic Equations322.0 Introduction2.1 Gauss-Jordan Elimination2.2 Gaussian Elimination with Backsubstitution2.3 LU Decomposition and Its Applications2.4 Tridiagonal and Band Diagonal Systems of Equations2.5 Iterative Improvement of a Solution to Linear Equations2.6 Singular Value Decomposition2.7 Sparse Linear Systems2.8 Vandermonde Matrices and Toeplitz Matrices2.9 Cholesky Decomposition2.10 QR Decomposition2.11 Is Matrix Inversion an N 3 Process?3236414350555971909698102Interpolation and Extrapolation3.0 Introduction3.1 Polynomial Interpolation and Extrapolation3.2 Rational Function Interpolation and Extrapolation3.3 Cubic Spline Interpolation3.4 How to Search an Ordered Table3.5 Coefficients of the Interpolating Polynomial3.6 Interpolation in Two or More Dimensionsv105105108111113117120123vi4ContentsIntegration of Functions4.0 Introduction4.1 Classical Formulas for Equally Spaced Abscissas4.2 Elementary Algorithms4.3 Romberg Integration4.4 Improper Integrals4.5 Gaussian Quadratures and Orthogonal Polynomials4.6 Multidimensional Integrals5Evaluation of Functions5.0 Introduction5.1 Series and Their Convergence5.2 Evaluation of Continued Fractions5.3 Polynomials and Rational Functions5.4 Complex Arithmetic5.5 Recurrence Relations and Clenshaw’s Recurrence Formula5.6 Quadratic and Cubic Equations5.7 Numerical Derivatives5.8 Chebyshev Approximation5.9 Derivatives or Integrals of a Chebyshev-approximated Function5.10 Polynomial Approximation from Chebyshev Coefficients5.11 Economization of Power Series5.12 Padé Approximants5.13 Rational Chebyshev Approximation5.14 Evaluation of Functions by Path Integration6Special Functions6.0 Introduction6.1 Gamma Function, Beta Function, Factorials, Binomial Coefficients6.2 Incomplete Gamma Function, Error Function, Chi-SquareProbability Function, Cumulative Poisson Function6.3 Exponential Integrals6.4 Incomplete Beta Function, Student’s Distribution, F-Distribution,Cumulative Binomial Distribution6.5 Bessel Functions of Integer Order6.6 Modified Bessel Functions of Integer Order6.7 Bessel Functions of Fractional Order, Airy Functions, SphericalBessel Functions6.8 Spherical Harmonics6.9 Fresnel Integrals, Cosine and Sine Integrals6.10 Dawson’s Integral6.11 Elliptic Integrals and Jacobian Elliptic Functions6.12 Hypergeometric Functions7Random Numbers7.0 Introduction7.1 Uniform Deviates129129130136140141147161165165165169173176178183186190195197198200204208212212213216222226230236240252255259261271274274275Contents7.2 Transformation Method: Exponential and Normal Deviates7.3 Rejection Method: Gamma, Poisson, Binomial Deviates7.4 Generation of Random Bits7.5 Random Sequences Based on Data Encryption7.6 Simple Monte Carlo Integration7.7 Quasi- (that is, Sub-) Random Sequences7.8 Adaptive and Recursive Monte Carlo Methods8Sorting8.0 Introduction8.1 Straight Insertion and Shell’s Method8.2 Quicksort8.3 Heapsort8.4 Indexing and Ranking8.5 Selecting the M th Largest8.6 Determination of Equivalence Classes9Root Finding and Nonlinear Sets of Equations9.0 Introduction9.1 Bracketing and Bisection9.2 Secant Method, False Position Method, and Ridders’ Method9.3 Van Wijngaarden–Dekker–Brent Method9.4 Newton-Raphson Method Using Derivative9.5 Roots of Polynomials9.6 Newton-Raphson Method for Nonlinear Systems of Equations9.7 Globally Convergent Methods for Nonlinear Systems of Equations10 Minimization or Maximization of Functions10.0 Introduction10.1 Golden Section Search in One Dimension10.2 Parabolic Interpolation and Brent’s Method in One Dimension10.3 One-Dimensional Search with First Derivatives10.4 Downhill Simplex Method in Multidimensions10.5 Direction Set (Powell’s) Methods in Multidimensions10.6 Conjugate Gradient Methods in Multidimensions10.7 Variable Metric Methods in Multidimensions10.8 Linear Programming and the Simplex Method10.9 Simulated Annealing Methods11 Eigensystems11.0 Introduction11.1 Jacobi Transformations of a Symmetric Matrix11.2 Reduction of a Symmetric Matrix to Tridiagonal Form:Givens and Householder Reductions11.3 Eigenvalues and Eigenvectors of a Tridiagonal Matrix11.4 Hermitian Matrices11.5 Reduction of a General Matrix to Hessenberg Formvii287290296300304309316329329330332336338341345347347350354359362369379383394394397402405408412420425430444456456463469475481482viiiContents11.6 The QR Algorithm for Real Hessenberg Matrices11.7 Improving Eigenvalues and/or Finding Eigenvectors byInverse Iteration12 Fast Fourier Transform12.0 Introduction12.1 Fourier Transform of Discretely Sampled Data12.2 Fast Fourier Transform (FFT)12.3 FFT of Real Functions, Sine and Cosine Transforms12.4 FFT in Two or More Dimensions12.5 Fourier Transforms of Real Data in Two and Three Dimensions12.6 External Storage or Memory-Local FFTs13 Fourier and Spectral Applications13.0 Introduction13.1 Convolution and Deconvolution Using the FFT13.2 Correlation and Autocorrelation Using the FFT13.3 Optimal (Wiener) Filtering with the FFT13.4 Power Spectrum Estimation Using the FFT13.5 Digital Filtering in the Time Domain13.6 Linear Prediction and Linear Predictive Coding13.7 Power Spectrum Estimation by the Maximum Entropy(All Poles) Method13.8 Spectral Analysis of Unevenly Sampled Data13.9 Computing Fourier Integrals Using the FFT13.10 Wavelet Transforms13.11 Numerical Use of the Sampling Theorem14 Statistical Description of Data14.0 Introduction14.1 Moments of a Distribution: Mean, Variance, Skewness,and So Forth14.2 Do Two Distributions Have the Same Means or Variances?14.3 Are Two Distributions Different?14.4 Contingency Table Analysis of Two Distributions14.5 Linear Correlation14.6 Nonparametric or Rank Correlation14.7 Do Two-Dimensional Distributions Differ?14.8 Savitzky-Golay Smoothing Filters15 Modeling of Data15.0 Introduction15.1 Least Squares as a Maximum Likelihood Estimator15.2 Fitting Data to a Straight Line15.3 Straight-Line Data with Errors in Both Coordinates15.4 General Linear Least Squares15.5 Nonlinear Models486493496496500504510521525532537537538545547549558564572575584591606609609610615620628636639645650656656657661666671681Contents15.6 Confidence Limits on Estimated Model Parameters15.7 Robust Estimation16 Integration of Ordinary Differential Equations16.0 Introduction16.1 Runge-Kutta Method16.2 Adaptive Stepsize Control for Runge-Kutta16.3 Modified Midpoint Method16.4 Richardson Extrapolation and the Bulirsch-Stoer Method16.5 Second-Order Conservative Equations16.6 Stiff Sets of Equations16.7 Multistep, Multivalue, and Predictor-Corrector Methods17 Two Point Boundary Value Problems17.0 Introduction17.1 The Shooting Method17.2 Shooting to a Fitting Point17.3 Relaxation Methods17.4 A Worked Example: Spheroidal Harmonics17.5 Automated Allocation of Mesh Points17.6 Handling Internal Boundary Conditions or Singular Points18 Integral Equations and Inverse Theory18.0 Introduction18.1 Fredholm Equations of the Second Kind18.2 Volterra Equations18.3 Integral Equations with Singular Kernels18.4 Inverse Problems and the Use of A Priori Information18.5 Linear Regularization Methods18.6 Backus-Gilbert Method18.7 Maximum Entropy Image Restoration19 Partial Differential Equations19.0 Introduction19.1 Flux-Conservative Initial Value Problems19.2 Diffusive Initial Value Problems19.3 Initial Value Problems in Multidimensions19.4 Fourier and Cyclic Reduction Methods for BoundaryValue Problems19.5 Relaxation Methods for Boundary Value Problems19.6 Multigrid Methods for Boundary Value Problems20 Less-Numerical Algorithms20.0 Introduction20.1 Diagnosing Machine Parameters20.2 Gray Codesix689699707707710714722724732734747753753757760762772783784788788791794797804808815818827827834847853857863871889889889894xContents20.3 Cyclic Redundancy and Other Checksums20.4 Huffman Coding and Compression of Data20.5 Arithmetic Coding20.6 Arithmetic at Arbitrary Precision896903910915References926Appendix A: Table of Prototype Declarations930Appendix B: Utility Routines940Appendix C: Complex Arithmetic948Index of Programs and Dependencies951General Index965Preface to the Second EditionOur aim in writing the original edition of Numerical Recipes was to provide abook that combined general discussion, analytical mathematics, algorithmics, andactual working programs.