octave-optim 1.3.0-1 / usr / share / octave / packages / optim-1.3.0 / minimize.m

## Copyright (C) 2002 Etienne Grossmann <etienne@egdn.net>
##
## 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/>.

## [x,v,nev,...] = minimize (f,args,...) - Minimize f
##
## ARGUMENTS
## f    : string  : Name of function. Must return a real value
## args : list or : List of arguments to f (by default, minimize the first)
##        matrix  : f's only argument
##
## RETURNED VALUES
## x   : matrix  : Local minimum of f. Let's suppose x is M-by-N.
## v   : real    : Value of f in x0
## nev : integer : Number of function evaluations 
##     or 1 x 2  : Number of function and derivative evaluations (if
##                 derivatives are used)
## 
##
## Extra arguments are either a succession of option-value pairs or a single
## list or struct of option-value pairs (for unary options, the value in the
## struct is ignored).
## 
## OPTIONS : DERIVATIVES   Derivatives may be used if one of these options
## ---------------------   uesd. Otherwise, the Nelder-Mean (see
##                         nelder_mead_min) method is used.
## 
## 'd2f', d2f     : Name of a function that returns the value of f, of its
##                  1st and 2nd derivatives : [fx,dfx,d2fx] = feval (d2f, x)
##                  where fx is a real number, dfx is 1x(M*N) and d2fx is
##                  (M*N)x(M*N). A Newton-like method (d2_min) will be used.
##
## 'hess'         : Use [fx,dfx,d2fx] = leval (f, args) to compute 1st and
##                  2nd derivatives, and use a Newton-like method (d2_min).
##
## 'd2i', d2i     : Name of a function that returns the value of f, of its
##                  1st and pseudo-inverse of second derivatives : 
##                  [fx,dfx,id2fx] = feval (d2i, x) where fx is a real
##                  number, dfx is 1x(M*N) and d2ix is (M*N)x(M*N).
##                  A Newton-like method will be used (see d2_min).
##
## 'ihess'        : Use [fx,dfx,id2fx] = leval (f, args) to compute 1st
##                  derivative and the pseudo-inverse of 2nd derivatives,
##                  and use a Newton-like method (d2_min).
##
##            NOTE : df, d2f or d2i take the same arguments as f.
## 
## 'order', n     : Use derivatives of order n. If the n'th order derivative
##                  is not specified by 'df', 'd2f' or 'd2i', it will be
##                  computed numerically. Currently, only order 1 works.
## 
## 'ndiff'        : Use a variable metric method (bfgs) using numerical
##                  differentiation.
##
## OPTIONS : STOPPING CRITERIA  Default is to use 'tol'
## ---------------------------
## 'ftol', ftol   : Stop search when value doesn't improve, as tested by
##
##              ftol > Deltaf/max(|f(x)|,1)
##
##                 where Deltaf is the decrease in f observed in the last
##                 iteration.                                 Default=10*eps
##
## 'utol', utol   : Stop search when updates are small, as tested by
##
##              tol > max { dx(i)/max(|x(i)|,1) | i in 1..N }
##
##                 where  dx is the change in the x that occured in the last
##                 iteration.
##
## 'dtol',dtol    : Stop search when derivatives are small, as tested by
##
##              dtol > max { df(i)*max(|x(i)|,1)/max(v,1) | i in 1..N }
##
##                 where x is the current minimum, v is func(x) and df is
##                 the derivative of f in x. This option is ignored if
##                 derivatives are not used in optimization.
##
## MISC. OPTIONS
## -------------
## 'maxev', m     : Maximum number of function evaluations             <inf>
##
## 'narg' , narg  : Position of the minimized argument in args           <1>
## 'isz'  , step  : Initial step size (only for 0 and 1st order method)  <1>
##                  Should correspond to expected distance to minimum
## 'verbose'      : Display messages during execution
##
## 'backend'      : Instead of performing the minimization itself, return
##                  [backend, control], the name and control argument of the
##                  backend used by minimize(). Minimimzation can then be
##                  obtained without the overhead of minimize by calling, if
##                  a 0 or 1st order method is used :
##
##              [x,v,nev] = feval (backend, args, control)
##                   
##                  or, if a 2nd order method is used :
##
##              [x,v,nev] = feval (backend, control.d2f, args, control)

function [x,v,nev,varargout] = minimize (f,args,varargin)

## Oldies
##
## 'df' , df      : Name of a function that returns the derivatives of f
##                  in x : dfx = feval (df, x) where dfx is 1x(M*N). A
##                  variable metric method (see bfgs) will be used.
##
## 'jac'          : Use [fx, dfx] = leval(f, args) to compute derivatives
##                  and use a variable metric method (bfgs).
##


# ####################################################################
# Read the options ###################################################
# ####################################################################
# Options with a value
op1 = "ftol utol dtol df d2f d2i order narg maxev isz";
# Boolean options 
op0 = "verbose backend jac hess ihess ndiff" ;

default = struct ("backend",0,"verbose",0,...
		    "df","",  "df", "","d2f","","d2i","",  ...
		    "hess", 0,  "ihess", 0,  "jac", 0,"ndiff", 0,  ...
		    "ftol" ,nan, "utol",nan, "dtol", nan,...
		    "order",nan, "narg",nan, "maxev",nan,...
		    "isz",  nan);

if nargin == 3			# Accomodation to struct and list optional
  tmp = varargin{1};

  if isstruct (tmp)
    opls = {};
    for [v,k] = tmp		# Treat separately unary and binary opts
      if findstr ([" ",k," "],op0)
	opls (end+1) = {k};        # append k 
      else
	opls (end+[1:2]) = {k, v}; # append k and v 
      end
    end
  elseif iscell (tmp)
    opls = tmp;
  else
    opls = {tmp};
  end
