For deterministic simulations of a path between an arbitrary initial condition and the steady state, the path is initialized with the solution obtained from a first order approximation of the model => The number of iterations required for the convergence of the Newton is less important, and on a simple growth model with CES technology I observe also a reduction of the error size.

git-svn-id: https://www.dynare.org/svn/dynare/dynare_v4@1703 ac1d8469-bf42-47a9-8791-bf33cf982152
time-shift
adjemian 2008-02-11 14:42:46 +00:00
parent 806c29cc4a
commit 48f23191fb
3 changed files with 58 additions and 27 deletions

View File

@ -31,7 +31,17 @@ function make_y_
oo_.endo_simul = [ys0_*ones(1,M_.maximum_lag) oo_.steady_state*ones(1,options_.periods+M_.maximum_lead)];
end
elseif size(oo_.endo_simul,2) < M_.maximum_lag+M_.maximum_lead+options_.periods
oo_.endo_simul = [oo_.endo_simul oo_.steady_state*ones(1,M_.maximum_lag+options_.periods+M_.maximum_lead-size(oo_.endo_simul,2),1)];
end
%% A linear approximation is used to initiate the solution.
oldopt = options_;
options_.order = 1;
dr = oo_.dr;
dr.ys = oo_.steady_state;
[dr,info]=dr1(dr,0);
exogenous_variables = zeros(M_.maximum_lag+options_.periods+M_.maximum_lead-size(oo_.endo_simul,2)+1,0);
y0 = oo_.endo_simul(:,1:M_.maximum_lag);
oo_.endo_simul=simult_(y0,dr,exogenous_variables,1);
options_ = oldopt;
end

View File

@ -22,8 +22,11 @@
function y_=simult_(y0,dr,ex_,iorder)
global M_ options_ it_
iter = size(ex_,1);
nx = size(dr.ghu,2);
if ~isempty(dr.ghu)
nx = size(dr.ghu,2);
end
y_ = zeros(size(y0,1),iter+M_.maximum_lag);
y_(:,1:M_.maximum_lag) = y0;
k1 = [M_.maximum_lag:-1:1];
k2 = dr.kstate(find(dr.kstate(:,2) <= M_.maximum_lag+1),[1 2]);
@ -43,21 +46,38 @@ global M_ options_ it_
end
if iorder == 1
for i = M_.maximum_lag+1: iter+M_.maximum_lag
tempx1 = y_(dr.order_var,k1);
tempx2 = tempx1-repmat(dr.ys(dr.order_var),1,M_.maximum_lag);
tempx = tempx2(k2);
if options_.simul_algo == 0
y_(dr.order_var,i) = dr.ys(dr.order_var)+dr.ghx*tempx+dr.ghu* ...
ex_(i-M_.maximum_lag,:)';
elseif options_.simul_algo == 1
it_ = i;
m = dr.ys(dr.order_var);
[y_(:,i), check] = dynare_solve('ff_simul1',y_(:,i-1),tempx1(k3), ...
m(o3:end),tempx(k4),o1,o2,o3,k6);
if ~isempty(dr.ghu)
for i = M_.maximum_lag+1: iter+M_.maximum_lag
tempx1 = y_(dr.order_var,k1);
tempx2 = tempx1-repmat(dr.ys(dr.order_var),1,M_.maximum_lag);
tempx = tempx2(k2);
if options_.simul_algo == 0
y_(dr.order_var,i) = dr.ys(dr.order_var)+dr.ghx*tempx+dr.ghu* ...
ex_(i-M_.maximum_lag,:)';
elseif options_.simul_algo == 1
it_ = i;
m = dr.ys(dr.order_var);
[y_(:,i), check] = dynare_solve('ff_simul1',y_(:,i-1),tempx1(k3), ...
m(o3:end),tempx(k4),o1,o2,o3,k6);
end
k1 = k1+1;
end
else
for i = M_.maximum_lag+1: iter+M_.maximum_lag
tempx1 = y_(dr.order_var,k1);
tempx2 = tempx1-repmat(dr.ys(dr.order_var),1,M_.maximum_lag);
tempx = tempx2(k2);
if options_.simul_algo == 0
y_(dr.order_var,i) = dr.ys(dr.order_var)+dr.ghx*tempx;
elseif options_.simul_algo == 1
it_ = i;
m = dr.ys(dr.order_var);
[y_(:,i), check] = dynare_solve('ff_simul1',y_(:,i-1),tempx1(k3), ...
m(o3:end),tempx(k4),o1,o2,o3,k6);
end
k1 = k1+1;
end
end
k1 = k1+1;
end
elseif iorder == 2
for i = M_.maximum_lag+1: iter+M_.maximum_lag
tempx1 = y_(dr.order_var,k1);
@ -81,4 +101,4 @@ global M_ options_ it_
end
end
% MJ 08/30/02 corrected bug at order 2
% MJ 08/30/02 corrected bug at order 2

View File

@ -35,10 +35,11 @@ function steady()
disp(sprintf('%s \t\t %g',endo_names(i,:),steady_state(i)));
end
if isempty(ys0_)
oo_.endo_simul(:,1:M_.maximum_lag) = oo_.steady_state * ones(1,M_.maximum_lag);
else
options_ =set_default_option(options_,'periods',1);
oo_.endo_simul(:,M_.maximum_lag+1:M_.maximum_lag+options_.periods+M_.maximum_lead) = oo_.steady_state * ones(1,options_.periods+M_.maximum_lead);
end
%%% Unless I'm wrong, this is (should be?) done in make_y_.m
% $$$ if isempty(ys0_)
% $$$ oo_.endo_simul(:,1:M_.maximum_lag) = oo_.steady_state * ones(1, M_.maximum_lag);
% $$$ else
% $$$ options_ =set_default_option(options_,'periods',1);
% $$$ oo_.endo_simul(:,M_.maximum_lag+1:M_.maximum_lag+options_.periods+ ...
% $$$ M_.maximum_lead) = oo_.steady_state * ones(1,options_.periods+M_.maximum_lead);
% $$$ end