%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%  Purpose
%
%       Driver for the minimization of a nonlinear least squares 
%       function.
%
%       See Example 2.1 in J. Nocedal and S. J. Wright, 
%       "Numerical Optimization", Springer-Verlag, 1999.
%
%  Functions called: func.m, gradient.m, hessian.m
%
%
%  Written by L. N. Vicente.
%  February 21, 2000.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Begin.
%
% Clear work space and screen.
  clf;
  clear;
%
% Set formats.
  format compact;
  format long;
%
% Initializations.
%
% fid   file identifier for output file. No output if fid < 0.
%       Use fid = 1 for standard output (usually screen).
  fid = 1;
%
  t = [ 1 2 3 4 5 6 7 8 9 10]';
% t = [ 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2.0]';
  x = [ 1 1 1 1 1 1 ]';
  r = x(1) + x(2)*exp(-(x(3)-t).^2/x(4)) + x(5)*cos(x(6).*t);
  y = r + 10^(-5);
%
% Plotting the observations.
  plot(t,y,'.')
  hold on;
%
% Initial points.
% x0 = x + 10^(0)*ones(6,1);
  x0 = x + 10^(-1)*ones(6,1);
%
  if 0
%
% Newton's method with a backtracking line search.
  c1    = 10^(-4);
  rho   = 0.5;
  iter  = 0; ifunc=0; ig=0; iH=0;
  x     = x0;
  [f]   = func(x,t,y);
  ifunc = ifunc+1;
  [g]   = gradient(x,t,y);
  ig    = ig+1;
  [H]   = hessian(x,t,y);
  iH    = iH+1;
  while ( norm(g) > 10^(-8) & cond(H) < 10^(10) & iter <= 99 )
          [H]       = hessian(x,t,y);
          iH        = iH+1;
          p         = - H\g;
          alpha     = 1;
          xalphap   = x + alpha*p;
          [falphap] = func(xalphap,t,y);
          ifunc     = ifunc+1;
          linearf   = f+c1*alpha*g'*p;
          while ( falphap > linearf )
%                 Backtrack.
                  alpha     = rho*alpha;
                  xalphap   = x + alpha*p;
                  [falphap] = func(xalphap,t,y);
                  ifunc     = ifunc+1;
                  linearf   = f+c1*alpha*g'*p;
          end
          x    = xalphap;
          f    = falphap;
          [g]  = gradient(x,t,y);
          ig   = ig+1;
          iter = iter + 1;
  end
%
  r      = x(1) + x(2)*exp(-(x(3)-t).^2/x(4)) + x(5)*cos(x(6).*t);
  normr2 = norm(y-r);
  normg  = norm(g);
  condH  = cond(H);
%
% Print numerical results.
  fprintf( fid,'\nNewton''s method with a backtracking line search \n');
  fprintf( fid,'\n|  iter  | fevals | gevals | Hevals |   norm(g)  |   cond(H)  |   normr2   | \n');
  format_str = '| %6d | %6d | %6d | %6d | %9.4e | %9.4e | %9.4e | \n \n';
  fprintf( fid, format_str, iter, ifunc, ig, iH, normg, condH, normr2);
%
% Plotting the values obtained by the model after minimization.
  plot(t,r,'x')
%
  else
%
% Newton's method.
  iter = 0; ifunc=0; ig=0; iH=0;
  x    = x0;
  [g]  = gradient(x,t,y);
  ig   = ig+1;
  [H]  = hessian(x,t,y);
  iH   = iH+1;
  while ( norm(g) > 10^(-8) & cond(H) < 10^(10) & iter <= 99 )
          x    = x - H\g;
          [g]  = gradient(x,t,y);
          ig   = ig+1;
          [H]  = hessian(x,t,y);
          iH   = iH+1;
          iter = iter + 1;
  end
%
  r      = x(1) + x(2)*exp(-(x(3)-t).^2/x(4)) + x(5)*cos(x(6).*t);
  normr2 = norm(y-r);
  normg  = norm(g);
  condH  = cond(H);
%
% Print numerical results.
  fprintf( fid,'\nNewton''s method \n');
  fprintf( fid,'\n|  iter  | fevals | gevals | Hevals |   norm(g)  |   cond(H)  |   normr2   | \n');
  format_str = '| %6d | %6d | %6d | %6d | %9.4e | %9.4e | %9.4e | \n \n';
  fprintf( fid, format_str, iter, ifunc, ig, iH, normg, condH, normr2);
%
% Plotting the values obtained by the model after minimization.
  plot(t,r,'x')
%
  end
%
% End of driver for the minimization of the nonlinear least squares function.
