adding missing occbin functions

time-shift
Michel Juillard 2017-10-05 14:19:24 +02:00
parent a0ea6e3a21
commit 53b1910201
8 changed files with 325 additions and 0 deletions

View File

@ -0,0 +1,19 @@
function [r,g1,g2,g3] = evaluate_model(z,x,M,ss)
ll = M.lead_lag_incidence';
y = z(find(ll(:)));
switch nargout
case 1
r = feval([M.fname '_dynamic'],y,x, ...
M.params, ss, 1);
case 2
[r,g1] = feval([M.fname '_dynamic'],y,x, ...
M.params, ss, 1);
case 3
[r,g1,g2] = feval([M.fname '_dynamic'],y,x, ...
M.params, ss, 1);
case 4
[r,g1,g2,g3] = feval([M.fname '_dynamic'],y,x, ...
M.params, ss, 1);
end

26
matlab/occbin/get_coef.m Normal file
View File

@ -0,0 +1,26 @@
function [coef_y,coef_u] = get_coef(jacobian,M)
ll = M.lead_lag_incidence;
endo_nbr = M.endo_nbr;
coef_y = zeros(endo_nbr,3*endo_nbr);
coef_u = zeros(endo_nbr,M.exo_nbr);
if M.maximum_lag > 0
[junk,c1,c2] = find(ll(1,:));
coef_y(:,c1) = jacobian(:,c2);
[junk,c1,c2] = find(ll(2,:));
coef_y(:,c1+endo_nbr) = jacobian(:,c2);
if M.maximum_lead > 0
[junk,c1,c2] = find(ll(3,:));
coef_y(:,c1+2*endo_nbr) = jacobian(:,c2);
end
else
[junk,c1,c2] = find(ll(1,:));
coef_y(:,c1+endo_nbr) = jacobian(:,c2);
if M.maximum_lead > 0
[junk,c1,c2] = find(ll(2,:));
coef_y(:,c1+2*endo_nbr) = jacobian(:,c2);
end
end
coef_u = jacobian(:,max(ll(end,:))+1:end);

View File

@ -0,0 +1,72 @@
function [lb,ub,eq_index] = get_complementarity_conditions(M,ramsey_policy)
% Copyright (C) 2014 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
ub = inf(M.endo_nbr,1);
lb = -ub;
eq_index = (1:M.endo_nbr)';
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
lb(rc{i}{1}) = eval(rc{i}{3});
case {'<','<='}
ub(rc{i}{1}) = eval(rc{i}{3});
otherwise
error('Wrong operator in get_complementarity_conditions')
end
eq_index(i) = 1;
end
end
end
etags = M.equations_tags;
for i=1:size(etags,1)
if strcmp(etags{i,2},'mcp')
str = etags{i,3};
kop = strfind(etags{i,3},'<');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)),cellstr(M.endo_names)));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i,3},b{1}]))
end
ub(k) = str2num(str(kop+1:end));
eq_index(etags{i,1}) = k;
eq_index(k) = etags{i,1};
else
kop = strfind(etags{i,3},'>');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)),cellstr(M.endo_names)));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i},b{1}]))
end
lb(k) = str2num(str(kop+1:end));
eq_index(etags{i,1}) = k;
eq_index(k) = etags{i,1};
else
error(sprintf(['Complementarity condition %s can''t be ' ...
'parsed'],etags{i,3}))
end
end
end
end

View File

@ -0,0 +1,75 @@
function [ivar,ieq,lb,ub] = get_occbin_complementarity_conditions(M,ramsey_policy)
% Copyright (C) 2015 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
nrow = 1;
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
lb(nrow) = eval(rc{i}{3});
case {'<','<='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
ub(nrow) = eval(rc{i}{3});
otherwise
error('Wrong operator in get_complementarity_conditions')
end
nrow = nrow + 1;
end
end
end
etags = M.equations_tags;
for i=1:size(etags,1)
if strcmp(etags{i,2},'mcp')
str = etags{i,3};
kop = strfind(etags{i,3},'<');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)),cellstr(M.endo_names)));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i,3},b{1}]))
end
ivar(nrow) = k;
ieq(nrow) = etags{i,1};
ub(nrow) = eval(str(kop+1:end));
else
kop = strfind(etags{i,3},'>');
if ~isempty(kop)
k = find(strcmp(strtrim(str(1:kop-1)),cellstr(M.endo_names)));
if isempty(k)
error(sprintf(['Complementarity condition %s: variable %s is ' ...
'not recognized',etags{i},b{1}]))
end
ivar(nrow) = k;
ieq(nrow) = etags{i,1};
lb(k) = eval(str(kop+1:end));
else
error(sprintf(['Complementarity condition %s can''t be ' ...
'parsed'],etags{i,3}))
end
end
nrow = nrow + 1;
end
end

