diff --git a/matlab/ep/extended_path.m b/matlab/ep/extended_path.m
index 1bcd697e23bd283db176fff1f288cc648f158cd7..bec06dc1d73634f0690de39c270fccf25db7bf71 100644
--- a/matlab/ep/extended_path.m
+++ b/matlab/ep/extended_path.m
@@ -81,6 +81,10 @@ pfm.i_upd = pfm.ny+(1:pfm.periods*pfm.ny);
% Set the algorithm for the perfect foresight solver
options_.stack_solve_algo = options_.ep.stack_solve_algo;
+% Set check_stability flag
+do_not_check_stability_flag = ~options_.ep.check_stability;
+
+
% Compute the first order reduced form if needed.
%
% REMARK. It is assumed that the user did run the same mod file with stoch_simul(order=1) and save
@@ -248,84 +252,90 @@ while (t<sample_size)
end
end
end
- % 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) '.'])
+ if do_not_check_stability_flag
+ % Exit from the while loop.
+ oo_.endo_simul = tmp;
+ break
+ else
+ % 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
+ 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;
else
- disp(['Time: ' int2str(t) '. I increase the number of periods to ' int2str(options_.periods) '.'])
+ % 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
- 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
+ % 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
- 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...'])
+ 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
end
+ % Exit from the while loop.
+ break
end
- % Exit from the while loop.
- break
- end
- end% if info.convergence
+ end% if info.convergence
+ end
end% while
if ~info.convergence% If exited from the while loop without achieving convergence, use an homotopic approach
[INFO,tmp] = homotopic_steps(.5,.01);
diff --git a/matlab/global_initialization.m b/matlab/global_initialization.m
index e4cc50cc36ba490f6e25cb25e92e243bb4500d62..3ae136f2c6d02b3e3657a1ac0bef45863afd3607 100644
--- a/matlab/global_initialization.m
+++ b/matlab/global_initialization.m
@@ -126,6 +126,8 @@ ep.maxit = 500;
ep.periods = 200;
% Default step for increasing the number of periods if needed
ep.step = 50;
+% Set check_stability flag
+ep.check_stability = 1;
% Define last periods used to test if the solution is stable with respect to an increase in the number of periods.
ep.lp = 5;
% Define first periods used to test if the solution is stable with respect to an increase in the number of periods.