diff --git a/matlab/+occbin/DSGE_smoother.m b/matlab/+occbin/DSGE_smoother.m
index 03d1c560baf7a2754f433eb90e499a98f0a4240d..281fb7224bd6ad7b56a0ade2d98a824df5975a60 100644
--- a/matlab/+occbin/DSGE_smoother.m
+++ b/matlab/+occbin/DSGE_smoother.m
@@ -106,6 +106,7 @@ opts_simul.maxit = options_.occbin.smoother.maxit;
opts_simul.waitbar = options_.occbin.smoother.waitbar;
opts_simul.periods = options_.occbin.smoother.periods;
opts_simul.check_ahead_periods = options_.occbin.smoother.check_ahead_periods;
+opts_simul.max_check_ahead_periods = options_.occbin.smoother.max_check_ahead_periods;
opts_simul.periodic_solution = options_.occbin.smoother.periodic_solution;
opts_simul.full_output = options_.occbin.smoother.full_output;
opts_simul.piecewise_only = options_.occbin.smoother.piecewise_only;
diff --git a/matlab/+occbin/IVF_core.m b/matlab/+occbin/IVF_core.m
index b9e96622bc8e8464a8e4c46d26282d0647579f6b..efc3f6e844012c95f724fec1f01043a006f2e5ae 100644
--- a/matlab/+occbin/IVF_core.m
+++ b/matlab/+occbin/IVF_core.m
@@ -57,6 +57,7 @@ opts_simul.maxit = options_.occbin.likelihood.maxit;
opts_simul.waitbar = false;
opts_simul.periods = options_.occbin.likelihood.periods;
opts_simul.check_ahead_periods = options_.occbin.likelihood.check_ahead_periods;
+opts_simul.max_check_ahead_periods = options_.occbin.likelihood.max_check_ahead_periods;
opts_simul.periodic_solution = options_.occbin.likelihood.periodic_solution;
opts_simul.restrict_state_space = options_.occbin.likelihood.restrict_state_space;
opts_simul.piecewise_only = 1;
diff --git a/matlab/+occbin/set_default_options.m b/matlab/+occbin/set_default_options.m
index 3d827be16ead618da19230b668be42134bbb4e62..17e74a323448935fcdb9a57d6bb23a3f2c59798d 100644
--- a/matlab/+occbin/set_default_options.m
+++ b/matlab/+occbin/set_default_options.m
@@ -76,6 +76,7 @@ if ismember(flag,{'likelihood','all'})
options_occbin_.likelihood.IVF_shock_observable_mapping = [];
options_occbin_.likelihood.maxit = 30; % this is for occbin solver algo
options_occbin_.likelihood.max_number_of_iterations = 10; % this is for occbin_kalman_update loop
+ options_occbin_.likelihood.max_check_ahead_periods=inf;
options_occbin_.likelihood.periods = 100;
options_occbin_.likelihood.check_ahead_periods=200;
options_occbin_.likelihood.periodic_solution=false;
@@ -170,11 +171,12 @@ if ismember(flag,{'simul','all'})
options_occbin_.simul.init_binding_indicator=false(0);
options_occbin_.simul.exo_pos=1:M_.exo_nbr;
options_occbin_.simul.maxit=30;
- options_occbin_.simul.max_periods=inf;
+ options_occbin_.simul.max_check_ahead_periods=inf;
options_occbin_.simul.periods=100;
options_occbin_.simul.check_ahead_periods=200;
options_occbin_.simul.periodic_solution=false;
options_occbin_.simul.piecewise_only = false;
+ options_occbin_.simul.reset_check_ahead_periods_in_new_period = false;
options_occbin_.simul.reset_regime_in_new_period = false;
options_occbin_.simul.restrict_state_space=false;
options_occbin_.simul.SHOCKS=zeros(1,M_.exo_nbr);
@@ -193,6 +195,7 @@ if ismember(flag,{'smoother','all'})
options_occbin_.smoother.inversion_filter = false;
options_occbin_.smoother.linear_smoother = true;
options_occbin_.smoother.maxit = 30; % this is for occbin solver algo
+ options_occbin_.smoother.max_check_ahead_periods=inf;
options_occbin_.smoother.max_number_of_iterations = 10; % this is for smoother loop
options_occbin_.smoother.periods = 100;
options_occbin_.smoother.check_ahead_periods=200;
diff --git a/matlab/+occbin/setup.m b/matlab/+occbin/setup.m
index 7c4749623f1776958e100a8db3d75cb0efda37b3..a050ee787b0b62fc03a48608cc989b2ca7766619 100644
--- a/matlab/+occbin/setup.m
+++ b/matlab/+occbin/setup.m
@@ -2,12 +2,12 @@ function [M_, options_] = setup(M_,options_, options_occbin_)
% function [M_, options_] = setup(M_, options_, options_occbin_)
% Sets up run of Occbin: creates shock matrix, sets options
%
-% INPUT:
+% INPUT:
% - M_ [structure] Matlab's structure describing the model
% - options_ [structure] Matlab's structure containing the options
% - options_occbin_ [structure] Matlab's structure containing Occbin options
%
-% OUTPUT:
+% OUTPUT:
% - M_ [structure] Matlab's structure describing the model
% - options_occbin_ [structure] Matlab's structure containing Occbin options
@@ -55,7 +55,7 @@ if isfield(M_,'surprise_shocks') && ~isempty(M_.surprise_shocks)
for ii = 1:length(M_.surprise_shocks)
ivar = M_.surprise_shocks(ii).exo_id;
temp(M_.surprise_shocks(ii).periods,ivar) = M_.surprise_shocks(ii).value;
- end
+ end
shock_index=~all(temp==0,1);
options_.occbin.simul.SHOCKS=temp(:,shock_index);
options_.occbin.simul.exo_pos=find(shock_index);
diff --git a/matlab/+occbin/solve_one_constraint.m b/matlab/+occbin/solve_one_constraint.m
index 0801f8c1ebfab61a3b8e4658b30fee5bd22454a2..f5f0c1446c79ff6edf9b1c829e8b98c90f2216ca 100644
--- a/matlab/+occbin/solve_one_constraint.m
+++ b/matlab/+occbin/solve_one_constraint.m
@@ -125,21 +125,29 @@ for shock_period = 1:n_shocks_periods
regime_change_this_iteration=true;
iter = 0;
+ nperiods_endogenously_increased = false;
guess_history_it = false;
if guess_history && (shock_period<=length(regime_history_guess)) %beyond guess regime history
guess_history_it = true;
end
is_periodic=false;
+ is_periodic_loop=false;
binding_indicator_history={};
max_err = NaN(max_iter,1);
+ regime_violates_constraint_in_expectation = false(max_iter,1);
- while (regime_change_this_iteration && iter<max_iter && ~is_periodic)
+ while (regime_change_this_iteration && iter<max_iter && ~is_periodic && ~is_periodic_loop)
iter = iter +1;
- if binding_indicator(end) && nperiods_0<opts_simul_.max_periods
+ if binding_indicator(end) && nperiods_0<opts_simul_.max_check_ahead_periods
binding_indicator = [binding_indicator; false ];
nperiods_0 = nperiods_0 + 1;
+ nperiods_endogenously_increased = true;
disp_verbose(['nperiods has been endogenously increased up to ' int2str(nperiods_0) '.'],opts_simul_.debug)
+ elseif nperiods_0>=opts_simul_.max_check_ahead_periods
+ % enforce endogenously increased nperiods to not violate max_check_ahead_periods
+ binding_indicator = binding_indicator(1:opts_simul_.max_check_ahead_periods+1);
+ binding_indicator(end)=false;
end
if length(binding_indicator)<(nperiods_0 + 1)
binding_indicator=[binding_indicator; false(nperiods_0 + 1-length(binding_indicator),1)];
@@ -164,7 +172,13 @@ for shock_period = 1:n_shocks_periods
% get the hypothesized piece wise linear solution
if shock_period==1 || shock_period>1 && any(data.shocks_sequence(shock_period,:))
if iter==1 && opts_simul_.reset_regime_in_new_period
- binding_indicator=false(size(binding_indicator));
+ if opts_simul_.reset_check_ahead_periods_in_new_period
+ % I re-set check ahead periods to initial value, when in previous period it was endogenously increased
+ nperiods_0 = max(opts_simul_.check_ahead_periods,n_periods-n_shocks_periods);
+ binding_indicator = false(nperiods_0+1,1);
+ else
+ binding_indicator=false(size(binding_indicator));
+ end
binding_indicator_history{iter}=binding_indicator;
% analyse violvec and isolate contiguous periods in the other regime.
[regime, regime_start, error_code_period]=occbin.map_regime(binding_indicator,opts_simul_.debug);
@@ -178,42 +192,81 @@ for shock_period = 1:n_shocks_periods
[binding, relax, err]=feval([M_.fname,'.occbin_difference'],zdatalinear_+repmat(dr_base.ys',size(zdatalinear_,1),1),M_.params,dr_base.ys);
+ if ~isinf(opts_simul_.max_check_ahead_periods) && opts_simul_.max_check_ahead_periods<length(binding_indicator)
+ end_periods = opts_simul_.max_check_ahead_periods;
+ last_indicator = false(length(binding_indicator)-end_periods,1);
+ else
+ end_periods = length(binding_indicator);
+ last_indicator = false(0);
+ end
% check if changes to the hypothesis of the duration for each
% regime
- if any(binding.constraint_1 & ~binding_indicator) || any(relax.constraint_1 & binding_indicator)
- err_viol = err.binding_constraint_1(binding.constraint_1 & ~binding_indicator);
- err_relax = err.relax_constraint_1(relax.constraint_1 & binding_indicator);
+ if any(binding.constraint_1(1:end_periods) & ~binding_indicator(1:end_periods)) || any(relax.constraint_1(1:end_periods) & binding_indicator(1:end_periods))
+
+ err_viol = err.binding_constraint_1(binding.constraint_1(1:end_periods) & ~binding_indicator(1:end_periods));
+ err_relax = err.relax_constraint_1(relax.constraint_1(1:end_periods) & binding_indicator(1:end_periods));
max_err(iter) = max(abs([err_viol;err_relax]));
regime_change_this_iteration = true;
+ regime_violates_constraint_in_expectation(iter) = any(binding.constraint_1(1:end_periods) & ~binding_indicator(1:end_periods));
else
regime_change_this_iteration = false;
max_err(iter) = 0;
end
+ % if max_check_ahead_periods<inf, enforce unconstrained regime for period larger than max_check_ahead_periods
+ binding_constraint_new=[binding.constraint_1(1:end_periods);last_indicator];
+ relaxed_constraint_new = [relax.constraint_1(1:end_periods);not(last_indicator)];
if curb_retrench % apply Gauss-Seidel idea of slowing down the change in the guess
% for the constraint -- only relax one
% period at a time starting from the last
% one when each of the constraints is true.
retrench = false(numel(binding_indicator),1);
- max_relax_constraint_1=find(relax.constraint_1 & binding_indicator,1,'last');
- if ~isempty(max_relax_constraint_1) && find(relax.constraint_1,1,'last')>=find(binding_indicator,1,'last')
+ max_relax_constraint_1=find(relax.constraint_1(1:end_periods) & binding_indicator(1:end_periods),1,'last');
+ if ~isempty(max_relax_constraint_1) && find(relax.constraint_1(1:end_periods),1,'last')>=find(binding_indicator(1:end_periods),1,'last')
retrench(max_relax_constraint_1) = true;
end
- binding_indicator = (binding_indicator | binding.constraint_1) & ~ retrench;
+ binding_indicator = (binding_indicator | binding_constraint_new) & ~ retrench;
else
- binding_indicator= (binding_indicator | binding.constraint_1) & ~(binding_indicator & relax.constraint_1);
+ binding_indicator = (binding_indicator | binding_constraint_new) & ~(binding_indicator & relaxed_constraint_new);
end
- if iter>1 && regime_change_this_iteration
+ if iter>1 && regime_change_this_iteration && ~nperiods_endogenously_increased
+ % check for periodic solution only if nperiods is not
+ % increased endogenously
+ % first check for infinite loop
+ is_periodic_loop = false(iter-1,1);
for kiter=1:iter-1
- vvv = [binding_indicator_history{kiter}; false(size(binding_indicator,1)- size(binding_indicator_history{kiter},1), 1)];
- is_periodic(kiter) = isequal(vvv, binding_indicator);
+ if size(binding_indicator,1)== size(binding_indicator_history{kiter},1)
+ % vvv = [binding_indicator_history{kiter}; false(size(binding_indicator,1)- size(binding_indicator_history{kiter},1), 1)];
+ % is_periodic(kiter) = isequal(vvv, binding_indicator);
+ is_periodic_loop(kiter) = isequal(binding_indicator_history{kiter}, binding_indicator);
+ else
+ is_periodic_loop(kiter) = false;
+ end
+ end
+ is_periodic_loop_all =is_periodic_loop;
+ is_periodic_loop = any(is_periodic_loop);
+ % only accept periodicity where regimes differ by one
+ % period!
+ is_periodic = false(iter-1,1);
+ for kiter=iter-1
+ if size(binding_indicator,1)== size(binding_indicator_history{kiter},1)
+ % vvv = [binding_indicator_history{kiter}; false(size(binding_indicator,1)- size(binding_indicator_history{kiter},1), 1)];
+ % is_periodic(kiter) = isequal(vvv, binding_indicator);
+ is_periodic(kiter) = isequal(binding_indicator_history{kiter}, binding_indicator) && length(find(binding_indicator_history{iter}-binding_indicator))==1;
+ else
+ is_periodic(kiter)=false;
+ end
end
is_periodic_all =is_periodic;
is_periodic = any(is_periodic);
if is_periodic && periodic_solution
- [merr,imerr]=min(max_err(find(is_periodic_all,1):end));
inx = find(is_periodic_all,1):iter;
+ inx1 = inx(find(~regime_violates_constraint_in_expectation(inx)));
+ if not(isempty(inx1))
+ inx=inx1;
+ end
+ [merr,imerr]=min(max_err(inx));
inx = inx(imerr);
binding_indicator=binding_indicator_history{inx};
if inx<iter
@@ -265,9 +318,14 @@ for shock_period = 1:n_shocks_periods
return
end
else
+ if is_periodic_loop
+ disp_verbose('Did not converge -- infinite loop of guess regimes.',opts_simul_.debug)
+ error_flag = 313;
+ else
disp_verbose('Did not converge -- increase maxit.',opts_simul_.debug)
+ error_flag = 311;
+ end
if opts_simul_.waitbar; dyn_waitbar_close(hh); end
- error_flag = 311;
return
end
else
diff --git a/matlab/+occbin/solve_two_constraints.m b/matlab/+occbin/solve_two_constraints.m
index 234aaa78b510dfc93cc14d6b6c4b1532fc5c621c..55fff162d7501edf1f93ba9a6b14cb180229febf 100644
--- a/matlab/+occbin/solve_two_constraints.m
+++ b/matlab/+occbin/solve_two_constraints.m
@@ -133,23 +133,35 @@ for shock_period = 1:n_shocks_periods
dyn_waitbar(shock_period/n_shocks_periods, hh, sprintf('Period %u of %u', shock_period,n_shocks_periods));
end
regime_change_this_iteration=true;
+ nperiods_endogenously_increased = false;
iter = 0;
guess_history_it = false;
if guess_history && (shock_period<=length(regime_history_guess)) %beyond guess regime history
guess_history_it = true;
end
is_periodic=false;
+ is_periodic_loop=false;
binding_indicator_history={};
max_err = NaN(max_iter,1);
+ regime_violates_constraint_in_expectation = false(max_iter,1);
- while (regime_change_this_iteration && iter<max_iter && ~is_periodic)
+ while (regime_change_this_iteration && iter<max_iter && ~is_periodic && ~is_periodic_loop)
iter = iter +1;
- if any(binding_indicator(end,:)) && nperiods_0<opts_simul_.max_periods
+ if any(binding_indicator(end,:)) && nperiods_0<opts_simul_.max_check_ahead_periods
binding_indicator = [binding_indicator; false(1,2)];
nperiods_0 = nperiods_0 + 1;
+ nperiods_endogenously_increased = true;
disp_verbose(['nperiods has been endogenously increased up to ' int2str(nperiods_0) '.'],opts_simul_.debug)
+ elseif nperiods_0>=opts_simul_.max_check_ahead_periods
+ % enforce endogenously increased nperiods to not violate max_check_ahead_periods
+ binding_indicator = binding_indicator(1:opts_simul_.max_check_ahead_periods+1,:);
+ binding_indicator(end,:)=false(1,2);
end
if size(binding_indicator,1)<(nperiods_0 + 1)
+ % to ensure the simulation is run for the required nperiods
+ % even beyond max_check_ahead_periods: the latter controls check ahead periods
+ % and NOT how many periods we simulate after we are back to
+ % unconstrained regime (nperiods_0)
binding_indicator=[binding_indicator; false(nperiods_0 + 1-size(binding_indicator,1),2)];
end
binding_indicator_history{iter}=binding_indicator;
@@ -178,7 +190,13 @@ for shock_period = 1:n_shocks_periods
regime_history(shock_period).regimestart2 = regime_start_2;
if shock_period==1 || shock_period>1 && any(data.shocks_sequence(shock_period,:)) % first period or shock happening
if iter==1 && opts_simul_.reset_regime_in_new_period
- binding_indicator=false(size(binding_indicator));
+ if opts_simul_.reset_check_ahead_periods_in_new_period
+ % I re-set check ahead periods to initial value, when in previous period it was endogenously increased
+ nperiods_0 = max(opts_simul_.check_ahead_periods,n_periods-n_shocks_periods);
+ binding_indicator = false(nperiods_0+1,2);
+ else
+ binding_indicator=false(size(binding_indicator));
+ end
binding_indicator_history{iter}=binding_indicator;
% analyse violvec and isolate contiguous periods in the other regime.
[regime_1, regime_start_1, error_code_period(1)]=occbin.map_regime(binding_indicator(:,1),opts_simul_.debug);
@@ -195,16 +213,25 @@ for shock_period = 1:n_shocks_periods
data.exo_pos,data.shocks_sequence(shock_period,:),endo_init, update_flag);
[binding, relax, err]=feval([M_.fname,'.occbin_difference'],zdatalinear_+repmat(dr.ys',size(zdatalinear_,1),1),M_.params,dr.ys);
- binding_constraint_new=[binding.constraint_1;binding.constraint_2];
- relaxed_constraint_new = [relax.constraint_1;relax.constraint_2];
+
+ if ~isinf(opts_simul_.max_check_ahead_periods) && opts_simul_.max_check_ahead_periods<length(binding_indicator)
+ end_periods = opts_simul_.max_check_ahead_periods;
+ last_indicator = false(length(binding_indicator)-end_periods,1);
+ else
+ end_periods = length(binding_indicator);
+ last_indicator = false(0);
+ end
+ binding_constraint_new=[binding.constraint_1(1:end_periods); binding.constraint_2(1:end_periods)];
+ relaxed_constraint_new = [relax.constraint_1(1:end_periods); relax.constraint_2(1:end_periods)];
+ my_binding_indicator = binding_indicator(1:end_periods,:);
- err_binding_constraint_new = [err.binding_constraint_1; err.binding_constraint_2];
- err_relaxed_constraint_new = [err.relax_constraint_1; err.relax_constraint_2];
+ err_binding_constraint_new = [err.binding_constraint_1(1:end_periods); err.binding_constraint_2(1:end_periods)];
+ err_relaxed_constraint_new = [err.relax_constraint_1(1:end_periods); err.relax_constraint_2(1:end_periods)];
% check if changes_
- if any(binding_constraint_new & ~binding_indicator(:)) || any(relaxed_constraint_new & binding_indicator(:))
- err_violation = err_binding_constraint_new(binding_constraint_new & ~binding_indicator(:));
- err_relax = err_relaxed_constraint_new(relaxed_constraint_new & binding_indicator(:));
+ if any(binding_constraint_new & ~my_binding_indicator(:)) || any(relaxed_constraint_new & my_binding_indicator(:))
+ err_violation = err_binding_constraint_new(binding_constraint_new & ~my_binding_indicator(:));
+ err_relax = err_relaxed_constraint_new(relaxed_constraint_new & my_binding_indicator(:));
max_err(iter) = max(abs([err_violation;err_relax]));
regime_change_this_iteration = true;
else
@@ -212,36 +239,64 @@ for shock_period = 1:n_shocks_periods
max_err(iter) = 0;
end
+ binding_constraint_new=[binding.constraint_1(1:end_periods);last_indicator; binding.constraint_2(1:end_periods);last_indicator];
+ relaxed_constraint_new = [relax.constraint_1(1:end_periods);not(last_indicator); relax.constraint_2(1:end_periods);not(last_indicator)];
if curb_retrench % apply Gauss-Seidel idea of slowing down the change in the guess
% for the constraint -- only relax one
% period at a time starting from the last
% one when each of the constraints is true.
retrench = false(numel(binding_indicator),1);
- max_relax_constraint_1=find(relax.constraint_1 & binding_indicator(:,1),1,'last');
- if ~isempty(max_relax_constraint_1) && find(relax.constraint_1,1,'last')>=find(binding_indicator(:,1),1,'last')
+ max_relax_constraint_1=find(relax.constraint_1(1:end_periods) & binding_indicator(1:end_periods,1),1,'last');
+ if ~isempty(max_relax_constraint_1) && find(relax.constraint_1(1:end_periods),1,'last')>=find(binding_indicator(1:end_periods,1),1,'last')
retrench(max_relax_constraint_1) = true;
end
- max_relax_constraint_2=find(relax.constraint_2 & binding_indicator(:,2),1,'last');
- if ~isempty(max_relax_constraint_2) && find(relax.constraint_2,1,'last')>=find(binding_indicator(:,2),1,'last')
+ max_relax_constraint_2=find(relax.constraint_2(1:end_periods) & binding_indicator(1:end_periods,2),1,'last');
+ if ~isempty(max_relax_constraint_2) && find(relax.constraint_2(1:end_periods),1,'last')>=find(binding_indicator(1:end_periods,2),1,'last')
retrench(max_relax_constraint_2+nperiods_0+1) = true;
end
binding_indicator = (binding_indicator(:) | binding_constraint_new) & ~ retrench;
else
- binding_indicator= (binding_indicator(:) | binding_constraint_new) & ~(binding_indicator(:) & relaxed_constraint_new);
+ binding_indicator = (binding_indicator(:) | binding_constraint_new) & ~(binding_indicator(:) & relaxed_constraint_new);
end
binding_indicator = reshape(binding_indicator,nperiods_0+1,2);
- if iter>1 && regime_change_this_iteration
- is_periodic=false(1,iter-1);
+ if iter>1 && regime_change_this_iteration && ~nperiods_endogenously_increased
+ % check for periodic solution only if nperiods is not
+ % increased endogenously
+ % first check for infinite loop
+ is_periodic_loop = false(iter-1,1);
for kiter=1:iter-1
- vvv = [binding_indicator_history{kiter}; false(size(binding_indicator,1)- size(binding_indicator_history{kiter},1), 2)];
- is_periodic(kiter) = isequal(vvv, binding_indicator);
+ if size(binding_indicator,1)== size(binding_indicator_history{kiter},1)
+ % vvv = [binding_indicator_history{kiter}; false(size(binding_indicator,1)- size(binding_indicator_history{kiter},1), 1)];
+ % is_periodic(kiter) = isequal(vvv, binding_indicator);
+ is_periodic_loop(kiter) = isequal(binding_indicator_history{kiter}, binding_indicator);
+ else
+ is_periodic_loop(kiter) = false;
+ end
+ end
+% is_periodic_loop_all =is_periodic_loop;
+ is_periodic_loop = any(is_periodic_loop);
+ % only accept periodicity where regimes differ by one
+ % period!
+ is_periodic=false(1,iter-1);
+ for kiter=iter-1
+ if size(binding_indicator,1)== size(binding_indicator_history{kiter},1)
+ % vvv = [binding_indicator_history{kiter}; false(size(binding_indicator,1)- size(binding_indicator_history{kiter},1), 1)];
+ % is_periodic(kiter) = isequal(vvv, binding_indicator);
+ is_periodic(kiter) = isequal(binding_indicator_history{kiter}, binding_indicator) && length(find(binding_indicator_history{iter}(:,1)-binding_indicator(:,1)))<=1 && length(find(binding_indicator_history{iter}(:,2)-binding_indicator(:,2)))<=1;
+ else
+ is_periodic(kiter)=false;
+ end
end
is_periodic_all = is_periodic;
is_periodic = any(is_periodic);
if is_periodic && periodic_solution
- [min_err,index_min_err]=min(max_err(find(is_periodic_all,1):end));
inx = find(is_periodic_all,1):iter;
+ inx1 = inx(find(~regime_violates_constraint_in_expectation(inx)));
+ if not(isempty(inx1))
+ inx=inx1;
+ end
+ [min_err,index_min_err]=min(max_err(inx));
inx = inx(index_min_err);
binding_indicator=binding_indicator_history{inx}; %select regime history with same result, but smallest error
if inx<iter %update if needed
@@ -295,8 +350,13 @@ for shock_period = 1:n_shocks_periods
return;
end
else
- disp_verbose('Did not converge -- increase maxit.',opts_simul_.debug)
- error_flag = 311;
+ if is_periodic_loop
+ disp_verbose('Did not converge -- infinite loop of guess regimes.',opts_simul_.debug)
+ error_flag = 313;
+ else
+ disp_verbose('Did not converge -- increase maxit.',opts_simul_.debug)
+ error_flag = 311;
+ end
if opts_simul_.waitbar; dyn_waitbar_close(hh); end
return;
end
diff --git a/matlab/+occbin/solver.m b/matlab/+occbin/solver.m
index 2d9be140532f7673d3b3c3682f71fb868ad4cd9a..9dd44f88fa68a76464de38df9d1c838a29f1166c 100644
--- a/matlab/+occbin/solver.m
+++ b/matlab/+occbin/solver.m
@@ -68,6 +68,12 @@ else
oo_.dr=sto_dr;
end
+if options_.occbin.simul.check_ahead_periods>options_.occbin.simul.max_check_ahead_periods
+ options_.occbin.simul.check_ahead_periods=options_.occbin.simul.max_check_ahead_periods;
+ disp(['occbin::options::' simul '_check_ahead_periods cannot exceed ' simul '_max_check_ahead_periods'])
+ disp(['occbin::options::' simul '_check_ahead_periods is re-set to be equal to ' simul '_max_check_ahead_periods'])
+end
+
if M_.occbin.constraint_nbr==1
[out, ss, error_flag ] = occbin.solve_one_constraint(M_,oo_.dr,options_.occbin.simul,solve_dr);
elseif M_.occbin.constraint_nbr==2
diff --git a/matlab/dsge_likelihood.m b/matlab/dsge_likelihood.m
index 72395f1667b5e8d5874c64baaff951c01809ee7c..e6fcd04f9d11e8868797be7b39b4f934b92cd62b 100644
--- a/matlab/dsge_likelihood.m
+++ b/matlab/dsge_likelihood.m
@@ -976,6 +976,7 @@ occbin_options.opts_simul.curb_retrench = DynareOptions.occbin.likelihood.curb_r
occbin_options.opts_simul.maxit = DynareOptions.occbin.likelihood.maxit;
occbin_options.opts_simul.periods = DynareOptions.occbin.likelihood.periods;
occbin_options.opts_simul.check_ahead_periods = DynareOptions.occbin.likelihood.check_ahead_periods;
+occbin_options.opts_simul.max_check_ahead_periods = DynareOptions.occbin.likelihood.max_check_ahead_periods;
occbin_options.opts_simul.periodic_solution = DynareOptions.occbin.likelihood.periodic_solution;
occbin_options.opts_simul.restrict_state_space = DynareOptions.occbin.likelihood.restrict_state_space;
diff --git a/matlab/get_error_message.m b/matlab/get_error_message.m
index 8b463f619b5b0e6480a88d4cb5e3d2de386cfea8..57063a7d4bb2dc005b92d14b2fb05318945250c7 100644
--- a/matlab/get_error_message.m
+++ b/matlab/get_error_message.m
@@ -188,6 +188,8 @@ switch info(1)
message = 'Occbin: Simulation did not converge, increase maxit or check_ahead_periods.';
case 312
message = 'Occbin: Constraint(s) are binding at the end of the sample.';
+ case 313
+ message = 'Occbin: Simulation did not converge -- infinite loop of guess regimes';
case 320
message = 'Piecewise linear Kalman filter: There was a problem in obtaining the likelihood.';
case 401
diff --git a/tests/occbin/model_irrcap_twoconstraints/dynrbc_common.inc b/tests/occbin/model_irrcap_twoconstraints/dynrbc_common.inc
index 768dac606c5471ac7826bd1cb92d069703e21c60..d294ce5560e547eebf7a2a343aa6560404081ae1 100644
--- a/tests/occbin/model_irrcap_twoconstraints/dynrbc_common.inc
+++ b/tests/occbin/model_irrcap_twoconstraints/dynrbc_common.inc
@@ -44,13 +44,13 @@ if option==1
elseif option==2
- nperiods = 100;
+ nperiods = 90;
randn('seed',1);
shockssequence = 1*randn(nperiods,1)*0.02 ;
shocks(surprise,overwrite);
var erra;
- periods 1:100;
+ periods 1:90;
values (shockssequence);
end;
@@ -60,7 +60,7 @@ end
occbin_setup(smoother_debug);
occbin_solver(simul_debug,simul_periodic_solution,simul_periods=200,simul_maxit=300,simul_curb_retrench,simul_check_ahead_periods=300);
-occbin_write_regimes(filename='test',periods=[1:100]);
+occbin_write_regimes(filename='test',periods=[1:90]);
%% Modify to plot IRFs
titlelist = char('c','lambdak','k','i','a');