Nash - Compact Numerical Methods for Computers (523163), страница 32
Текст из файла (страница 32)
Theonly genuine polynomial root-finding problem I have encountered in practice isthe internal rate of return (example 12.4). However, accountants and economistshave very good ideas about where they would like the root to be found, so I havenot tried to develop general methods for finding all the roots of a polynomial, forinstance by methods such as those discussed by Jenkins and Traub (1975). Someexperiments I have performed with S G Nash (unpublished) on the use of matrixeigenvalue algorithms applied to the companion matrices of polynomials werenot very encouraging as to accuracy or speed, even though we had expected suchmethods to be slow.13.2.
THE LINEAR SEARCH PROBLEMThe linear search problem can be stated:minimise S(b) with respect to b.(13.1)However, it is not usual to leave the domain of b unrestricted. In all casesconsidered here, b is real and will often be confined to some interval [u, v ].If S(b) is differentiable, this problem can be approached by applying a rootfinding algorithm toS'(b) = dS(b)/ db.148(13.2)One-dimensional problems149Since local maxima also zero the derivative of the function S(b), such solutionswill have to be checked either by ensuring the second derivative S”(b) is positiveor equivalently by examining the values of the function at points near thesupposed minimum.When the derivative of S is not available or is expensive to compute, a methodfor minimisation along a line is needed which depends only on function values.Obviously, any method which evaluates the function at a finite number of points isnot guaranteed to detect deep, narrow wells in the function.
Therefore, someassumption must be made about the function on the interval [u,v ]. Here it will beassumed that the function S(b) is unimodal in [u, v], that is, that there is only onestationary value (either a maximum or a minimum) in the interval. In the casewhere the stationary value is a maximum, the minimum of the function will beeither at u or at v.Given that S(b) is unimodal in [u, v], a grid or equal-interval search could beused to decrease the size of the interval. For instance, the function values could becomputed for each of the pointsb j = u+ jhj = 0, 1, 2, .
. . , n(13.3)whereh = (v - u)/n.(13.4)If the smallest of these values is S(bk), then the minimum lies in the interval[bk-l ,bk+1]. Equation (13.3) can be used to compute the endpoints of thisinterval. If S(b0) or S(bn) is the smallest value, then the interval is [u , b 1] or[bn- 1, u], though for simplicity this could be left out of a program. The search cannow be repeated, noting, of course, that the function values at the endpoints ofthe interval have already been computed.Algorithm 16. Grid search along a lineprocedure gridsrch( var 1bound, ubound : real; {the lower andupper bounds to the interval to be tested}nint : integer; {the number of grid intervals}var finin: real; {the lowest function value found}var minarg: integer; {the grid point in the set0,1,...,nint at which the minimum function valuewas found}var changarg: integer {the grid point in the set1,2,...,nint which is nearest the upper boundubound such that there has been a sign changein the function values f(lbound+(changarg-1)*h)and f(lbound+changarg*h) where h is the step== (ubound - lbound)/nint} );{alg16.pas == one-dimensional grid search over function valuesThis version halts execution if the function is not computable at allthe grid points.
Note that it is not equivalent to the version in thefirst edition of Compact Numerical Methods.Copyright 1988 J.C.Nash}150Compact numerical methods for computersAlgorithm 16. Grid search along a line (cont.)varj : integer;h, p, t : real;notcomp : boolean;beginwriteln(‘alg16.pas-- one-dimensional grid search’);{STEP 0 via the procedure call}writeln(‘In gridsrch 1bound=‘,lbound,’ubound=‘,ubound);notcomp:=false; {function must be called with notcomp false orroot will be displayed -- STEP 1}t:=fnld(lbound, notcomp); {compute function at lower bound and setpointer k to lowest function value found so far i.e.
the 0’th}writeln(’lbf(‘,lbound,‘)=‘,t);if notcomp then halt;fmin:=t; {to save the function value}minarg:=0; {so far this is the lowest value found}changarg:=0; {as a safety setting of this value}h:=(ubound-lbound)/nint;for j:=l to nint do {STEP 2}{Note: we increase the number of steps so that the upper bound ubound isnow a function argument.
Warning: because the argument is nowcalculated, we may get a value slightly different from ubound informing (lbound+nint*h). }beginp:=fn1d(lbound+j*h, notcomp); {STEP 3}write(’f(‘,lbound+j*h,‘)=‘,p);if notcomp then halt;if p<fmin then {STEP 4}begin {STEP 5}fmin:=p; minarg:=j;end; {if p<fmin}if p*t<=0 thenbeginwriteln(‘ *** sign change ***‘);changarg:=j; {to save change point argument}endelsebeginwriteln; {no action since sign change}end;t:=p; {to save latest function value} {STEP 6}end; {loop on j}writeln( ‘Minimum so far is f( ‘,lbound+minarg*h,‘)=‘,fmin);if changarg>0 thenbeginwritehr(‘Sign change observed last in interval ‘);writeln(‘[‘,lbound+(changarg-l)*h,‘,‘,lbound+changarg*h,‘]‘);endelsewriteln(‘Apparently no sign change in [ ‘,lbound,‘,‘,ubound,‘]‘);end; {alg16.pas == gridsrch}One-dimensional problems151A grid or equal-interval search is not a particularly efficient method for findinga minimum.
However, it is very simple to program and provides a set of points forplotting the function along the line. Even a crude plot, such as that producedusing a line printer or teletypewriter, can be used to gain a visual appreciation ofthe function studied. It is in this fashion that a grid search can be used mosteffectively, not to find the minimum of a function but to carry out a preliminaryinvestigation of it. Furthermore, if the grid search is used in this way, there is nonecessity that the function be unimodal.Consider now the case where instead of the interval [u, v], only a startingposition u and a step h are given. If S(u + h) < S(u), the step will be termed asuccess, otherwise it will be called a failure.
This is the basis of Rosenbrock’s(1960) algorithm for minimising a function of several variables. Here, however,the goal will be to find points at which the function value is reduced. A procedurewhich does this can be repeated until a minimum is found.Suppose S(u + h) < S(u), then an obvious tactic is to replace u by (u + h) andtry again, perhaps increasing the step somewhat.
If, however, S(u + h) > S(u),then either the step is too large or it is in the wrong direction. In this case apossibility is to reduce the size of h and change its sign. (If only the size isincorrect, two operations of this sort will find a lower function value.)The above paragraph defines the success-failure algorithm except for thefactors to be used to increase or decrease the size of h . Notice that so far neitherthe grid search nor the success-failure algorithm take into account the functionvalues apart from the fact that one function value is larger than another.
Howmuch larger it is, which indicates how rapidly the function is changing with respectto h, is not considered. Using the function values can, however, lead to moreefficient minimisation procedures. The idea is to use the function values and theirassociated points (as well as the derivatives of the function at these points ifavailable) to generate some simple function, call it I(b), which interpolates thefunction S(b) between the points at which S has been computed. The mostpopular choices for I(b) are polynomials, in particular of degree 2 (parabolicapproximation) or degree 3 (cubic approximation).
Consider the case where S( b)is known at three points b 0, b1 and b 2. Then using a parabola to interpolate thefunction requiresfor j=0,1,2(13.5)which provides three linear equations for the three unknowns A, B and C. Oncethese are found, it is simple to set the derivative of the interpolant to zerodI( b) / db = 2Cb + B = 0(13.6)to find the value of b which minimises the interpolating polynomial I(b). Thispresumes that the parabola has a minimum at the stationary point, and upon thispresumption hang many ifs, ands and buts, and no doubt miles of computer listingof useless results.