dyn_ramsey_static.m 7.15 KB
Newer Older
1
function [steady_state,params,check] = dyn_ramsey_static(ys_init,M,options_,oo)
2

3
4
% function  [steady_state,params,check] = dyn_ramsey_static_(ys_init,M,options_,oo)
% Computes the steady state for optimal policy
5
6
%
% INPUTS
7
8
9
10
%    ys_init:       vector of endogenous variables or instruments
%    M:             Dynare model structure
%    options:       Dynare options structure
%    oo:            Dynare results structure
11
%
12
% OUTPUTS
13
14
15
%    steady_state:  steady state value
%    params:        parameters at steady state, potentially updated by
%                   steady_state file
16
%    check:         error indicator, 0 if everything is OK
17
18
19
20
%
% SPECIAL REQUIREMENTS
%    none

21
% Copyright (C) 2003-2018 Dynare Team
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
%
% 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/>.
37
38
39
40


params = M.params;
check = 0;
41
options_.steadystate.nocheck = 1; %locally disable checking because Lagrange multipliers are not accounted for in evaluate_steady_state_file
42
                                  % dyn_ramsey_static_1 is a subfunction
43
44
nl_func = @(x) dyn_ramsey_static_1(x,M,options_,oo);

45
% check_static_model is a subfunction
46
47
if check_static_model(ys_init,M,options_,oo) && ~options_.steadystate_flag
    steady_state = ys_init;
48
49
    return
elseif options_.steadystate_flag
50
    k_inst = [];
MichelJuillard's avatar
MichelJuillard committed
51
    inst_nbr = size(options_.instruments,1);
52
    for i = 1:inst_nbr
53
        k_inst = [k_inst; strmatch(options_.instruments{i}, M.endo_names, 'exact')];
54
55
    end
    if inst_nbr == 1
56
        %solve for instrument, using univariate solver, starting at initial value for instrument
57
58
59
60
61
62
        [inst_val, info1]= csolve(nl_func,ys_init(k_inst),'',options_.solve_tolf,options_.ramsey.maxit);
        if info1==1 || info1==3
            check=81;
        end
        if info1==4
            check=87;
63
        end
64
    else
65
66
67
68
69
70
        %solve for instrument, using multivariate solver, starting at
        %initial value for instrument
        opt = options_;
        opt.jacobian_flag = 0;
        [inst_val,info1] = dynare_solve(nl_func,ys_init(k_inst), ...
                                        opt);
71
72
73
        if info1~=0
            check=81;
        end
74
    end
75
    ys_init(k_inst) = inst_val;
76
    exo_ss = [oo.exo_steady_state oo.exo_det_steady_state];
77
    [xx,params] = evaluate_steady_state_file(ys_init,exo_ss,M,options_,~options_.steadystate.nocheck); %run steady state file again to update parameters
78
    [~,~,steady_state] = nl_func(inst_val); %compute and return steady state
79
else
80
    n_var = M.orig_endo_nbr;
81
    xx = oo.steady_state(1:n_var);
82
83
    opt = options_;
    opt.jacobian_flag = 0;
84
85
86
87
    [xx,info1] = dynare_solve(nl_func,xx,opt);
    if info1~=0
        check=81;
    end
88
    [~,~,steady_state] = nl_func(xx);
89
90
91
92
end



93
function [resids,rJ,steady_state] = dyn_ramsey_static_1(x,M,options_,oo)
94
95
96
resids = [];
rJ = [];
mult = [];
97
98

% recovering usefull fields
99
params = M.params;
100
endo_nbr = M.endo_nbr;
101
endo_names = M.endo_names;
102
orig_endo_nbr = M.orig_endo_nbr;
103
aux_vars_type = [M.aux_vars.type];
104
orig_endo_aux_nbr = orig_endo_nbr + min(find(aux_vars_type == 6)) - 1;
105
orig_eq_nbr = M.orig_eq_nbr;
106
inst_nbr = orig_endo_aux_nbr - orig_eq_nbr;
107
% indices of Lagrange multipliers
108
fname = M.fname;
109
110
111
112
113
114


if options_.steadystate_flag
    k_inst = [];
    instruments = options_.instruments;
    for i = 1:size(instruments,1)
115
        k_inst = [k_inst; strmatch(instruments{i}, endo_names, 'exact')];
116
    end
