removed "global options_" in dynare_solve.m and made it an

argument. Modified all calling functions.
time-shift
Michel Juillard 2015-05-25 10:06:13 +02:00
parent d5cb78bbea
commit f03839cf12
7 changed files with 53 additions and 49 deletions

View File

@ -50,9 +50,11 @@ elseif options_.steadystate_flag
M.endo_names,'exact')]; M.endo_names,'exact')];
end end
if inst_nbr == 1 if inst_nbr == 1
inst_val = csolve(nl_func,ys_init(k_inst),'',options_.solve_tolf,100); %solve for instrument, using univariate solver, starting at initial value for instrument %solve for instrument, using univariate solver, starting at initial value for instrument
inst_val = csolve(nl_func,ys_init(k_inst),'',options_.solve_tolf,100);
else else
[inst_val,info1] = dynare_solve(nl_func,ys_init(k_inst),0); %solve for instrument, using multivariate solver, starting at initial value for instrument %solve for instrument, using multivariate solver, starting at initial value for instrument
[inst_val,info1] = dynare_solve(nl_func,ys_init(k_inst),options_);
end end
ys_init(k_inst) = inst_val; ys_init(k_inst) = inst_val;
exo_ss = [oo.exo_steady_state oo.exo_det_steady_state]; exo_ss = [oo.exo_steady_state oo.exo_det_steady_state];
@ -61,7 +63,7 @@ elseif options_.steadystate_flag
else else
n_var = M.orig_endo_nbr; n_var = M.orig_endo_nbr;
xx = oo.steady_state(1:n_var); xx = oo.steady_state(1:n_var);
[xx,check] = dynare_solve(nl_func,xx,0); [xx,check] = dynare_solve(nl_func,xx,options_);
[junk,junk,steady_state] = nl_func(xx); [junk,junk,steady_state] = nl_func(xx);
end end

View File

