Skip to content
Snippets Groups Projects
Commit 43f46f28 authored by Stéphane Adjemian's avatar Stéphane Adjemian
Browse files

Removed stability test over the last periods of the perfect foresight solution.

parent 64ebd1d0
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
matlab/ep/extended_path.m
+ 73
− 83
View file @ 43f46f28
......@@ -248,94 +248,84 @@ while (t<sample_size)
end
end
end
% Test if periods is big enough. The variable delta measures the maximum absolute variation during
% the last periods of the simulated path. This variation has to be close to zero (because the
% economy is assumed to be in the steady state at the end of the simulated path).
if info.convergence
delta = max(max(abs(tmp(:,end-options_.ep.lp:end)-tmp(:,end-options_.ep.lp-1:end-1))));
end
if info.convergence && ~increase_periods && delta<options_.dynatol.x
% Exit from the while loop (the number of periods is big enough).
break
else
% Increase the number of periods.
options_.periods = options_.periods + options_.ep.step;
pfm.periods = options_.periods;
pfm.i_upd = pfm.ny+(1:pfm.periods*pfm.ny);
% Increment the counter.
increase_periods = increase_periods + 1;
if verbosity
if t<10
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
elseif t<100
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
elseif t<1000
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
else
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
end
end
if info.convergence
% If the previous call to the perfect foresight model solver exited
% announcing that the routine converged, adapt the size of oo_.endo_simul
% and oo_.exo_simul.
oo_.endo_simul = [ tmp , repmat(oo_.steady_state,1,options_.ep.step) ];
oo_.exo_simul = [ oo_.exo_simul ; zeros(options_.ep.step,size(shocks,2)) ];
tmp_old = tmp;
% Test if periods is big enough.
% Increase the number of periods.
options_.periods = options_.periods + options_.ep.step;
pfm.periods = options_.periods;
pfm.i_upd = pfm.ny+(1:pfm.periods*pfm.ny);
% Increment the counter.
increase_periods = increase_periods + 1;
if verbosity
if t<10
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
elseif t<100
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
elseif t<1000
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
else
% If the previous call to the perfect foresight model solver exited
% announcing that the routine did not converge, then tmp=1... Maybe
% should change that, because in some circonstances it may usefull
% to know where the routine did stop, even if convergence was not
% achieved.
oo_.endo_simul = [ oo_.endo_simul , repmat(oo_.steady_state,1,options_.ep.step) ];
oo_.exo_simul = [ oo_.exo_simul ; zeros(options_.ep.step,size(shocks,2)) ];
disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
end
% Solve the perfect foresight model with an increased number of periods.
[flag,tmp] = bytecode('dynamic');
if flag
[flag,tmp] = solve_perfect_foresight_model(oo_.endo_simul,oo_.exo_simul,pfm);
end
if info.convergence
% If the previous call to the perfect foresight model solver exited
% announcing that the routine converged, adapt the size of oo_.endo_simul
% and oo_.exo_simul.
oo_.endo_simul = [ tmp , repmat(oo_.steady_state,1,options_.ep.step) ];
oo_.exo_simul = [ oo_.exo_simul ; zeros(options_.ep.step,size(shocks,2)) ];
tmp_old = tmp;
else
% If the previous call to the perfect foresight model solver exited
% announcing that the routine did not converge, then tmp=1... Maybe
% should change that, because in some circonstances it may usefull
% to know where the routine did stop, even if convergence was not
% achieved.
oo_.endo_simul = [ oo_.endo_simul , repmat(oo_.steady_state,1,options_.ep.step) ];
oo_.exo_simul = [ oo_.exo_simul ; zeros(options_.ep.step,size(shocks,2)) ];
end
% Solve the perfect foresight model with an increased number of periods.
[flag,tmp] = bytecode('dynamic');
if flag
[flag,tmp] = solve_perfect_foresight_model(oo_.endo_simul,oo_.exo_simul,pfm);
end
info.convergence = ~flag;
if info.convergence
% If the solver achieved convergence, check that simulated paths did not
% change during the first periods.
% Compute the maximum deviation between old path and new path over the
% first periods
delta = max(max(abs(tmp(:,2:options_.ep.fp)-tmp_old(:,2:options_.ep.fp))));
if delta<options_.dynatol.x
% If the maximum deviation is close enough to zero, reset the number
% of periods to ep.periods
options_.periods = options_.ep.periods;
pfm.periods = options_.periods;
pfm.i_upd = pfm.ny+(1:pfm.periods*pfm.ny);
% Cut oo_.exo_simul and oo_.endo_simul consistently with the resetted
% number of periods and exit from the while loop.
oo_.exo_simul = oo_.exo_simul(1:(options_.periods+2),:);
oo_.endo_simul = oo_.endo_simul(:,1:(options_.periods+2));
break
end
info.convergence = ~flag;
if info.convergence
% If the solver achieved convergence, check that simulated paths did not
% change during the first periods.
% Compute the maximum deviation between old path and new path over the
% first periods
delta = max(max(abs(tmp(:,2:options_.ep.fp)-tmp_old(:,2:options_.ep.fp))));
if delta<options_.dynatol.x
% If the maximum deviation is close enough to zero, reset the number
% of periods to ep.periods
options_.periods = options_.ep.periods;
pfm.periods = options_.periods;
pfm.i_upd = pfm.ny+(1:pfm.periods*pfm.ny);
% Cut oo_.exo_simul and oo_.endo_simul consistently with the resetted
% number of periods and exit from the while loop.
oo_.exo_simul = oo_.exo_simul(1:(options_.periods+2),:);
oo_.endo_simul = oo_.endo_simul(:,1:(options_.periods+2));
break
end
else
% The solver did not converge... Try to solve the model again with a bigger
% number of periods, except if the number of periods has been increased more
% than 10 times.
if increase_periods==10;
if verbosity
if t<10
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
elseif t<100
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
elseif t<1000
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
else
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
end
else
% The solver did not converge... Try to solve the model again with a bigger
% number of periods, except if the number of periods has been increased more
% than 10 times.
if increase_periods==10;
if verbosity
if t<10
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
elseif t<100
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
elseif t<1000
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
else
disp(['Time: ' int2str(t) '. Even with ' int2str(options_.periods) ', I am not able to solve the perfect foresight model. Use homotopy instead...'])
end
% Exit from the while loop.
break
end
end% if info.convergence
end% ~increase_periods && delta<options_.dynatol.x
% Exit from the while loop.
break
end
end% if info.convergence
end% while
if ~info.convergence% If exited from the while loop without achieving convergence, use an homotopic approach
[INFO,tmp] = homotopic_steps(.5,.01);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment