Mercurial > hg > octave-lyh

comparison scripts/optimization/fzero.m @ 17182:c3c1ebfaa7dc

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
maint: Use common indentation for switch statement. * scripts/general/interp1.m, scripts/geometry/delaunay.m, scripts/help/__unimplemented__.m, scripts/image/cmunique.m, scripts/miscellaneous/edit.m, scripts/optimization/fzero.m, scripts/optimization/sqp.m, scripts/plot/__gnuplot_drawnow__.m, scripts/plot/hidden.m, scripts/plot/legend.m, scripts/plot/print.m, scripts/plot/printd.m, scripts/plot/private/__contour__.m, scripts/plot/private/__fltk_print__.m, scripts/plot/private/__gnuplot_print__.m, scripts/plot/private/__go_draw_axes__.m, scripts/plot/private/__print_parse_opts__.m, scripts/signal/periodogram.m, scripts/sparse/bicg.m, test/slice.tst, test/switch.tst: Use common indentation for switch statement.
author Rik <rik@octave.org>
date Sun, 04 Aug 2013 15:11:34 -0700
parents e7a059a9a644
children bc924baa2c4e
comparison
equal deleted inserted replaced
17181:3a23cbde59d5 17182:c3c1ebfaa7dc
194 d = e = u; 194 d = e = u;
195 fd = fe = fu; 195 fd = fe = fu;
196 mba = mu*(b - a); 196 mba = mu*(b - a);
197 while (niter < maxiter && nfev < maxfev) 197 while (niter < maxiter && nfev < maxfev)
198 switch (itype) 198 switch (itype)
199 case 1 199 case 1
200 ## The initial test. 200 ## The initial test.
201 if (b - a <= 2*(2 * abs (u) * eps + tolx)) 201 if (b - a <= 2*(2 * abs (u) * eps + tolx))
202 x = u; fval = fu; 202 x = u; fval = fu;
203 info = 1; 203 info = 1;
204 break; 204 break;
205 endif 205 endif
206 if (abs (fa) <= 1e3*abs (fb) && abs (fb) <= 1e3*abs (fa)) 206 if (abs (fa) <= 1e3*abs (fb) && abs (fb) <= 1e3*abs (fa))
207 ## Secant step. 207 ## Secant step.
208 c = u - (a - b) / (fa - fb) * fu; 208 c = u - (a - b) / (fa - fb) * fu;
209 else 209 else
210 ## Bisection step.
211 c = 0.5*(a + b);
212 endif
213 d = u; fd = fu;
214 itype = 5;
215 case {2, 3}
216 l = length (unique ([fa, fb, fd, fe]));
217 if (l == 4)
218 ## Inverse cubic interpolation.
219 q11 = (d - e) * fd / (fe - fd);
220 q21 = (b - d) * fb / (fd - fb);
221 q31 = (a - b) * fa / (fb - fa);
222 d21 = (b - d) * fd / (fd - fb);
223 d31 = (a - b) * fb / (fb - fa);
224 q22 = (d21 - q11) * fb / (fe - fb);
225 q32 = (d31 - q21) * fa / (fd - fa);
226 d32 = (d31 - q21) * fd / (fd - fa);
227 q33 = (d32 - q22) * fa / (fe - fa);
228 c = a + q31 + q32 + q33;
229 endif
230 if (l < 4 || sign (c - a) * sign (c - b) > 0)
231 ## Quadratic interpolation + newton.
232 a0 = fa;
233 a1 = (fb - fa)/(b - a);
234 a2 = ((fd - fb)/(d - b) - a1) / (d - a);
235 ## Modification 1: this is simpler and does not seem to be worse.
236 c = a - a0/a1;
237 if (a2 != 0)
238 c = a - a0/a1;
239 for i = 1:itype
240 pc = a0 + (a1 + a2*(c - b))*(c - a);
241 pdc = a1 + a2*(2*c - a - b);
242 if (pdc == 0)
243 c = a - a0/a1;
244 break;
245 endif
246 c -= pc/pdc;
247 endfor
248 endif
249 endif
250 itype += 1;
251 case 4
252 ## Double secant step.
253 c = u - 2*(b - a)/(fb - fa)*fu;
254 ## Bisect if too far.
255 if (abs (c - u) > 0.5*(b - a))
256 c = 0.5 * (b + a);
257 endif
258 itype = 5;
259 case 5
210 ## Bisection step. 260 ## Bisection step.
211 c = 0.5*(a + b);
212 endif
213 d = u; fd = fu;
214 itype = 5;
215 case {2, 3}
216 l = length (unique ([fa, fb, fd, fe]));
217 if (l == 4)
218 ## Inverse cubic interpolation.
219 q11 = (d - e) * fd / (fe - fd);
220 q21 = (b - d) * fb / (fd - fb);
221 q31 = (a - b) * fa / (fb - fa);
222 d21 = (b - d) * fd / (fd - fb);
223 d31 = (a - b) * fb / (fb - fa);
224 q22 = (d21 - q11) * fb / (fe - fb);
225 q32 = (d31 - q21) * fa / (fd - fa);
226 d32 = (d31 - q21) * fd / (fd - fa);
227 q33 = (d32 - q22) * fa / (fe - fa);
228 c = a + q31 + q32 + q33;
229 endif
230 if (l < 4 || sign (c - a) * sign (c - b) > 0)
231 ## Quadratic interpolation + newton.
232 a0 = fa;
233 a1 = (fb - fa)/(b - a);
234 a2 = ((fd - fb)/(d - b) - a1) / (d - a);
235 ## Modification 1: this is simpler and does not seem to be worse.
236 c = a - a0/a1;
237 if (a2 != 0)
238 c = a - a0/a1;
239 for i = 1:itype
240 pc = a0 + (a1 + a2*(c - b))*(c - a);
241 pdc = a1 + a2*(2*c - a - b);
242 if (pdc == 0)
243 c = a - a0/a1;
244 break;
245 endif
246 c -= pc/pdc;
247 endfor
248 endif
249 endif
250 itype += 1;
251 case 4
252 ## Double secant step.
253 c = u - 2*(b - a)/(fb - fa)*fu;
254 ## Bisect if too far.
255 if (abs (c - u) > 0.5*(b - a))
256 c = 0.5 * (b + a); 261 c = 0.5 * (b + a);
257 endif 262 itype = 2;
258 itype = 5;
259 case 5
260 ## Bisection step.
261 c = 0.5 * (b + a);
262 itype = 2;
263 endswitch 263 endswitch
264 264
265 ## Don't let c come too close to a or b. 265 ## Don't let c come too close to a or b.
266 delta = 2*0.7*(2 * abs (u) * eps + tolx); 266 delta = 2*0.7*(2 * abs (u) * eps + tolx);
267 if ((b - a) <= 2*delta) 267 if ((b - a) <= 2*delta)