c3-3 (779473), страница 2

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 2)

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).for (i=2;i<=n-1;i++) {This is the decomposition loop of the tridiagonal alsig=(x[i]-x[i-1])/(x[i+1]-x[i-1]);gorithm. y2 and u are used for temp=sig*y2[i-1]+2.0;porary storage of the decomposedy2[i]=(sig-1.0)/p;factors.u[i]=(y[i+1]-y[i])/(x[i+1]-x[i]) - (y[i]-y[i-1])/(x[i]-x[i-1]);u[i]=(6.0*u[i]/(x[i+1]-x[i-1])-sig*u[i-1])/p;}if (ypn > 0.99e30)The upper boundary condition is set either to beqn=un=0.0;“natural”else {or else to have a specified first derivative.qn=0.5;un=(3.0/(x[n]-x[n-1]))*(ypn-(y[n]-y[n-1])/(x[n]-x[n-1]));}y2[n]=(un-qn*u[n-1])/(qn*y2[n-1]+1.0);for (k=n-1;k>=1;k--)This is the backsubstitution loop of the tridiagonaly2[k]=y2[k]*y2[k+1]+u[k];algorithm.free_vector(u,1,n-1);3.4 How to Search an Ordered Table1173.4 How to Search an Ordered Tablevoid locate(float xx[], unsigned long n, float x, unsigned long *j)Given an array xx[1..n], and given a value x, returns a value j such that x is between xx[j]and xx[j+1].

xx must be monotonic, either increasing or decreasing. j=0 or j=n is returnedto indicate that x is out of range.{unsigned long ju,jm,jl;int ascnd;}jl=0;Initialize lowerju=n+1;and upper limits.ascnd=(xx[n] >= xx[1]);while (ju-jl > 1) {If we are not yet done,jm=(ju+jl) >> 1;compute a midpoint,if (x >= xx[jm] == ascnd)jl=jm;and replace either the lower limitelseju=jm;or the upper limit, as appropriate.}Repeat until the test condition is satisfied.if (x == xx[1]) *j=1;Then set the outputelse if(x == xx[n]) *j=n-1;else *j=jl;and return.A unit-offset array xx is assumed.

To use locate with a zero-offset array,remember to subtract 1 from the address of xx, and also from the returned value j.Search with Correlated ValuesSometimes you will be in the situation of searching a large table many times,and with nearly identical abscissas on consecutive searches. For example, youmay be generating a function that is used on the right-hand side of a differentialequation: Most differential-equation integrators, as we shall see in Chapter 16, callSample 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).Suppose that you have decided to use some particular interpolation scheme,such as fourth-order polynomial interpolation, to compute a function f(x) from aset of tabulated xi ’s and fi ’s. Then you will need a fast way of finding your placein the table of xi ’s, given some particular value x at which the function evaluationis desired.

This problem is not properly one of numerical analysis, but it occurs sooften in practice that it would be negligent of us to ignore it.Formally, the problem is this: Given an array of abscissas xx[j], j=1, 2, . . . ,n,with the elements either monotonically increasing or monotonically decreasing, andgiven a number x, find an integer j such that x lies between xx[j] and xx[j+1].For this task, let us define fictitious array elements xx[0] and xx[n+1] equal toplus or minus infinity (in whichever order is consistent with the monotonicity of thetable).

Then j will always be between 0 and n, inclusive; a value of 0 indicates“off-scale” at one end of the table, n indicates off-scale at the other end.In most cases, when all is said and done, it is hard to do better than bisection,which will find the right place in the table in about log2 n tries. We already did usebisection in the spline evaluation routine splint of the preceding section, so youmight glance back at that. Standing by itself, a bisection routine looks like this:.

Характеристики

Список файлов книги