@ -1,12 +1,11 @@
function [x,info] = dynare_solve(func,x,jacobian_flag,varargin) function [x,info] = dynare_solve(func,x,options,varargin)
% function [x,info] = dynare_solve(func,x,jacobian_flag,varargin) % function [x,info] = dynare_solve(func,x,options,varargin)
% proposes different solvers % proposes different solvers
% %
% INPUTS % INPUTS
% func: name of the function to be solved % func: name of the function to be solved
% x: guess values % x: guess values
% jacobian_flag=1: jacobian given by the 'func' function % options: struct of Dynare options
% jacobian_flag=0: jacobian obtained numerically
% varargin: list of arguments following jacobian_flag % varargin: list of arguments following jacobian_flag
% %
% OUTPUTS % OUTPUTS
@ -16,7 +15,7 @@ function [x,info] = dynare_solve(func,x,jacobian_flag,varargin)
% SPECIAL REQUIREMENTS % SPECIAL REQUIREMENTS
% none % none
% Copyright (C) 2001-2014 Dynare Team % Copyright (C) 2001-2015 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -33,14 +32,16 @@ function [x,info] = dynare_solve(func,x,jacobian_flag,varargin)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>. % along with Dynare. If not, see <http://www.gnu.org/licenses/>.
global options_ % jacobian_flag=1: jacobian given by the 'func' function
% jacobian_flag=0: jacobian obtained numerically
jacobian_flag = options.jacobian_flag;
% Set tolerance parameter depending the the caller function. % Set tolerance parameter depending the the caller function.
stack = dbstack; stack = dbstack;
if strcmp(stack(2).file,'simulation_core.m') if strcmp(stack(2).file,'simulation_core.m')
tolf = options_.dynatol.f; tolf = options.dynatol.f;
else else
tolf = options_.solve_tolf; tolf = options.solve_tolf;
end end
info = 0; info = 0;
@ -79,7 +80,7 @@ if max(abs(fvec)) < tolf
return ; return ;
end end
if options_.solve_algo == 0 if options.solve_algo == 0
if ~isoctave if ~isoctave
if ~user_has_matlab_license('optimization_toolbox') if ~user_has_matlab_license('optimization_toolbox')
error('You can''t use solve_algo=0 since you don''t have MATLAB''s Optimization Toolbox') error('You can''t use solve_algo=0 since you don''t have MATLAB''s Optimization Toolbox')
@ -87,7 +88,7 @@ if options_.solve_algo == 0
end end
options=optimset('fsolve'); options=optimset('fsolve');
options.MaxFunEvals = 50000; options.MaxFunEvals = 50000;
options.MaxIter = options_.steady.maxit; options.MaxIter = options.steady.maxit;
options.TolFun = tolf; options.TolFun = tolf;
options.Display = 'iter'; options.Display = 'iter';
if jacobian_flag if jacobian_flag
@ -115,17 +116,17 @@ if options_.solve_algo == 0
else else
info = 1; info = 1;
end end
elseif options_.solve_algo == 1 elseif options.solve_algo == 1
[x,info]=solve1(func,x,1:nn,1:nn,jacobian_flag,options_.gstep, ... [x,info]=solve1(func,x,1:nn,1:nn,jacobian_flag,options.gstep, ...
tolf,options_.solve_tolx, ... tolf,options.solve_tolx, ...
options_.steady.maxit,options_.debug,varargin{:}); options.steady.maxit,options.debug,varargin{:});
elseif options_.solve_algo == 9 elseif options.solve_algo == 9
[x,info]=trust_region(func,x,1:nn,1:nn,jacobian_flag,options_.gstep, ... [x,info]=trust_region(func,x,1:nn,1:nn,jacobian_flag,options.gstep, ...
tolf,options_.solve_tolx, ... tolf,options.solve_tolx, ...
options_.steady.maxit,options_.debug,varargin{:}); options.steady.maxit,options.debug,varargin{:});
elseif options_.solve_algo == 2 || options_.solve_algo == 4 elseif options.solve_algo == 2 || options.solve_algo == 4
if options_.solve_algo == 2 if options.solve_algo == 2
solver = @solve1; solver = @solve1;
else else
solver = @trust_region; solver = @trust_region;
@ -133,7 +134,7 @@ elseif options_.solve_algo == 2 || options_.solve_algo == 4
if ~jacobian_flag if ~jacobian_flag
fjac = zeros(nn,nn) ; fjac = zeros(nn,nn) ;
dh = max(abs(x),options_.gstep(1)*ones(nn,1))*eps^(1/3); dh = max(abs(x),options.gstep(1)*ones(nn,1))*eps^(1/3);
for j = 1:nn for j = 1:nn
xdh = x ; xdh = x ;
xdh(j) = xdh(j)+dh(j) ; xdh(j) = xdh(j)+dh(j) ;
@ -143,18 +144,18 @@ elseif options_.solve_algo == 2 || options_.solve_algo == 4
[j1,j2,r,s] = dmperm(fjac); [j1,j2,r,s] = dmperm(fjac);
if options_.debug if options.debug
disp(['DYNARE_SOLVE (solve_algo=2|4): number of blocks = ' num2str(length(r))]); disp(['DYNARE_SOLVE (solve_algo=2|4): number of blocks = ' num2str(length(r))]);
end end
for i=length(r)-1:-1:1 for i=length(r)-1:-1:1
if options_.debug if options.debug
disp(['DYNARE_SOLVE (solve_algo=2|4): solving block ' num2str(i) ', of size ' num2str(r(i+1)-r(i)) ]); disp(['DYNARE_SOLVE (solve_algo=2|4): solving block ' num2str(i) ', of size ' num2str(r(i+1)-r(i)) ]);
end end
[x,info]=solver(func,x,j1(r(i):r(i+1)-1),j2(r(i):r(i+1)-1),jacobian_flag, ... [x,info]=solver(func,x,j1(r(i):r(i+1)-1),j2(r(i):r(i+1)-1),jacobian_flag, ...
options_.gstep, ... options.gstep, ...
tolf,options_.solve_tolx, ... tolf,options.solve_tolx, ...
options_.steady.maxit,options_.debug,varargin{:}); options.steady.maxit,options.debug,varargin{:});
if info if info
return return
end end
@ -162,17 +163,17 @@ elseif options_.solve_algo == 2 || options_.solve_algo == 4
fvec = feval(func,x,varargin{:}); fvec = feval(func,x,varargin{:});
if max(abs(fvec)) > tolf if max(abs(fvec)) > tolf
[x,info]=solver(func,x,1:nn,1:nn,jacobian_flag, ... [x,info]=solver(func,x,1:nn,1:nn,jacobian_flag, ...
options_.gstep, tolf,options_.solve_tolx, ... options.gstep, tolf,options.solve_tolx, ...
options_.steady.maxit,options_.debug,varargin{:}); options.steady.maxit,options.debug,varargin{:});
end end
elseif options_.solve_algo == 3 elseif options.solve_algo == 3
if jacobian_flag if jacobian_flag
[x,info] = csolve(func,x,func,1e-6,500,varargin{:}); [x,info] = csolve(func,x,func,1e-6,500,varargin{:});
else else
[x,info] = csolve(func,x,[],1e-6,500,varargin{:}); [x,info] = csolve(func,x,[],1e-6,500,varargin{:});
end end
elseif options_.solve_algo == 10 elseif options.solve_algo == 10
olmmcp = options_.lmmcp; olmmcp = options.lmmcp;
[x,fval,exitflag] = lmmcp(func,x,olmmcp.lb,olmmcp.ub,olmmcp,varargin{:}); [x,fval,exitflag] = lmmcp(func,x,olmmcp.lb,olmmcp.ub,olmmcp,varargin{:});
if exitflag == 1 if exitflag == 1
info = 0; info = 0;

View File

