stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Efficiency changes. 

+ Use scalar product when possible, ie replace sum(x.*x) by x'*x,
+ Add parenthesis around scalar products, ie replace a*x'*x by a*(x'*x),
+ Removed global options_ in lnsrch1 (only used for solve_tolx).

Note that the scalar products could be replaced by a loop for small
problems. If the number of unknow is smaller than 70 the loop is
faster (10 times faster if the number of unknowns is less than 10).
time-shift
Stéphane Adjemian (Hermes) 2016-02-26 22:06:49 +01:00
parent 87955c61d0
commit 837c812c17
4 changed files with 16 additions and 24 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

24
matlab/lnsrch1.m
View File

@ -1,5 +1,5 @@
function [x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,varargin)
% function [x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,varargin)
function [x,f,fvec,check]=lnsrch1(xold, fold, g, p, stpmax, func, j1, j2, tolx, varargin)
% function [x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,tolx,varargin)
% Computes the optimal step by minimizing the residual sum of squares
%
% INPUTS
@ -11,6 +11,7 @@ function [x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,varargin)
% func: name of the function
% j1: equations index to be solved
% j2: unknowns index
% tolx: tolerance parameter
% varargin: list of arguments following j2
%
% OUTPUTS
@ -23,7 +24,7 @@ function [x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,varargin)
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2001-2010 Dynare Team
% Copyright (C) 2001-2016 Dynare Team
%
% This file is part of Dynare.
%
@ -40,15 +41,13 @@ function [x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,varargin)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global options_
alf = 1e-4 ;
tolx = options_.solve_tolx;
alam = 1;
x = xold;
nn = length(j2);
summ = sqrt(sum(p.*p)) ;
summ = sqrt(p'*p);
if ~isfinite(summ)
eq_number_string=[];
for ii=1:length(j1)-1
@ -73,7 +72,7 @@ if ~isfinite(summ)
end
if summ > stpmax
p=p.*stpmax/summ ;
p = p*stpmax/summ ;
end
slope = g'*p ;
@ -90,19 +89,16 @@ while 1
check = 1 ;
return
end
x(j2) = xold(j2) + (alam*p) ;
fvec = feval(func,x,varargin{:}) ;
fvec = fvec(j1);
f = 0.5*fvec'*fvec ;
f = 0.5*(fvec'*fvec) ;
if any(isnan(fvec))
alam = alam/2 ;
alam2 = alam ;
f2 = f ;
fold2 = fold ;
else
if f <= fold+alf*alam*slope
check = 0;
break ;
@ -124,15 +120,11 @@ while 1
else
tmplam = (-b+sqrt(disc))/(3*a) ;
end
end
if tmplam > 0.5*alam
tmplam = 0.5*alam;
end
end
alam2 = alam ;
f2 = f ;
fold2 = fold ;

6
matlab/solve1.m
View File

@ -23,7 +23,7 @@ function [x,check] = solve1(func,x,j1,j2,jacobian_flag,gstep,tolf,tolx,maxit,deb
% SPECIAL REQUIREMENTS
% none
% Copyright (C) 2001-2014 Dynare Team
% Copyright (C) 2001-2016 Dynare Team
%
% This file is part of Dynare.
%
@ -62,7 +62,7 @@ if ~isempty(i)
return;
end
f = 0.5*fvec'*fvec ;
f = 0.5*(fvec'*fvec) ;
if max(abs(fvec)) < tolf
return ;
@ -108,7 +108,7 @@ for its = 1:maxit
xold = x ;
fold = f ;
[x,f,fvec,check]=lnsrch1(xold,fold,g,p,stpmax,func,j1,j2,varargin{:});
[x, f, fvec, check] = lnsrch1(xold, fold, g, p, stpmax, func, j1, j2, tolx, varargin{:});
if debug
disp([its f])

6
matlab/solve_one_boundary.m
View File

@ -55,7 +55,7 @@ function [y, info] = solve_one_boundary(fname, y, x, params, steady_state, ...
% none.
%
% Copyright (C) 1996-2015 Dynare Team
% Copyright (C) 1996-2016 Dynare Team
%
% This file is part of Dynare.
%
@ -247,10 +247,10 @@ for it_=start:incr:finish
if is_dynamic
[ya,f,r,check]=lnsrch1(y(it_,:)',f,g,p,stpmax, ...
'lnsrch1_wrapper_one_boundary',nn, ...
y_index_eq, y_index_eq, fname, y, x, params, steady_state, it_);
y_index_eq, options.solve_tolx, y_index_eq, fname, y, x, params, steady_state, it_);
dx = ya' - y(it_, :);
else
[ya,f,r,check]=lnsrch1(y,f,g,p,stpmax,fname,nn,y_index_eq,x, ...
[ya,f,r,check]=lnsrch1(y,f,g,p,stpmax,fname,nn,y_index_eq, options.solve_tolx, x, ...
params, steady_state,0);
dx = ya - y(y_index_eq);
end

4
matlab/solve_two_boundaries.m
View File

@ -46,7 +46,7 @@ function [y, oo]= solve_two_boundaries(fname, y, x, params, steady_state, y_inde
% none.
%
% Copyright (C) 1996-2015 Dynare Team
% Copyright (C) 1996-2016 Dynare Team
%
% This file is part of Dynare.
%
@ -307,7 +307,7 @@ while ~(cvg==1 || iter>maxit_)
g = (ra'*g1a)';
f = 0.5*ra'*ra;
p = -g1a\ra;
[yn,f,ra,check]=lnsrch1(ya,f,g,p,stpmax,'lnsrch1_wrapper_two_boundaries',nn,nn, fname, y, y_index,x, params, steady_state, periods, y_kmin, Blck_size,options.periods);
[yn,f,ra,check]=lnsrch1(ya,f,g,p,stpmax,'lnsrch1_wrapper_two_boundaries',nn,nn, options.solve_tolx, fname, y, y_index,x, params, steady_state, periods, y_kmin, Blck_size,options.periods);
dx = ya - yn;
y(1+y_kmin:periods+y_kmin,y_index)=reshape(yn',length(y_index),periods)';
end