else
  opls = varargin;
end
ops = read_options (opls,...
		    "op0",op0, "op1",op1, "default",default);

backend=ops.backend; verbose=ops.verbose; 
df=ops.df; d2f=ops.d2f; d2i=ops.d2i; 
hess=ops.hess; ihess=ops.ihess; jac=ops.jac; 
ftol=ops.ftol; utol=ops.utol; dtol=ops.dtol;
order=ops.order; narg=ops.narg; maxev=ops.maxev; 
isz=ops.isz; ndiff=ops.ndiff;

if length (df), error ("Option 'df' doesn't exist any more. Sorry.\n");end
if jac, error ("Option 'jac' doesn't exist any more. Sorry.\n");end

				# Basic coherence checks #############

ws = "";			# Warning string
es = "";			# Error string

				# Warn if more than 1 differential is given
if !!length (df) + !!length (d2f) + !!length (d2i) + jac + hess + ihess + ...
      ndiff > 1
				# Order of preference of 
  if length (d2i), ws = [ws,"d2i='",d2i,"', "]; end
  if length (d2f), ws = [ws,"d2f='",d2f,"', "]; end
  if length (df),  ws = [ws,"df='",df,"', "]; end
  if hess       ,  ws = [ws,"hess, "]; end
  if ihess      ,  ws = [ws,"ihess, "]; end
  if jac        ,  ws = [ws,"jac, "]; end
  if ndiff      ,  ws = [ws,"ndiff, "]; end
  ws = ws(1:length(ws)-2);
  ws = ["Options ",ws," were passed. Only one will be used\n"]
end

				# Check that enough args are passed to call
				# f(), unless backend is specified, in which
				# case I don't need to call f()
if ! isnan (narg) && ! backend
  if narg > 1
    es = [es,sprintf("narg=%i, but a single argument was passed\n",narg)];
  end
end

if length (ws), warn (ws); end
if length (es), error (es); end	# EOF Basic coherence checks #########


op = 0;				# Set if any option is passed and should be
				# passed to backend

if ! isnan (ftol)   , ctls.ftol    = ftol;  op = 1; end
if ! isnan (utol)   , ctls.utol    = utol;  op = 1; end
if ! isnan (dtol)   , ctls.dtol    = dtol;  op = 1; end
if ! isnan (maxev)  , ctls.maxev   = maxev; op = 1; end
if ! isnan (narg)   , ctls.narg    = narg;  op = 1; end
if ! isnan (isz)    , ctls.isz     = isz;   op = 1; end
if         verbose  , ctls.verbose = 1;     op = 1; end

				# defaults That are used in this function :
if isnan (narg), narg = 1; end

				# ##########################################
				# Choose one optimization method ###########
				# Choose according to available derivatives 
if     ihess, d2f = f;  ctls.id2f = 1; op = 1;
elseif hess,  d2f = f;
end
  

if     length (d2i), method = "d2_min"; ctls.id2f = 1; op = 1; d2f = d2i;
elseif length (d2f), method = "d2_min";
### elseif length (df) , method = "bfgsmin"; ctls.df  = df; op = 1;
### elseif jac         , method = "bfgsmin"; ctls.jac = 1 ; op = 1;
  ## else                 method = "nelder_mead_min";
  ## end
				# Choose method because ndiff is passed ####
elseif ndiff       , method = "bfgsmin";

				# Choose method by specifying order ########
elseif ! isnan (order)

  if     order == 0, method = "nelder_mead_min";
  elseif order == 1
    method = "bfgsmin";

  elseif order == 2
    if ! (length (d2f) || length (d2i))
      error ("minimize(): 'order' is 2, but 2nd differential is missing\n");
    end
  else
    error ("minimize(): 'order' option only implemented for order<=2\n");
  end
else				# Default is nelder_mead_min
  method = "nelder_mead_min";
end				# EOF choose method ########################

if verbose
  printf ("minimize(): Using '%s' as back-end\n",method);
end

				# More checks ##############################
ws = "";
if !isnan (isz) && strcmp (method,"d2_min")
  ws = [ws,"option 'isz' is passed to method that doesn't use it"];
end
if length (ws), warn (ws); end
				# EOF More checks ##########################

if     strcmp (method, "d2_min"), all_args = {f, d2f, args};
elseif strcmp (method, "bfgsmin"),all_args = {f, args};
else                              all_args = {f, args};
end
				# Eventually add ctls to argument list
if op, all_args{end+1} = ctls; end

if ! backend			# Call the backend ###################
  if strcmp (method, "d2_min"),
    [x,v,nev,h] = d2_min(all_args{:});
				# Eventually return inverse of Hessian
    if nargout > 3, varargout{1} = h; vr_val_cnt=2; end 
  elseif strcmp (method, "bfgsmin")
    nev = nan;
    if !iscell(args), args = {args}; end
    if isnan (ftol), ftol = 1e-12; end # Use bfgsmin's defaults
    if isnan (utol), utol = 1e-6; end
    if isnan (dtol), dtol = 1e-5; end
    if isnan (maxev), maxev = inf; end
    [x, v, okcv] = bfgsmin (f, args, {maxev,verbose,1,narg,0,ftol,utol,dtol});
  else
    [x,v,nev] = feval (method, all_args{:});
  end

else				# Don't call backend, just return its name
				# and arguments. 

  x = method;
  if op, v = ctls; else v = []; end
end