octave_packages/optim-1.2.0/test_min_2.m

   1 ## Copyright (C) 2002 Etienne Grossmann <etienne@egdn.net>
   2 ## Copyright (C) 2004 Michael Creel <michael.creel@uab.es>
   3 ##
   4 ## This program is free software; you can redistribute it and/or modify it under
   5 ## the terms of the GNU General Public License as published by the Free Software
   6 ## Foundation; either version 3 of the License, or (at your option) any later
   7 ## version.
   8 ##
   9 ## This program is distributed in the hope that it will be useful, but WITHOUT
  10 ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  11 ## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
  12 ## details.
  13 ##
  14 ## You should have received a copy of the GNU General Public License along with
  15 ## this program; if not, see <http://www.gnu.org/licenses/>.
  16
  17 ## test_min_2                   - Test that bfgs works
  18 ##
  19 ## Defines some simple functions and verifies that calling
  20 ##
  21 ## bfgs on them returns the correct minimum.
  22 ##
  23 ## Sets 'ok' to 1 if success, 0 otherwise
  24
  25 if ! exist ("optim_func"), optim_func = "bfgsmin"; end
  26
  27 ok = 1;
  28
  29 if ! exist ("verbose"), verbose = 0; end
  30
  31 P = 15;
  32 R = 20;                 # must have R >= P
  33
  34
  35 global obsmat ;
  36 ## Make test_min_2 reproducible by using fixed obsmat
  37 ## obsmat = randn(R,P) ;
  38 obsmat = zeros (R,P);
  39 obsmat(sub2ind([R,P],1:R,1+rem(0:R-1,P))) = 1:R;
  40
  41 global truep ;
  42
  43 ## Make test_min_2 reproducible by using fixed starting point
  44 ## truep = randn(P,1) ;
  45 ## xinit = randn(P,1) ;
  46 truep = rem (1:P, P/4)';
  47 xinit = truep + 2*(1:P)'/(P);
  48
  49 global obses ;
  50 obses = obsmat*truep ;
  51
  52
  53 function v = ff(x)
  54   global obsmat;
  55   global obses;
  56   v = mean ((obses - obsmat*x).^2) + 1 ;
  57 endfunction
  58
  59
  60 function dv = dff(x)
  61   global obsmat;
  62   global obses;
  63   er = -obses + obsmat*x ;
  64   dv = 2*er'*obsmat / rows(obses) ;
  65   ## dv = 2*er'*obsmat ;
  66 endfunction
  67
  68 ##       dt = mytic()
  69 ##
  70 ## Returns the cputime since last call to 'mytic'.
  71
  72 function dt = mytic()
  73    persistent last_mytic = 0 ;
  74    [t,u,s] = cputime() ;
  75    dt = t - last_mytic ;
  76    last_mytic = t ;
  77 endfunction
  78
  79
  80 if verbose
  81   printf ("\n   Testing %s on a quadratic problem\n\n", optim_func);
  82
  83   printf (["     Set 'optim_func' to the name of the optimization\n",\
  84            "     function you want to test (must have same synopsis\n",\
  85            "     as 'bfgs')\n\n"]);
  86
  87   printf ("  Nparams = P = %i,  Nobses = R = %i\n",P,R);
  88   fflush (stdout);
  89 end
  90
  91 ctl.df = "dff";
  92 ctl.ftol = eps;
  93 ctl.dtol = 1e-7;
  94 mytic() ;
  95 if strcmp(optim_func,"bfgsmin")
  96         ctl = {-1,2,1,1};
  97         xinit2 = {xinit};
  98 else xinit2 = xinit;
  99 endif
 100 ## [xlev,vlev,nlev] = feval(optim_func, "ff", "dff", xinit) ;
 101 [xlev,vlev,nlev] = feval(optim_func, "ff", xinit2, ctl) ;
 102 tlev = mytic() ;
 103
 104
 105 if max (abs(xlev-truep)) > 1e-4,
 106   if verbose
 107     printf ("Error is too big : %8.3g\n", max (abs (xlev-truep)));
 108   end
 109   ok = 0;
 110 elseif verbose,  printf ("ok 1\n");
 111 end
 112
 113 if verbose,
 114   printf ("  Costs :     init=%8.3g, final=%8.3g, best=%8.3g\n",\
 115           ff(xinit), vlev, ff(truep));
 116 end
 117 if verbose
 118     printf ( "   time : %8.3g\n",tlev);
 119 end
 120 if verbose && ok
 121   printf ( "All tests ok (there's just one test)\n");
 122 end
 123