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).

test_mtron
%
% Test Routine for the mtron gateway function.
%
% F(u) is the finite difference discretization of the p-Laplacian
%           \int_0^1  |u|^p - u  dx
% test_mtron minimizes F(u) subject to the constraints u(0)=u(1)=0
%
%

function test_mtron

% Number of nodes for the finite difference discretization of the
% p-Laplacian functional
N = 31;
% sobolev index for p-Laplacian
p = 6;

% TRON parameters
delta = 0.1;            % initial trust region radius
frtol = 0.0;            % relative and absolute tolerance for termination
fatol = 0.0;            % due to change of F
cgtol = 1e-1;           % tolerance for the CG iteration. This value
% is suggested by the authors of TRON
fmin = -1e30;           % lower bound for F (usually it is best to make
% it a trivial lower bound

gtol = 1e-8;            % gtol is a tolerance setting not used by TRON
% but in this routine. The iteration terminates
% if the norm of the gradient is less then
% gtol.

% generate the initial condition
x = linspace(0, 1, N)';
u0 = sin(pi*x);

% the bound constraints. On interior nodes, they are (essentially)
% -\infty \leq u_j \leq \infty. on the boundary nodes the constraint
% can be used to define the boundary condition u_0 = u_N = 0
ul = -1e30 * ones(N, 1);
uu = -ul;
ul(1) = 0.0; uu(1) = 0.0;
ul(N) = 0.0; uu(N) = 0.0;

% plot initial condition
plot(x, u0); hold on;

% define the objective function
% note that plap implements the entire (sparse) hessian matrix, even
% though mtron requires the striclty lower triangular part and the
% diagonal separately. This is done to demonstrate how to convert
% efficiently.
ofun = @(uu)(plap(uu, p));

% evaluate functional at initial condition and output interation info
[F, G, H] = ofun(u0);
normG = norm(G(2:(N-1)), inf);
disp([' F = ', num2str(F), '    ||G|| = ', num2str(normG)]);

% initialize MTRON
% the hessian is always passed in two parts, the strictly lower
% triangular part is passed as a sparse matrix and the diagonal
% as a full vector.
[task, u, info] = mtron('INIT', u0, ul, uu, fmin, ...
F, G, tril(H, -1), full(diag(H)), ...
delta, frtol, fatol, cgtol);

% initialize iteration counter
nit = 0;

% start actual loop. loop as long as the task returned by mtron
% is neither 'CONVERGENCE' nor 'WARNING'.
while ( (task(1) ~= 'C') && (task(1) ~= 'W') )

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

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

% if the task is 'NEWX', we plot the current iterate, evaluate the
% gradient and the function and display the iteration info.
elseif (task(1) == 'N')
nit = nit + 1;
plot(x, u);
[F, G] = ofun(u);
% (don't take first and last entry in norm since those are
%  restricted by the bound constraints.)
normG = norm(G(2:(N-1)), inf);
disp([num2str(nit), ': F = ', num2str(F), ...
'    ||G|| = ', num2str(normG)]);
title([' ||G|| = ', num2str(normG)]);
drawnow;

% if the norm of G is less than a given tolerance
if (normG < 1e-6)
break;
end

% call TRON again. Since the current task is 'NEWX', no new
% information has to be passed on.
[task, u, info] = mtron('STEP');
end

end

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

% output termination reason.
if (task(1) == 'C')
disp('Termination due to change in function value less then tol.');