Mercurial > mxe-octave
changeset 4942:805494784837
of-optim: Cherry-picked upstream fix (bug #55325).
* src/of-optim-2-deprecated.patch: new file
* dist-files.mk: add ref to of-optim-2-deprecated.patch
author | Markus Mützel <markus.muetzel@gmx.de> |
---|---|
date | Sun, 20 Jan 2019 17:59:50 +0100 |
parents | 1fa65a9facb7 |
children | c2953d95fbc3 |
files | dist-files.mk src/of-optim-2-deprecated.patch |
diffstat | 2 files changed, 574 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- a/dist-files.mk Sat Jan 26 10:19:17 2019 -0500 +++ b/dist-files.mk Sun Jan 20 17:59:50 2019 +0100 @@ -515,6 +515,7 @@ of-odepkg-3-deprecated.patch \ of-odepkg.mk \ of-optim-1-fixes.patch \ + of-optim-2-deprecated.patch \ of-optim.mk \ of-optiminterp-1-dev-fixes.patch \ of-optiminterp.mk \
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/src/of-optim-2-deprecated.patch Sun Jan 20 17:59:50 2019 +0100 @@ -0,0 +1,573 @@ +# HG changeset patch +# User Olaf Till <i7tiol@t-online.de> +# Date 1547218871 -3600 +# Fri Jan 11 16:01:11 2019 +0100 +# Node ID 5fcbf411d5c2a4b076e3ababeb6126125243e49d +# Parent 158b61118c7557f549656bf76ddc144c70bd3aa5 +remove deprecated samin.cc, still available as backend to nonlin_min + +diff -r 158b61118c75 -r 5fcbf411d5c2 src/Makefile.in +--- a/src/Makefile.in Fri Jan 11 15:29:45 2019 +0100 ++++ b/src/Makefile.in Fri Jan 11 16:01:11 2019 +0100 +@@ -46,7 +46,7 @@ + TEXIFILE := $(addsuffix .texi,$(basename $(INFOFILE))) + TXIFILE := $(addsuffix .txi,$(basename $(INFOFILE))) + HTMLDIR := ../doc/html/ +-OCTFILES := __bfgsmin.oct numgradient.oct numhessian.oct samin.oct __disna_optim__.oct ++OCTFILES := __bfgsmin.oct numgradient.oct numhessian.oct __disna_optim__.oct + + OCTSOURCEFILES := $(addsuffix .cc,$(basename $(OCTFILES))) + DSFILES := $(addsuffix .docstrings,$(OCTSOURCEFILES)) +diff -r 158b61118c75 -r 5fcbf411d5c2 src/samin.cc +--- a/src/samin.cc Fri Jan 11 15:29:45 2019 +0100 ++++ /dev/null Thu Jan 01 00:00:00 1970 +0000 +@@ -1,549 +0,0 @@ +-// Copyright (C) 2004, 2006 Michael Creel <michael.creel@uab.es> +-// +-// This program is free software; you can redistribute it and/or modify it under +-// the terms of the GNU General Public License as published by the Free Software +-// Foundation; either version 3 of the License, or (at your option) any later +-// version. +-// +-// This program is distributed in the hope that it will be useful, but WITHOUT +-// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +-// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more +-// details. +-// +-// You should have received a copy of the GNU General Public License along with +-// this program; if not, see <http://www.gnu.org/licenses/>. +- +-// References: +-// +-// The code follows the article +-// Goffe, William L. (1996) "SIMANN: A Global Optimization Algorithm +-// using Simulated Annealing " Studies in Nonlinear Dynamics & Econometrics +-// Oct96, Vol. 1 Issue 3. +-// +-// The code uses the same names for control variables, +-// for the most part. A notable difference is that the initial +-// temperature is found automatically to ensure that the active +-// bounds when the temperature begins to reduce cover the entire +-// parameter space (defined as a n-dimensional +-// rectangle that is the Cartesian product of the +-// (lb_i, ub_i), i = 1,2,..n +-// +-// Also of note: +-// Corana et. al., (1987) "Minimizing Multimodal Functions of Continuous +-// Variables with the "Simulated Annealing" Algorithm", +-// ACM Transactions on Mathematical Software, V. 13, N. 3. +-// +-// Goffe, et. al. (1994) "Global Optimization of Statistical Functions +-// with Simulated Annealing", Journal of Econometrics, +-// V. 60, N. 1/2. +- +-#include <oct.h> +-#include <octave/parse.h> +-#include <octave/Cell.h> +-#include <octave/lo-mappers.h> +-#include <octave/oct-rand.h> +-#include <float.h> +- +-#include "error-helpers.h" +- +-// define argument checks +-static bool any_bad_argument(const octave_value_list& args) +-{ +- +- // objective function name is a string? +- if (!args(0).is_string()) +- { +- _p_error("samin: first argument must be string holding objective function name"); +- return true; +- } +- +- // are function arguments contained in a cell? +- if (!args(1).OV_ISCELL ()) +- { +- _p_error("samin: second argument must cell array of function arguments"); +- return true; +- } +- +- // is control a cell? +- Cell control; +- bool err; +- SET_ERR (control =args(2).cell_value(), err); +- if (err) +- { +- _p_error("samin: third argument must cell array of algorithm controls"); +- return true; +- } +- +- // does control have proper number of elements? +- if (!(control.numel () == 11)) +- { +- _p_error("samin: third argument must be a cell array with 11 elements"); +- return true; +- } +- +- // now check type of each element of control +- if (!(control(0).is_real_matrix()) && !(control(0).is_real_scalar())) +- { +- _p_error("samin: 1st element of controls must be LB: a vector of lower bounds"); +- return true; +- } +- +- if ((control(0).is_real_matrix()) && (control(0).columns() != 1)) +- { +- _p_error("samin: 1st element of controls must be LB: a vector of lower bounds"); +- return true; +- } +- +- if (!(control(1).is_real_matrix()) && !(control(1).is_real_scalar())) +- { +- _p_error("samin: 1st element of controls must be UB: a vector of lower bounds"); +- return true; +- } +- +- if ((control(1).is_real_matrix()) && (control(1).columns() != 1)) +- { +- _p_error("samin: 2nd element of controls must be UB: a vector of lower bounds"); +- return true; +- } +- +- int tmp; +- SET_ERR (tmp = control(2).int_value(), err); +- if (err || tmp < 1) +- { +- _p_error("samin: 3rd element of controls must be NT: positive integer\n\ +-loops per temperature reduction"); +- return true; +- } +- +- SET_ERR (tmp = control(3).int_value(), err); +- if (err || tmp < 1) +- { +- _p_error("samin: 4th element of controls must be NS: positive integer\n\ +-loops per stepsize adjustment"); +- return true; +- } +- +- double tmp2; +- SET_ERR (tmp2 = control(4).double_value(), err); +- if (err || tmp < 0) +- { +- _p_error("samin: 5th element of controls must be RT:\n\ +-temperature reduction factor, RT > 0"); +- return true; +- } +- +- SET_ERR (tmp2 = control(5).double_value(), err); +- if (err || tmp < 0) +- { +- _p_error("samin: 6th element of controls must be integer MAXEVALS > 0 "); +- return true; +- } +- +- SET_ERR (tmp = control(6).int_value(), err); +- if (err || tmp < 0) +- { +- _p_error("samin: 7th element of controls must be NEPS: positive integer\n\ +-number of final obj. values that must be within EPS of eachother "); +- return true; +- } +- +- SET_ERR (tmp2 = control(7).double_value(), err); +- if (err || tmp2 < 0) +- { +- _p_error("samin: 8th element of controls must must be FUNCTOL (> 0)\n\ +-used to compare the last NEPS obj values for convergence test"); +- return true; +- } +- +- SET_ERR (tmp2 = control(8).double_value(), err); +- if (err || tmp2 < 0) +- { +- _p_error("samin: 9th element of controls must must be PARAMTOL (> 0)\n\ +-used to compare the last NEPS obj values for convergence test"); +- return true; +- } +- +- SET_ERR (tmp = control(9).int_value(), err); +- if (err || tmp < 0 || tmp > 2) +- { +- _p_error("samin: 9th element of controls must be VERBOSITY (0, 1, or 2)"); +- return true; +- } +- +- SET_ERR (tmp = control(10).int_value(), err); +- if (err || tmp < 0) +- { +- _p_error("samin: 10th element of controls must be MINARG (integer)\n\ +- position of argument to minimize wrt"); +- return true; +- } +- +- // make sure that minarg points to an existing element +- if ((tmp > args(1).length ()) || (tmp < 1)) +- { +- _p_error("bfgsmin: 4th argument must be a positive integer that indicates \n\ +-which of the elements of the second argument is the one minimization is over"); +- return true; +- } +- +- return false; +-} +- +-//-------------- The annealing algorithm -------------- +-DEFUN_DLD(samin, args, , "samin: simulated annealing minimization of a function. See samin_example.m\n\ +-\n\ +-samin will be removed from a future version of the optim package.\n\ +-Equivalent functionality is now in the samin backend of nonlin_min.\n\ +-\n\ +-usage: [x, obj, convergence, details] = samin(\"f\", {args}, {control})\n\ +-\n\ +-Arguments:\n\ +-* \"f\": function name (string)\n\ +-* {args}: a cell array that holds all arguments of the function,\n\ +-* {control}: a cell array with 11 elements\n\ +- * LB - vector of lower bounds\n\ +- * UB - vector of upper bounds\n\ +- * nt - integer: # of iterations between temperature reductions\n\ +- * ns - integer: # of iterations between bounds adjustments\n\ +- * rt - (0 < rt <1): temperature reduction factor\n\ +- * maxevals - integer: limit on function evaluations\n\ +- * neps - integer: number of values final result is compared to\n\ +- * functol - (> 0): the required tolerance level for function value\n\ +- comparisons\n\ +- * paramtol - (> 0): the required tolerance level for parameters\n\ +- * verbosity - scalar: 0, 1, or 2.\n\ +- * 0 = no screen output\n\ +- * 1 = only final results to screen\n\ +- * 2 = summary every temperature change\n\ +- * minarg - integer: which of function args is minimization over?\n\ +-\n\ +-Returns:\n\ +-* x: the minimizer\n\ +-* obj: the value of f() at x\n\ +-* convergence:\n\ +- 0 if no convergence within maxevals function evaluations\n\ +- 1 if normal convergence to a point interior to the parameter space\n\ +- 2 if convergence to point very near bounds of parameter space\n\ +- (suggest re-running with looser bounds)\n\ +-* details: a px3 matrix. p is the number of times improvements were found.\n\ +- The columns record information at the time an improvement was found\n\ +- * first: cumulative number of function evaluations\n\ +- * second: temperature\n\ +- * third: function value\n\ +-\n\ +-Example: see samin_example\n\ +-") +-{ +- static bool warned = false; +- if (! warned) +- { +- warned = true; +- warning_with_id ("Octave:deprecated-function", +- "samin will be removed from a future version of the optim package, equivalent functionality is now in the samin backend of nonlin_min"); +- } +- +- int nargin = args.length(); +- if (!(nargin == 3)) +- { +- error("samin: you must supply 3 arguments"); +- return octave_value_list(); +- } +- +- // check the arguments +- if (any_bad_argument (args)) +- { +- error ("error in samin"); +- return octave_value_list(); +- } +- +- std::string obj_fn (args(0).string_value()); +- Cell f_args_cell = args(1).cell_value (); // args to obj fn come in as a cell to allow arbitrary number +- Cell control (args(2).cell_value()); +- +- octave_value_list f_args; +- octave_value_list f_return; // holder for feval returns +- +- int m, i, j, k, h, n, nacc, func_evals; +- int nup, nrej, nnew, ndown, lnobds; +- int converge, test, coverage_ok; +- +- // user provided controls +- const ColumnVector lb (control(0).column_vector_value()); +- const ColumnVector ub (control(1).column_vector_value()); +- const int nt (control(2).int_value()); +- const int ns (control(3).int_value()); +- const double rt (control(4).double_value()); +- const int maxevals (control(5).int_value()); +- const int neps (control(6).int_value()); +- const double functol (control(7).double_value()); +- const double paramtol (control(8).double_value()); +- const int verbosity (control(9).int_value()); +- const int minarg (control(10).int_value()); +- +- // type checking for minimization parameter done here, since we don't know minarg +- // until now +- if (!(f_args_cell(minarg - 1).is_real_matrix() || (f_args_cell(minarg - 1).is_real_scalar()))) +- { +- error("samin: minimization must be with respect to a column vector"); +- return octave_value_list(); +- } +- if ((f_args_cell(minarg - 1).is_real_matrix()) && (f_args_cell(minarg - 1).columns() != 1)) +- { +- error("samin: minimization must be with respect to a column vector"); +- return octave_value_list(); +- } +- +- double f, fp, p, fopt, rand_draw, ratio, t; +- +- Matrix details(1,3); // record function evaluations, temperatures and function values +- RowVector info(3); +- +- // copy cell contents over to octave_value_list to use feval() +- k = f_args_cell.numel (); +- f_args(k); // resize only once +- for (i = 0; i<k; i++) f_args(i) = f_args_cell(i); +- +- ColumnVector x = f_args(minarg - 1).column_vector_value(); +- ColumnVector bounds = ub - lb; +- n = x.rows(); +- ColumnVector xopt = x; +- ColumnVector xp(n); +- ColumnVector nacp(n); +- +- // Set initial values +- nacc = 0; // total accepted trials +- t = 1000.0; // temperature - will initially rise or fall to cover parameter space. Then it will fall +- converge = 0; // convergence indicator 0 (failure), 1 (normal success), or 2 (convergence but near bounds) +- coverage_ok = 0; // has parameter space been covered? When turns to 1, temperature starts to fall +- // most recent values, to compare to when checking convergend +- ColumnVector fstar(neps,1); +- fstar.fill(DBL_MAX); +- octave_rand::distribution("uniform"); // we'll be using draws from U(0,1) +- +- // check for out-of-bounds starting values +- for(i = 0; i < n; i++) +- { +- if(( x(i) > ub(i)) || (x(i) < lb(i))) +- { +- error("samin: initial parameter %d out of bounds", i+1); +- return octave_value_list(); +- } +- } +- +- // Initial obj_value +- f_return = OCTAVE__FEVAL (obj_fn, f_args); +- f = f_return(0).double_value(); +- fopt = f; // give it something to compare to +- func_evals = 0; // total function evaluations (limited by maxeval) +- details(0,0) = func_evals; +- details(0,1) = t; +- details(0,2) = fopt; +- +- // main loop, first increase temperature until parameter space covered, then reduce until convergence +- while(converge==0) +- { +- // statistics to report at each temp change, set back to zero +- nup = 0; +- nrej = 0; +- nnew = 0; +- ndown = 0; +- lnobds = 0; +- +- // repeat nt times then adjust temperature +- for(m = 0;m < nt;m++) +- { +- // repeat ns times, then adjust bounds +- for(j = 0;j < ns;j++) +- { +- // generate new point by taking last and adding a random value +- // to each of elements, in turn +- for(h = 0;h < n;h++) +- { +- // new Sept 2011, if bounds are same, skip the search for that vbl. Allows restrictions without complicated programming +- if (lb(h) != ub(h)) +- { +- xp = x; +- rand_draw = octave_rand::scalar(); +- xp(h) = x(h) + (2.0 * rand_draw - 1.0) * bounds(h); +- if ((xp(h) < lb(h)) || (xp(h) > ub(h))) +- { +- rand_draw = octave_rand::scalar(); // change 07-Nov-2007: avoid correlation with hitting bounds +- xp(h) = lb(h) + (ub(h) - lb(h)) * rand_draw; +- lnobds = lnobds + 1; +- } +- // Evaluate function at new point +- f_args(minarg - 1) = xp; +- f_return = OCTAVE__FEVAL (obj_fn, f_args); +- fp = f_return(0).double_value(); +- func_evals = func_evals + 1; +- // Accept the new point if the function value decreases +- if(fp <= f) +- { +- x = xp; +- f = fp; +- nacc = nacc + 1; // total number of acceptances +- nacp(h) = nacp(h) + 1; // acceptances for this parameter +- nup = nup + 1; +- // If lower than any other point, record as new optimum +- if(fp < fopt) +- { +- xopt = xp; +- fopt = fp; +- nnew = nnew + 1; +- info(0) = func_evals; +- info(1) = t; +- info(2) = fp; +- details = details.stack(info); +- } +- } +- // If the point is higher, use the Metropolis criteria to decide on +- // acceptance or rejection. +- else +- { +- p = exp(-(fp - f) / t); +- rand_draw = octave_rand::scalar(); +- if(rand_draw < p) +- { +- x = xp; +- f = fp; +- nacc = nacc + 1; +- nacp(h) = nacp(h) + 1; +- ndown = ndown + 1; +- } +- else nrej = nrej + 1; +- } +- } +- // If maxevals exceeded, terminate the algorithm +- if (func_evals >= maxevals) +- { +- if (verbosity >= 1) +- { +- printf("\n================================================\n"); +- printf("SAMIN results\n"); +- printf("NO CONVERGENCE: MAXEVALS exceeded\n"); +- printf("================================================\n"); +- printf("Convergence tolerances: Func. tol. %e Param. tol. %e\n", functol, paramtol); +- printf("Obj. fn. value %f\n\n", fopt); +- printf(" parameter search width\n"); +- for(i = 0; i < n; i++) printf("%20f%20f\n", xopt(i), bounds(i)); +- } +- f_return(3) = details; +- f_return(2) = 0; +- f_return(1) = fopt; +- f_return(0) = xopt; +- return octave_value_list(f_return); +- } +- } +- } +- // Adjust bounds so that approximately half of all evaluations are accepted +- test = 0; +- for(i = 0;i < n;i++) +- { +- if (lb(i) != ub(i)) +- { +- ratio = nacp(i) / ns; +- if(ratio > 0.6) bounds(i) = bounds(i) * (1.0 + 2.0 * (ratio - 0.6) / 0.4); +- else if (ratio < .4) +- bounds(i) = bounds(i) / (1.0 + 2.0 * ((0.4 - ratio) / 0.4)); +- // keep within initial bounds +- if(bounds(i) >= (ub(i) - lb(i))) +- { +- bounds(i) = ub(i) - lb(i); +- test = test + 1; +- } +- } +- else +- test = test + 1; // make sure coverage check passes for the fixed parameters +- } +- nacp.fill(0.0); +- // check if we cover parameter space. if we have yet to do so +- if (!coverage_ok) coverage_ok = (test == n); +- } +- // intermediate output, if desired +- if(verbosity == 2) +- { +- printf("\nsamin: intermediate results before next temperature change\n"); +- printf("\ntemperature %e", t); +- printf("\ncurrent best function value %f", fopt); +- printf("\ntotal evaluations so far %d", func_evals); +- printf("\ntotal moves since last temperature reduction %d", nup + ndown + nrej); +- printf("\ndownhill %d", nup); +- printf("\naccepted uphill %d", ndown); +- printf("\nrejected uphill %d", nrej); +- printf("\nout of bounds trials %d", lnobds); +- printf("\nnew minima this temperature %d", nnew); +- printf("\n\n parameter search width\n"); +- for(i = 0; i < n; i++) printf("%20f%20f\n", xopt(i), bounds(i)); +- printf("\n"); +- } +- // Check for convergence, if we have covered the parameter space +- if (coverage_ok) +- { +- // last value close enough to last neps values? +- fstar(0) = f; +- test = 0; +- for (i = 1; i < neps; i++) test = test + fabs(f - fstar(i)) > functol; +- test = (test > 0); // if different from zero, function conv. has failed +- // last value close enough to overall best? +- if (((fopt - f) <= functol) && (!test)) +- { +- // check for bound narrow enough for parameter convergence +- for (i = 0;i < n;i++) +- { +- if (bounds(i) > paramtol) +- { +- converge = 0; // no conv. if bounds too wide +- break; +- } +- else +- converge = 1; +- } +- } +- // check if optimal point is near boundary of parameter space, and change convergence message if so +- if (converge && lnobds > 0) converge = 2; +- // Like to see the final results? +- if (converge > 0) +- { +- if (verbosity >= 1) +- { +- printf("\n================================================\n"); +- printf("SAMIN results\n\n"); +- if (converge == 1) printf("==> Normal convergence <==\n\n"); +- if (converge == 2) +- { +- printf("==> WARNING <==: Last point satisfies convergence criteria,\n"); +- printf("but is near boundary of parameter space.\n"); +- printf("%d out of %d evaluations were out-of-bounds in the last round.\n", lnobds, (nup+ndown+nrej)); +- printf("Expand bounds and re-run, unless this is a constrained minimization.\n\n"); +- } +- printf("Convergence tolerances:\nFunction: %e\nParameters: %e\n", functol, paramtol); +- printf("\nObjective function value at minimum: %f\n\n", fopt); +- printf(" parameter search width\n"); +- for(i = 0; i < n; i++) printf("%20f%20f\n", xopt(i), bounds(i)); +- printf("================================================\n"); +- } +- f_return(3) = details; +- f_return(2) = converge; +- f_return(1) = fopt; +- f_return(0) = xopt; +- return f_return; // this breaks out, if we get here +- } +- // Reduce temperature, record current function value in the +- // list of last "neps" values, and loop again +- t = rt * t; +- for(i = neps-1; i > 0; i--) fstar(i) = fstar(i-1); +- f = fopt; +- x = xopt; +- } +- else +- { +- // coverage not ok - increase temperature quickly to expand search area, to cover parameter space +- t = t*t; +- for(i = neps-1; i > 0; i--) fstar(i) = fstar(i-1); +- f = fopt; +- x = xopt; +- } +- } +- // silence compiler warning +- return octave_value_list (); +-}