diff --git a/matlab/+occbin/solve_one_constraint.m b/matlab/+occbin/solve_one_constraint.m
index 27d01a7204148a977600be3fbf6b135cc45beb63..5b710cff5a467e49b573044a2295ca3c4c9343a5 100644
--- a/matlab/+occbin/solve_one_constraint.m
+++ b/matlab/+occbin/solve_one_constraint.m
@@ -125,20 +125,24 @@ 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
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_periods
% enforce endogenously increased nperiods to not violate max_periods
@@ -201,6 +205,7 @@ for shock_period = 1:n_shocks_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;
@@ -220,16 +225,43 @@ for shock_period = 1:n_shocks_periods
binding_indicator= (binding_indicator | binding.constraint_1) & ~(binding_indicator & relax.constraint_1);
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
@@ -281,9 +313,14 @@ for shock_period = 1:n_shocks_periods
return
end
else
+ if is_periodic_loop
+ disp_verbose('Did not converge -- infinite loop of regimes.',opts_simul_.debug)
+ error_flag = 313;
+ else
disp_verbose('Did not converge -- increase maxit.',opts_simul_.debug)
- if opts_simul_.waitbar; dyn_waitbar_close(hh); end
error_flag = 311;
+ end
+ if opts_simul_.waitbar; dyn_waitbar_close(hh); end
return
end
else
diff --git a/matlab/+occbin/solve_two_constraints.m b/matlab/+occbin/solve_two_constraints.m
index 4e448b241f89dc72c8d9fafd11bd05a1c3e59a1d..15c8ae1c2a819ff9c9d8daf72d364430f3800c42 100644
--- a/matlab/+occbin/solve_two_constraints.m
+++ b/matlab/+occbin/solve_two_constraints.m
@@ -133,20 +133,24 @@ 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
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_periods
% enforce endogenously increased nperiods to not violate max_periods
@@ -248,17 +252,43 @@ for shock_period = 1:n_shocks_periods
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
@@ -312,8 +342,13 @@ for shock_period = 1:n_shocks_periods
return;
end
else
+ if is_periodic_loop
+ disp_verbose('Did not converge -- infinite loop of 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/get_error_message.m b/matlab/get_error_message.m
index 1013184d69aeabf1c1434238082a4a718742bfee..cb236cfd2bfbb347cb32387317b25eae0b9848de 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 regimes';
case 320
message = 'Piecewise linear Kalman filter: There was a problem in obtaining the likelihood.';
otherwise