10.3 One-Dimensional Search with First Derivatives 407 the minimum function value is returned as dbrent,the returned function value int iter,ok1,ok2; Will be used as flags for whether pro- float a,b,d,d1,d2,du,dv,dw,dx,e=0.0; posed steps are acceptable or not. float fu,fv,fw,fx,olde,tol1,tol2,u,ul,u2,v,w,x,xm; Comments following will point out only differences from the routine brent.Read that routine first. a=(ax<cx？ax:cx); b=(ax>cx？ax:cx); x=w=v=bx; fw-=fv=fx=(*f)(x): dw=dv=dx=(*df)(x); All our housekeeping chores are dou- for (iter=1;iter<=ITMAX;iter++){ bled by the necessity of moving xm=0.5*(a+b); derivative values around as well 83 tol1-tol*fabs(x)+ZEPS: as function values. 鱼 granted for 18881992 to12=2.0*to11; if (fabs(x-xm) <(to12-0.5*(b-a)) 1.800 *xmin=x; return fx; if (fabs(e)>tol1){ from NUMERICAL RECIPES IN d1=2.0*(b-a): Initialize these d's to an out-of-bracket d2=d1; value. (North to make if (dw !dx)d1=(w-x)*dx/(dx-dw); Secant method with one point. if (dv !dx)d2=(v-x)*dx/(dx-dv); And the other. Which of these two estimates of d shall we take?We will insist that they be within Ameri computer, one paper /Cambridge University Press. THE the bracket,and on the side pointed to by the derivative at x: ART u1=x+d1; u2=x+d2; 是 ok1=(a-u1)*(u1-b)>0.0&dx*d1<=0.0; ok2=(a-u2)*(u2-b)>0.0dx*d2<=0.0; copy for their Programs olde=e; Movement on the step before last. e=d; strictly prohibited if (ok1 Il ok2){ Take only an acceptable d,and if if (ok1&&ok2) both are acceptable,then take d=(fabs(d1)<fabs(d2)?d1 d2); the smallest one. else if (ok1) d=d1; else ded2 if (fabs(d)<=fabs(0.5*olde)){ u=x+d; if (u-a tol2 II b-u tol2) d=SIGN(tol1,xm-x); else Bisect,not golden section. Numerical Recipes books or 1988-1992 by Numerical Recipes OF SCIENTIFIC COMPUTING(ISBN 0-521 d=0.5*(e=(dx>=0.0？a-x:b-x)); -431085 Decide which segment by the sign of the derivative. else d=0.5*(e=(dx>=0.0？a-x:b-x)); (outside North America) Software. else d=0.5*(e=(dx>=0.0？a-x:b-x)); if (fabs(d)>=toll){ machine u=x+d; fu=(*f)(u) else u=x+SIGN(tol1,d); fu=(*f)(u); if (fu>fx){ If the minimum step in the downhill 米Xxm1n=X; direction takes us uphill,then return fx; we are done10.3 One-Dimensional Search with First Derivatives 407 Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). the minimum function value is returned as dbrent, the returned function value. { int iter,ok1,ok2; Will be used as flags for whether pro￾float a,b,d,d1,d2,du,dv,dw,dx,e=0.0; posed steps are acceptable or not. float fu,fv,fw,fx,olde,tol1,tol2,u,u1,u2,v,w,x,xm; Comments following will point out only differences from the routine brent. Read that routine first. a=(ax < cx ? ax : cx); b=(ax > cx ? ax : cx); x=w=v=bx; fw=fv=fx=(*f)(x); dw=dv=dx=(*df)(x); All our housekeeping chores are dou￾bled by the necessity of moving derivative values around as well as function values. for (iter=1;iter<=ITMAX;iter++) { xm=0.5*(a+b); tol1=tol*fabs(x)+ZEPS; tol2=2.0*tol1; if (fabs(x-xm) <= (tol2-0.5*(b-a))) { *xmin=x; return fx; } if (fabs(e) > tol1) { d1=2.0*(b-a); Initialize these d’s to an out-of-bracket d2=d1; value. if (dw != dx) d1=(w-x)*dx/(dx-dw); Secant method with one point. if (dv != dx) d2=(v-x)*dx/(dx-dv); And the other. Which of these two estimates of d shall we take? We will insist that they be within the bracket, and on the side pointed to by the derivative at x: u1=x+d1; u2=x+d2; ok1 = (a-u1)*(u1-b) > 0.0 && dx*d1 <= 0.0; ok2 = (a-u2)*(u2-b) > 0.0 && dx*d2 <= 0.0; olde=e; Movement on the step before last. e=d; if (ok1 || ok2) { Take only an acceptable d, and if both are acceptable, then take the smallest one. if (ok1 && ok2) d=(fabs(d1) < fabs(d2) ? d1 : d2); else if (ok1) d=d1; else d=d2; if (fabs(d) <= fabs(0.5*olde)) { u=x+d; if (u-a < tol2 || b-u < tol2) d=SIGN(tol1,xm-x); } else { Bisect, not golden section. d=0.5*(e=(dx >= 0.0 ? a-x : b-x)); Decide which segment by the sign of the derivative. } } else { d=0.5*(e=(dx >= 0.0 ? a-x : b-x)); } } else { d=0.5*(e=(dx >= 0.0 ? a-x : b-x)); } if (fabs(d) >= tol1) { u=x+d; fu=(*f)(u); } else { u=x+SIGN(tol1,d); fu=(*f)(u); if (fu > fx) { If the minimum step in the downhill direction takes us uphill, then we are done. *xmin=x; return fx;