@ -1,5 +1,5 @@
function [x,info] = dynare_solve_block_or_bytecode(y, exo, params, options, M) function [x,info] = dynare_solve_block_or_bytecode(y, exo, params, options, M)
% Copyright (C) 2010-2012 Dynare Team % Copyright (C) 2010-2015 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -26,7 +26,7 @@ if options.block && ~options.bytecode
if options.solve_algo <= 4 if options.solve_algo <= 4
[y, check] = dynare_solve('block_mfs_steadystate', ... [y, check] = dynare_solve('block_mfs_steadystate', ...
ss(M.block_structure_stat.block(b).variable), ... ss(M.block_structure_stat.block(b).variable), ...
options.jacobian_flag, b, ss, exo, params, M); options, b, ss, exo, params, M);
if check ~= 0 if check ~= 0
% error(['STEADY: convergence % error(['STEADY: convergence
% problems in block ' int2str(b)]) % problems in block ' int2str(b)])
@ -62,18 +62,21 @@ elseif options.bytecode
if M.block_structure_stat.block(b).Simulation_Type ~= 1 && ... if M.block_structure_stat.block(b).Simulation_Type ~= 1 && ...
M.block_structure_stat.block(b).Simulation_Type ~= 2 M.block_structure_stat.block(b).Simulation_Type ~= 2
[y, check] = dynare_solve('block_bytecode_mfs_steadystate', ... [y, check] = dynare_solve('block_bytecode_mfs_steadystate', ...
x(M.block_structure_stat.block(b).variable), ... x(M.block_structure_stat ...
options.jacobian_flag, b, x, exo, params, M); .block(b).variable), ...
options, b, x, exo, params, M);
if check if check
% error(['STEADY: convergence % error(['STEADY: convergence problems in block '
% problems in block ' int2str(b)]) % int2str(b)])
info = 1; info = 1;
return return
end end
x(M.block_structure_stat.block(b).variable) = y; x(M.block_structure_stat.block(b).variable) = y;
else else
[chk, nulldev, nulldev1, x] = bytecode( x, exo, params, ... [chk, nulldev, nulldev1, x] = bytecode( x, exo, params, ...
x, 1, x, 'evaluate', 'static', ['block = ' int2str(b)]); x, 1, x, 'evaluate', ...
'static', ['block ' ...
'= '] int2str(b)]);
if chk if chk
info = 1; info = 1;
return return
@ -81,9 +84,8 @@ elseif options.bytecode
end; end;
end end
else else
[x, check] = dynare_solve('bytecode_steadystate', ... [x, check] = dynare_solve('bytecode_steadystate', y, ...
y, ... options, exo, params);
options.jacobian_flag, exo, params);
if check if check
% error('STEADY: convergence problems') % error('STEADY: convergence problems')
info = 1; info = 1;

View File

@ -146,7 +146,7 @@ y = repmat(steady_state,block_nbr,1);
old_options = options_; old_options = options_;
options_.solve_algo = options_.ep.solve_algo; options_.solve_algo = options_.ep.solve_algo;
options_.steady.maxit = options_.ep.maxit; options_.steady.maxit = options_.ep.maxit;
[y,info] = dynare_solve(@ep_problem_2,y,1,exo_simul,pfm); [y,info] = dynare_solve(@ep_problem_2,y,options_,exo_simul,pfm);
options_ = old_options; options_ = old_options;
if info if info
flag = 1; flag = 1;

View File

@ -117,8 +117,7 @@ function [ys,params,info] = evaluate_steady_state_file(ys_init,exo_ss,M,options)
end end
[ys,check] = dynare_solve('restricted_steadystate',... [ys,check] = dynare_solve('restricted_steadystate',...
ys(indv),... ys(indv),...
options.jacobian_flag, ... options, exo_ss,indv);
exo_ss,indv);
end end

View File

@ -89,7 +89,7 @@ else
% nonlinear models % nonlinear models
if max(abs(feval(fh,dr.ys,[oo_.exo_steady_state; ... if max(abs(feval(fh,dr.ys,[oo_.exo_steady_state; ...
oo_.exo_det_steady_state], M_.params))) > options_.dynatol.f oo_.exo_det_steady_state], M_.params))) > options_.dynatol.f
[dr.ys,check1] = dynare_solve(fh,dr.ys,options_.jacobian_flag,... [dr.ys,check1] = dynare_solve(fh,dr.ys,options_,...
[oo_.exo_steady_state; ... [oo_.exo_steady_state; ...
oo_.exo_det_steady_state], M_.params); oo_.exo_det_steady_state], M_.params);
end end

View File

@ -87,7 +87,7 @@ if options_.ramsey_policy && options_.ACES_solver == 0
end end
old_solve_algo = options_.solve_algo; old_solve_algo = options_.solve_algo;
% options_.solve_algo = 1; % options_.solve_algo = 1;
oo_.steady_state = dynare_solve('ramsey_static',oo_.steady_state,0,M_,options_,oo_,it_); oo_.steady_state = dynare_solve('ramsey_static',oo_.steady_state,options_,M_,options_,oo_,it_);
options_.solve_algo = old_solve_algo; options_.solve_algo = old_solve_algo;
[junk,junk,multbar] = ramsey_static(oo_.steady_state,M_,options_,oo_,it_); [junk,junk,multbar] = ramsey_static(oo_.steady_state,M_,options_,oo_,it_);
[jacobia_,M_] = ramsey_dynamic(oo_.steady_state,multbar,M_,options_,oo_,it_); [jacobia_,M_] = ramsey_dynamic(oo_.steady_state,multbar,M_,options_,oo_,it_);