117
118
    ys_init=zeros(size(oo.steady_state)); %create starting vector for steady state computation as only instrument value is handed over
    ys_init(k_inst) = x; %set instrument, the only value required for steady state computation, to current value
Houtan Bastani's avatar
Houtan Bastani committed
119
    [x,params,check] = evaluate_steady_state_file(ys_init,... %returned x now has size endo_nbr as opposed to input size of n_instruments
120
                                                  [oo.exo_steady_state; ...
121
                        oo.exo_det_steady_state], ...
122
                                                  M,options_,~options_.steadystate.nocheck);
123
124
125
    if any(imag(x(1:M.orig_endo_nbr))) %return with penalty
        resids=1+sum(abs(imag(x(1:M.orig_endo_nbr)))); %return with penalty
        steady_state=NaN(endo_nbr,1);
126
        return
127
128
    end

129
end
MichelJuillard's avatar
MichelJuillard committed
130

131
132
133
134
135
xx = zeros(endo_nbr,1); %initialize steady state vector
xx(1:M.orig_endo_nbr) = x(1:M.orig_endo_nbr); %set values of original endogenous variables based on steady state file or initial value

% setting steady state of auxiliary variables that depends on original endogenous variables
if any([M.aux_vars.type] ~= 6) %auxiliary variables other than multipliers
136
    needs_set_auxiliary_variables = 1;
137
    if M.set_auxiliary_variables
138
        fh = str2func([M.fname '.set_auxiliary_variables']);
139
140
141
142
143
144
145
        s_a_v_func = @(z) fh(z,...
            [oo.exo_steady_state,...
            oo.exo_det_steady_state],...
            params);
    else
        s_a_v_func = z;
    end
146
    xx = s_a_v_func(xx);
147
148
149
else
    needs_set_auxiliary_variables = 0;
end
150
151

% value and Jacobian of objective function
152
ex = zeros(1,M.exo_nbr);
153
[U,Uy,Uyy] = feval([fname '.objective.static'],x,ex, params);
154
155
Uyy = reshape(Uyy,endo_nbr,endo_nbr);

156
157
% set multipliers and auxiliary variables that
% depends on multipliers to 0 to compute residuals
158
if (options_.bytecode)
159
160
161
    [chck, res, junk] = bytecode('static',xx,[oo.exo_steady_state oo.exo_det_steady_state], ...
                                 params, 'evaluate');
    fJ = junk.g1;
162
else
163
    [res,fJ] = feval([fname '.static'],xx,[oo.exo_steady_state oo.exo_det_steady_state], ...
164
                     params);
165
end
166
167
168
% index of multipliers and corresponding equations
% the auxiliary variables before the Lagrange multipliers are treated
% as ordinary endogenous variables
169
aux_eq = [1:orig_endo_aux_nbr, orig_endo_aux_nbr+orig_eq_nbr+1:size(fJ,1)];
170
171
A = fJ(1:orig_endo_aux_nbr,orig_endo_nbr+find(aux_vars_type==6));
y = res(1:orig_endo_aux_nbr);
172
mult = -A\y;
173

174
resids1 = y+A*mult;
MichelJuillard's avatar
MichelJuillard committed
175
176
177
178
if inst_nbr == 1
    r1 = sqrt(resids1'*resids1);
else
    [q,r,e] = qr([A y]');
179
180
    k = size(A,1)+(1-inst_nbr:0);
    r1 = r(end,k)';
MichelJuillard's avatar
MichelJuillard committed
181
end
182
if options_.steadystate_flag
MichelJuillard's avatar
MichelJuillard committed
183
    resids = r1;
184
else
185
    resids = [res(orig_endo_nbr+(1:orig_endo_nbr-inst_nbr)); r1];
186
187
end
rJ = [];
188
189
190
191
192
193
194
195
196
if needs_set_auxiliary_variables
    steady_state = s_a_v_func([xx(1:orig_endo_aux_nbr); mult]);
else
    steady_state = [xx(1:orig_endo_aux_nbr); mult];
end

function result = check_static_model(ys,M,options_,oo)
result = false;
if (options_.bytecode)
197
198
    [chck, res, ~] = bytecode('static',ys,[oo.exo_steady_state oo.exo_det_steady_state], ...
                              M.params, 'evaluate');
199
else
200
    res = feval([M.fname '.static'],ys,[oo.exo_steady_state oo.exo_det_steady_state], ...
201
202
203
204
205
                M.params);
end
if norm(res) < options_.solve_tolf
    result = true;
end