Code covered by the BSD License

### Highlights from MTRON

• ...% ITRON finds a constrained minimum of a function of several variables
• plap(u, p)% PLAP implements a finite difference discretizaion of the p-Laplacian
• test_mtron% Test Routine for the mtron gateway function.
• make_mtron.m% make_mtron
• mtron.m% MTRON is a routine that interfaces the DTRON routine of the
• test_itron.m% TEST_ITRON is a script that demonstrates the use of the box-constraint
• View all files

# MTRON

### Christoph Ortner (view profile)

02 May 2007 (Updated )

MTRON is a Matlab wrapper for the Fortran software TRON (a large scale trust region Newton method).

...
%
% ITRON finds a constrained minimum of a function of several variables
%
%   ITRON attempts to solve problems of the form
%        min F(x)    subject to xl <= x <= xu  (box constraints)
%
%   Calling convention (required arguments):
%
%   [x, fval, exitflag, output] = itron(fun, x0, xl, xu, fmin)
%
%   The following additional arguments are optional. Replace by [] if
%   default values should be used.
%
%   [x, fval, exitflag, output] = itron(fun, x0, xl, xu, fmin, ...
%                                 delta, gtol, maxnit, frtol, fatol, ...
%                                 cgtol, display, output_fcn)
%
%   Required Arguments:
%      fun    : function handle. Takes one argument (x) and returns
%               [F(x), G(x), H(x)], where G(x) \in R^n is the gradient
%               of F at x and H(x) is sparse and is the hessian of F at x.
%               (note that the hessian MUST be supplied)
%      x0     :  initial guess, n-column vector
%      xl, xu : lower and upper bound constraints for x. If you don't
%               have bounds, just use +- 1e300.
%      fmin   : lower bound for F(x). If TRON/mtron finds an iterate for
%               which F(x) < fmin, it terminates with a WARNING message.
%               if you don't have a lower bound, just use -1e300.
%
%   Optional Arguments:
%      delta  : initial trust region radius.
%               [Default : 1.0]
%      gtol   : at the end of each iteration, itron, determines which
%               variables are free (= constraints are inactive) and
%               computed the inf-norm of the gradient for those variables
%               only. If the norm is less than gtol, itron terminates
%               [Default : 1e-8]
%      maxnit : Maximum number of iterations
%               [Default : 100]
%      frtol  : TRON terminates if the relative change in the objective
%               function is less than options.frtol.
%               Default : 0.0
%      fatol  : TRON terminates if the absolute change in the objective
%               function value is less than options.fatol
%               Default : 0.0
%      cgtol  : termination tolerance for PCG iteration.
%               Default : 1e-1 (as suggested by authors of TRON)
%     display : Sets the display level: 0 = no display, 1 = display only
%               after termination, 2 = at every iteration, 3 = full detail
%               Default is 2.
%  output_fcn : set to a function which can display output of the
%               optimization process at each iteration.
%
%   Output Arguments:
%           x : final iterate. This is always returned, even if ITRON quits
%               with and error. The user should therefore always call with
%               at least three output arguments.
%        fval : objective function value at final iterate
%    exitflag : itron was succesful if exitflag > 0, unseccesful otherwise.
%               0 : unknown error; 1 : termination as gradient tolerance
%               was reached; 2 : termination as change in objective
%               function was smaller than a tolerance; -1 : maximum number
%               of iterations reached; -2 : TRON issued a warning.
%               optimization process. currently only output.message is
%               implemented
%
%

%
% by Christoph Ortner
% First Version: 28/11/2006
% Changes:
%   01/12/2006
%     - extended DISPLAY options
%     - print delta at each iteration (still a bug it seems)
%   19/03/2007
%     - improved computation of norm of gradient using the indfree array
%   23/04/2007
%     - trust region radius is now correctly displayed
%     - new parameters: display and output_fcn.
%   02/05/2007
%     - display number of CG iterations
%     - better formatting of output
%

%
% TODO: - xtol,
%       - first computation of normG is invalid!!!!!
%       - statistics : number of F, G, H evaluations
%

function [x, fval, exitflag, output, G, H] = ...
itron(fun, x0, xl, xu, fmin, delta, ...
gtol, maxnit, frtol, fatol, cgtol, ...
display, output_fcn)

% 0 : no display
% 1 : final display
% 2 : at every iteration
% 3 : also  display every message.
DEFAULT_DISPLAY = 2;

%% plotting frequence
PLOTFREQ = 20;
plotcnt = 0;

%
% supply optional variables
%
if (nargin < 5)
error('mtron:itron:inputmin', 'itron : At least 5 input variables required.');
end
if (nargin > 13)
error('mtron:itron:inputmax', 'itron : At most 11 input variables.');
end