View File

@ -0,0 +1,88 @@
function [i_base,i_alt,c_base,c_alt] = get_occbin_constraints(M,steady_state,ramsey_policy)
% Copyright (C) 2015 Dynare Team
%
% This file is part of Dynare.
%
% Dynare is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% Dynare is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
nrow = 1;
if ramsey_policy
if isfield(M,'ramsey_model_constraints')
rc = M.ramsey_model_constraints;
for i = 1:length(rc)
switch rc{i}{2}
case {'>','>='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
lb(nrow) = eval(rc{i}{3});
case {'<','<='}
ivar(nrow) = rc{i}{1};
ieq(nrow) = rc{i}{1};
ub(nrow) = eval(rc{i}{3});
otherwise
error('Wrong operator in get_complementarity_conditions')
end
nrow = nrow + 1;
end
end
end
i_base = {};
i_alt = {};
etags = M.equations_tags;
m = 1;
base = true;
for i=1:size(etags,1)
[iv,boundary,operator] = parse_constraint(etags{i,3},M.endo_names,M.params,M.param_names);
if strcmp(etags{i,2},'OCCBIN')
if base
i_alt{m} = 1:M.eq_nbr;
i_alt{m}(etags{i,1}) = [];
c_base{m,1} = etags{i,1};
c_base{m,2} = iv;
c_base{m,3} = boundary - steady_state(iv);
c_base{m,4} = operator;
base = false;
else
i_base{m} = 1:M.eq_nbr;
i_base{m}(etags{i,1}) = [];
c_alt{m,1} = etags{i,1};
c_alt{m,2} = iv;
c_alt{m,3} = boundary - steady_state(iv);
c_alt{m,4} = operator;
base = true;
m = m + 1;
end
end
end
if ~base
error('OCCBIN: constraints must come by pair')
end
function [iv,boundary,operator] = parse_constraint(str,endo_names,params,param_names)
delim = {'<=','>=','<','>'};
[c,operator] = strsplit(str,delim);
operator = operator{1};
iv = strmatch(strtrim(c{1}),endo_names);
% try for a number
boundary = str2num(strtrim(c{2}));
% if not a number try for a parameter name
if isempty(boundary)
k = strmatch(strtrim(c{2}),param_names);
if isempty(k)
error(['OCCBIN: illegal constraint ' str]);
end
boundary = params(k);
end

View File

@ -0,0 +1,10 @@
function r = get_residuals(ivar,lb,ub,M,oo)
ss = oo.steady_state;
for i = 1:length(ivar)
% only one is different from zero
ss(ivar(i)) = lb(i) + ub(i);
end
oo.steady_state = ss;
r = evaluate_model(M,oo);

19
matlab/occbin/occbin.m Normal file
View File

@ -0,0 +1,19 @@
function [endo_simul,endo_simul_no_constraint,status] = occbin(M,oo,options)
% function oo=occbin(M,oo,options) solves linear models with occasionally
% binding constraints using OCCBIN by L. Guerrieri
status = 1;
constraint_nbr = sum(strcmp(upper(M.equations_tags(:,2)),'OCCBIN'))/2;
switch(constraint_nbr)
case 1
[zdatalinear_ zdatapiecewise_ zdatass_ oobase_ ] = ...
solve_one_constraint(M,oo,options);
case 2
[zdatalinear_ zdatapiecewise_ zdatass_ oobase_ ] = ...
solve_two_constraints(M,oo,options);
otherwise
error('OCCBIN can only handle two constraints in a model')
end
endo_simul = zdatapiecewise_';
endo_simul_no_constraint = zdatalinear_';

View File

@ -0,0 +1,16 @@
function b=test_constraint(x,constr)
b = zeros(size(x,1),size(constr,1));
for i=1:size(constr,1)
switch constr{i,4}
case '<'
b(:,i) = ~(x(:,constr{i,2}) < constr{i,3});
case '<='
b(:,i) = ~(x(:,constr{i,2}) <= constr{i,3});
case '>'
b(:,i) = ~(x(:,constr{i,2}) > constr{i,3});
case '>='
b(:,i) = ~(x(:,constr{i,2}) >= constr{i,3});
otherwise
error('OCCBIN: wrong inequality sign')
end
end