1 ## Copyright (C) 2002 Etienne Grossmann <etienne@egdn.net>
2 ## Copyright (C) 2004 Michael Creel <michael.creel@uab.es>
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
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
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/>.
17 ## test_min_2 - Test that bfgs works
19 ## Defines some simple functions and verifies that calling
21 ## bfgs on them returns the correct minimum.
23 ## Sets 'ok' to 1 if success, 0 otherwise
25 if ! exist ("optim_func"), optim_func = "bfgsmin"; end
29 if ! exist ("verbose"), verbose = 0; end
32 R = 20; # must have R >= P
36 ## Make test_min_2 reproducible by using fixed obsmat
37 ## obsmat = randn(R,P) ;
39 obsmat(sub2ind([R,P],1:R,1+rem(0:R-1,P))) = 1:R;
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);
50 obses = obsmat*truep ;
56 v = mean ((obses - obsmat*x).^2) + 1 ;
63 er = -obses + obsmat*x ;
64 dv = 2*er'*obsmat / rows(obses) ;
65 ## dv = 2*er'*obsmat ;
70 ## Returns the cputime since last call to 'mytic'.
73 persistent last_mytic = 0 ;
81 printf ("\n Testing %s on a quadratic problem\n\n", optim_func);
83 printf ([" Set 'optim_func' to the name of the optimization\n",\
84 " function you want to test (must have same synopsis\n",\
87 printf (" Nparams = P = %i, Nobses = R = %i\n",P,R);
95 if strcmp(optim_func,"bfgsmin")
100 ## [xlev,vlev,nlev] = feval(optim_func, "ff", "dff", xinit) ;
101 [xlev,vlev,nlev] = feval(optim_func, "ff", xinit2, ctl) ;
105 if max (abs(xlev-truep)) > 1e-4,
107 printf ("Error is too big : %8.3g\n", max (abs (xlev-truep)));
110 elseif verbose, printf ("ok 1\n");
114 printf (" Costs : init=%8.3g, final=%8.3g, best=%8.3g\n",\
115 ff(xinit), vlev, ff(truep));
118 printf ( " time : %8.3g\n",tlev);
121 printf ( "All tests ok (there's just one test)\n");