% make all variables which were not passed empty arrays
if (nargin < 6), delta = []; end
if (nargin < 7), gtol = []; end
if (nargin < 8), maxnit = []; end
if (nargin < 9), frtol = []; end
if (nargin < 10), fatol = []; end
if (nargin < 11), cgtol = []; end
if (nargin < 12), display = []; end
if (nargin < 13), output_fcn = []; end

% give default values to all variables which are empty arrays (i.e.
% either not passed, or passed as empty arrays
if isempty(delta), delta = 1.0; end
if isempty(gtol), gtol = 1e-8; end
if isempty(maxnit), maxnit = 100; end
if isempty(frtol), frtol = 1e-12; end
if isempty(fatol), fatol =  1e-12; end
if isempty(cgtol), cgtol = 1e-1; end
if isempty(display), display = DEFAULT_DISPLAY; end

% set the actual display variable
DISPLAY = display;
% set iteration counter to 0.
nit = 0;
% set exitflag to 0
exitflag = 0;
% set the counter for PCG iterations.
ncgit = 0;

% evaluate functional at initial condition and output interation info
[F, G, H] = fun(x0);
normG = norm(G, inf); % compute_norm(G);
if (DISPLAY >= 2)
output = sprintf('    n          F            ||G||           delta       #PCG');
disp(output);
output = sprintf('--------------------------------------------------------------');
disp(output);
output = sprintf('%6d    %e   %e   %e   %4d', nit, F, normG, delta, 0);
disp(output);
end

% initialize MTRON
[task, x, info] = mtron('INIT', x0, xl, xu, fmin, ...
F, G, tril(H, -1), full(diag(H)), ...
delta, frtol, fatol, cgtol);

% start actual loop. loop as long as the task returned by mtron
% is neither 'CONVERGENCE' nor 'WARNING'.

if (DISPLAY >= 3)
end

% if the task is 'F', it means we have to evaluate F(u) and return it
F = fun(x);
[task, x, info] = mtron('STEP', F);
delta = info.delta;

% if the task is 'GH' we have to evaluate the gradient and the
% hessian and return them.
[F, G, H] = fun(x);
[task, x, info] = mtron('STEP', G, tril(H, -1), full(diag(H)));
delta = info.delta;

% if the task is 'NEWX', then mtron just returned a new iterate,
% i.e. x changed
% increase iteration counter
nit = nit + 1;
% evaluate F, G
[F, G] = fun(x);
normG = compute_norm(G);
% output iteration info
if (DISPLAY >= 2)
output = sprintf('%6d    %e   %e   %e   %4d', nit, F, normG, ...
info.delta, info.ncgit - ncgit);
disp(output);
end
% update number of pCG iterations
ncgit = info.ncgit;

% if the norm of G is less than a given tolerance
if (normG < gtol)
exitflag = 1;
break;
end

% pass the data to an external routine to plot the progress
if ~isempty(output_fcn)
plotcnt = plotcnt + 1;
if (plotcnt >= PLOTFREQ)
output_fcn(x, F);
plotcnt = 0;
end
end

% call TRON again. Since the current task is 'NEWX', no new
% information has to be passed on.
delta = info.delta;
end

% check iteration count
if (nit >= maxnit)
exitflag = -1;
break;
end

end

% cleanup the memory used by mtron.
mtron('CLEAN');

% delete the output variable which we now need to use for
% something else.
clear output

% write F, G and H output.
fval = F;
if (nargout >= 6)
[F, G, H] = fun(x);
end

% output termination reason.
% 1. ||G|| small
if (exitflag == 1)
if (DISPLAY >= 1)
disp('itron: Termination due to ||G|| less than tolerance.');
end
output.message = '||G||_inf less than tolerance.';
% 2. maxnit reached
elseif (exitflag == -1)
if (DISPLAY >= 1)
disp('itron: Maximum number of iterations reached.');
end
output.message = 'Maximum number of iterations reached.';
% 3. broke loop through task message:
else
% 3.1  change in function value less than tol.
if (DISPLAY >= 1)
disp('itron: Change in function value less then tol.');
end
exitflag = 2;
% 3.2 warning message.
if (DISPLAY >= 1)
disp('itron: Termination of TRON due to error.');
end
exitflag = -2;
end
end

%
% computes the norm over all indices where
% x does not satisfy one of the bound constraints.
%
function normG = compute_norm(G)

% get a list of indices which are free.
indfree = mtron('FREE');
indfree = indfree(find(indfree));
% compute norm over those indices only.
normG = norm(G(indfree), inf);