diff options
Diffstat (limited to 'src')
| -rw-r--r-- | src/fit.c | 58 |
1 files changed, 58 insertions, 0 deletions
@@ -5691,6 +5691,64 @@ close: memcpy(xopt,x,sizeof(x)); } + if (time_elapsed >= maxtime) goto end; + + /* Now, minimize with SBPLX. + * + * The reason we minimize with a different algorithm is because occasionally + * the BOBYQA minimizer seems to get stuck somewhere which doesn't appear to + * be a local minima (this was noticed when fitting 300 MeV muons). + * + * Although SBPLX is *much* slower, the idea here is that if we really did + * already find a local minimum, then SBPLX should be relatively fast. + * + * FIXME: Ideally we would just figure out and fix whatever is causing + * BOBYQA to get stuck which would make things faster and be much better + * since it would fix the initial multiple optimizations we do in order to + * determine which particle id goes in which direction. */ + nlopt_destroy(opt); + opt = nlopt_create(NLOPT_LN_SBPLX, 4+3*n); + nlopt_set_min_objective(opt,nlopt_nll2,&fpars); + nlopt_set_lower_bounds(opt, lb); + nlopt_set_upper_bounds(opt, ub); + nlopt_set_initial_step(opt, ss); + nlopt_set_ftol_abs(opt, 1e-5); + nlopt_set_xtol_rel(opt, 1e-2); + nlopt_set_maxeval(opt, 1000); + nlopt_set_maxtime(opt, maxtime-time_elapsed); + + memcpy(x,xopt,sizeof(x)); + + gettimeofday(&tv_start, NULL); + nlopt_optimize(opt,x,&fval); + gettimeofday(&tv_stop, NULL); + + time_elapsed = tv_stop.tv_sec - tv_start.tv_sec + (tv_stop.tv_usec - tv_start.tv_usec)/1e6; + + if (time_elapsed > maxtime) goto end; + + if (stop) goto stop; + + do { + *fmin = fval; + memcpy(xopt,x,sizeof(x)); + + nlopt_set_maxtime(opt, maxtime-time_elapsed); + + gettimeofday(&tv_start, NULL); + nlopt_optimize(opt,x,&fval); + gettimeofday(&tv_stop, NULL); + + time_elapsed += tv_stop.tv_sec - tv_start.tv_sec + (tv_stop.tv_usec - tv_start.tv_usec)/1e6; + + if (stop) goto stop; + } while (fval < *fmin && fabs(fval-*fmin) > 1e-2 && time_elapsed < maxtime); + + if (fval < *fmin) { + *fmin = fval; + memcpy(xopt,x,sizeof(x)); + } + end: free(result); |