From 444d0a6b349b64773fe6429368f2e54b49d70681 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org> Date: Wed, 15 May 2024 18:23:45 +0200 Subject: [PATCH] Preprocessor: refactor output for complementarity conditions --- matlab/display_static_residuals.m | 14 ++-- matlab/ep/extended_path_initialization.m | 5 +- matlab/evaluate_steady_state.m | 4 +- matlab/lmmcp/get_complementarity_conditions.m | 84 ------------------- .../perfect_foresight_solver_core.m | 10 +-- .../private/initialize_stacked_problem.m | 6 +- .../solve_stacked_problem.m | 6 +- preprocessor | 2 +- .../Gali_2015_chapter_3_nonlinear.mod | 8 +- 9 files changed, 25 insertions(+), 114 deletions(-) delete mode 100644 matlab/lmmcp/get_complementarity_conditions.m diff --git a/matlab/display_static_residuals.m b/matlab/display_static_residuals.m index d6b62e81e5..331d4a67a8 100644 --- a/matlab/display_static_residuals.m +++ b/matlab/display_static_residuals.m @@ -15,7 +15,7 @@ function z = display_static_residuals(M_, options_, oo_,options_resid_) % SPECIAL REQUIREMENTS % none -% Copyright © 2001-2023 Dynare Team +% Copyright © 2001-2024 Dynare Team % % This file is part of Dynare. % @@ -60,15 +60,15 @@ z = evaluate_static_model(oo_.steady_state, [oo_.exo_steady_state; ... oo_.exo_det_steady_state], ... M_.params, M_, options_); +[lb, ub] = feval(sprintf('%s.static_complementarity_conditions', M_.fname), M_.params); if ismember(options_.solve_algo,[10,11]) - [lb,ub,eq_index] = get_complementarity_conditions(M_,options_.ramsey_policy); eq_to_check=find(isfinite(lb) | isfinite(ub)); eq_to_ignore=eq_to_check(oo_.steady_state(eq_to_check,:)<=lb(eq_to_check)+eps | oo_.steady_state(eq_to_check,:)>=ub(eq_to_check)-eps); - z(eq_index(eq_to_ignore))=0; - disp_string=' (accounting for MCP tags)'; + z(M_.static_mcp_equations_reordering(eq_to_ignore))=0; + disp_string=' (accounting for complementarity conditions)'; else - if istag && ~isempty(strmatch('mcp',M_.equations_tags(:,2),'exact')) - disp_string=' (ignoring MCP tags)'; + if any(isfinite(lb) | isfinite(ub)) + disp_string=' (ignoring complementarity conditions)'; else disp_string=''; end @@ -124,4 +124,4 @@ if nargout == 0 disp('All residuals are zero') end skipline(2) -end \ No newline at end of file +end diff --git a/matlab/ep/extended_path_initialization.m b/matlab/ep/extended_path_initialization.m index 3b9f02477c..362788dec9 100644 --- a/matlab/ep/extended_path_initialization.m +++ b/matlab/ep/extended_path_initialization.m @@ -16,7 +16,7 @@ function [initial_conditions, innovations, pfm, ep, verbosity, options_, oo_] = % % SPECIAL REQUIREMENTS -% Copyright © 2016-2023 Dynare Team +% Copyright © 2016-2024 Dynare Team % % This file is part of Dynare. % @@ -126,7 +126,8 @@ else end % set boundaries if mcp -[lb,ub,pfm.eq_index] = get_complementarity_conditions(M_, options_.ramsey_policy); +[lb, ub] = feval(sprintf('%s.dynamic_complementarity_conditions', M_.fname), M_.params); +pfm.eq_index = M_.dynamic_mcp_equations_reordering; if options_.ep.solve_algo == 10 options_.lmmcp.lb = repmat(lb,block_nbr,1); options_.lmmcp.ub = repmat(ub,block_nbr,1); diff --git a/matlab/evaluate_steady_state.m b/matlab/evaluate_steady_state.m index 6e78515fa2..047bdc71de 100644 --- a/matlab/evaluate_steady_state.m +++ b/matlab/evaluate_steady_state.m @@ -280,7 +280,7 @@ elseif ~options_.bytecode && ~options_.block if ~options_.linear % non linear model if ismember(options_.solve_algo,[10,11]) - [lb,ub,eq_index] = get_complementarity_conditions(M_,options_.ramsey_policy); + [lb, ub] = feval(sprintf('%s.static_complementarity_conditions', M_.fname), M_.params); if options_.solve_algo == 10 options_.lmmcp.lb = lb; options_.lmmcp.ub = ub; @@ -293,7 +293,7 @@ elseif ~options_.bytecode && ~options_.block options_.steady.maxit, options_.solve_tolf, options_.solve_tolx, ... options_, exo_ss, params,... M_.endo_nbr, static_resid, static_g1, ... - M_.static_g1_sparse_rowval, M_.static_g1_sparse_colval, M_.static_g1_sparse_colptr, eq_index); + M_.static_g1_sparse_rowval, M_.static_g1_sparse_colval, M_.static_g1_sparse_colptr, M_.static_mcp_equations_reordering); else [ys, check, fvec, fjac, errorcode] = dynare_solve(@static_problem, ys_init, ... options_.steady.maxit, options_.solve_tolf, options_.solve_tolx, ... diff --git a/matlab/lmmcp/get_complementarity_conditions.m b/matlab/lmmcp/get_complementarity_conditions.m deleted file mode 100644 index 9aed55f4d3..0000000000 --- a/matlab/lmmcp/get_complementarity_conditions.m +++ /dev/null @@ -1,84 +0,0 @@ -function [lb,ub,eq_index] = get_complementarity_conditions(M_,ramsey_policy) -% [lb,ub,eq_index] = get_complementarity_conditions(M_,ramsey_policy) -% INPUTS -% - M [struct] contains a description of the model. -% - ramsey_policy [boolean] indicator whether a Ramsey problem is considered -% OUTPUTS -% - lb [double] endo_nbr*1 array of lower bounds for -% endogenous variables -% - ub [double] endo_nbr*1 array of upper bounds for -% endogenous variables -% - eq_index [struct] endo_nbr*1 index vector describing residual mapping resulting -% from complementarity setup used in -% perfect_foresight_mcp_problem.m - -% Copyright © 2014-2024 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 <https://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 - end - end -end - -etags = M_.equations_tags; -for i=1:size(etags,1) - if strcmp(etags{i,2},'mcp') - eq_nbr = etags{i,1}; - str = etags{i,3}; - kop = strfind(etags{i,3},'<'); - if ~isempty(kop) - k = find(strcmp(strtrim(str(1:kop-1)), M_.endo_names)); %get variable index with restriction - if isempty(k) - error('Complementarity condition %s: variable %s is not recognized', etags{i,3}, strtrim(str(1:kop-1))) - end - ub(k) = str2num(str(kop+1:end)); - swapped_eq_nbr = find(eq_index == eq_nbr); - assert(length(swapped_eq_nbr) == 1); - eq_index(swapped_eq_nbr) = eq_index(k); - eq_index(k) = eq_nbr; - else - kop = strfind(etags{i,3},'>'); - if ~isempty(kop) - k = find(strcmp(strtrim(str(1:kop-1)), M_.endo_names)); %get variable index with restriction - if isempty(k) - error('Complementarity condition %s: variable %s is not recognized', etags{i,3}, strtrim(str(1:kop-1))) - end - lb(k) = str2num(str(kop+1:end)); - swapped_eq_nbr = find(eq_index == eq_nbr); - assert(length(swapped_eq_nbr) == 1); - eq_index(swapped_eq_nbr) = eq_index(k); - eq_index(k) = eq_nbr; - else - error('Complementarity condition %s can''t be parsed',etags{i,3}) - end - end - end -end diff --git a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m index 2c7202f6a1..fdeaa60e4b 100644 --- a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m +++ b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m @@ -112,14 +112,8 @@ else [y, success, maxerror, iter] = sim1_lbj(y, exo_simul, steady_state, M_, options_); case 7 if options_.linear_approximation - if isequal(options_.solve_algo, 10) - if options_.ramsey_policy && isfield(M_,'ramsey_model_constraints') && ~isempty(M_.ramsey_model_constraints) - warning('Due to ramsey_constraints you should not specify your model as model(linear)!') - elseif options_.lmmcp.status - warning('Due to lmmcp option, you should not specify your model as model(linear)!') - else - warning('It would be more efficient to set option solve_algo equal to 0!') - end + if options_.solve_algo == 10 || options_.solve_algo == 11 + warning('Since you are requesting an MCP solver, you should not specify your model as model(linear)!') end [y, success] = solve_stacked_linear_problem(y, exo_simul, steady_state, exo_steady_state, M_, options_); else diff --git a/matlab/perfect-foresight-models/private/initialize_stacked_problem.m b/matlab/perfect-foresight-models/private/initialize_stacked_problem.m index 333805e41d..5b6c00b6c3 100644 --- a/matlab/perfect-foresight-models/private/initialize_stacked_problem.m +++ b/matlab/perfect-foresight-models/private/initialize_stacked_problem.m @@ -30,7 +30,7 @@ function [options_, y0, yT, z, i_cols, i_cols_J1, i_cols_T, i_cols_j, i_cols_1, % - i_cols_J0 [double] indices of contemporaneous variables appearing in M_.lead_lag_incidence (relevant in problems with periods=1) % - dynamicmodel [handle] function handle to _dynamic-file -% Copyright © 2015-2023 Dynare Team +% Copyright © 2015-2024 Dynare Team % % This file is part of Dynare. % @@ -50,7 +50,7 @@ function [options_, y0, yT, z, i_cols, i_cols_J1, i_cols_T, i_cols_j, i_cols_1, periods = options_.periods; if (options_.solve_algo == 10) if ~isfield(options_.lmmcp,'lb') - [lb,ub,pfm.eq_index] = get_complementarity_conditions(M_,options_.ramsey_policy); + [lb, ub] = feval(sprintf('%s.dynamic_complementarity_conditions', M_.fname), M_.params); if options_.linear_approximation lb = lb - steadystate_y; ub = ub - steadystate_y; @@ -60,7 +60,7 @@ if (options_.solve_algo == 10) end elseif (options_.solve_algo == 11) if ~isfield(options_.mcppath,'lb') - [lb,ub,pfm.eq_index] = get_complementarity_conditions(M_,options_.ramsey_policy); + [lb, ub] = feval(sprintf('%s.dynamic_complementarity_conditions', M_.fname), M_.params); if options_.linear_approximation lb = lb - steadystate_y; ub = ub - steadystate_y; diff --git a/matlab/perfect-foresight-models/solve_stacked_problem.m b/matlab/perfect-foresight-models/solve_stacked_problem.m index 2a5d0bbcd4..cd5996343e 100644 --- a/matlab/perfect-foresight-models/solve_stacked_problem.m +++ b/matlab/perfect-foresight-models/solve_stacked_problem.m @@ -14,7 +14,7 @@ function [endogenousvariables, success, maxerror] = solve_stacked_problem(endoge % - success [logical] Whether a solution was found % - maxerror [double] 1-norm of the residual -% Copyright © 2015-2023 Dynare Team +% Copyright © 2015-2024 Dynare Team % % This file is part of Dynare. % @@ -35,7 +35,7 @@ function [endogenousvariables, success, maxerror] = solve_stacked_problem(endoge initialize_stacked_problem(endogenousvariables, options_, M_, steadystate); if (options_.solve_algo == 10 || options_.solve_algo == 11)% mixed complementarity problem - [lb,ub,eq_index] = get_complementarity_conditions(M_,options_.ramsey_policy); + [lb, ub] = feval(sprintf('%s.dynamic_complementarity_conditions', M_.fname), M_.params); if options_.linear_approximation lb = lb - steadystate_y; ub = ub - steadystate_y; @@ -54,7 +54,7 @@ if (options_.solve_algo == 10 || options_.solve_algo == 11)% mixed complementari exogenousvariables, M_.params, steadystate, ... M_.maximum_lag, options_.periods, M_.endo_nbr, i_cols, ... i_cols_J1, i_cols_1, i_cols_T, i_cols_j, i_cols_0, i_cols_J0, ... - eq_index); + M_.dynamic_mcp_equations_reordering); eq_to_ignore=find(isfinite(lb) | isfinite(ub)); else diff --git a/preprocessor b/preprocessor index 097f706c55..dc1ec15fc6 160000 --- a/preprocessor +++ b/preprocessor @@ -1 +1 @@ -Subproject commit 097f706c55b56d63a26781b38d8c90d48d7b2589 +Subproject commit dc1ec15fc6fec59cc4c70fdb719bf731d7a57c9c diff --git a/tests/discretionary_policy/Gali_2015_chapter_3_nonlinear.mod b/tests/discretionary_policy/Gali_2015_chapter_3_nonlinear.mod index ed443c95a2..cccd971004 100644 --- a/tests/discretionary_policy/Gali_2015_chapter_3_nonlinear.mod +++ b/tests/discretionary_policy/Gali_2015_chapter_3_nonlinear.mod @@ -19,8 +19,8 @@ */ /* - * Copyright © 2016-20 Johannes Pfeifer - * Copyright © 2020 Dynare Team + * Copyright © 2016-2020 Johannes Pfeifer + * Copyright © 2020-2024 Dynare Team * * This is free software: you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by @@ -162,7 +162,7 @@ model; log_A=log(A); [name='Definition log preference'] log_Z=log(Z); - [mcp='a'] + y_hat=log(Y)-STEADY_STATE(log(Y)); pi=log(Pi)-STEADY_STATE(log(Pi)); @@ -228,4 +228,4 @@ end if max(max(abs([temp.oo_.irfs.y_eps_a; temp.oo_.irfs.w_real_eps_a; temp.oo_.irfs.n_eps_a; temp.oo_.irfs.pi_ann_eps_a]-... [oo_.irfs.log_y_eps_a; oo_.irfs.log_W_real_eps_a; oo_.irfs.log_N_eps_a; oo_.irfs.pi_ann_eps_a])))>1e-6 error('Policy is different') -end \ No newline at end of file +end -- GitLab