diff --git a/matlab/perfect-foresight-models/perfect_foresight_solver.m b/matlab/perfect-foresight-models/perfect_foresight_solver.m index cc4f2d60a0a5067a61b482abfed15bf1a2f2659a..f55b6e3e2fa82cbad9452e97f54cdfe7e3f36678 100644 --- a/matlab/perfect-foresight-models/perfect_foresight_solver.m +++ b/matlab/perfect-foresight-models/perfect_foresight_solver.m @@ -23,7 +23,7 @@ function [oo_, ts]=perfect_foresight_solver(M_, options_, oo_, no_error_if_learn % SPECIAL REQUIREMENTS % none -% Copyright © 1996-2024 Dynare Team +% Copyright © 1996-2025 Dynare Team % % This file is part of Dynare. % @@ -158,13 +158,19 @@ if completed_share == 1 if options_.simul.endval_steady oo_.steady_state = steady_state; end - % NB: no need to modify oo_.exo_simul and oo_.exo_steady_state, since we simulated 100% of the shock + % NB: no need to modify oo_.exo_simul and oo_.exo_steady_state, since we simulated 100% of the shock, except if controlled paths + if ~isempty(M_.perfect_foresight_controlled_paths) + oo_.exo_simul = exo_simul; + end if ~options_.noprint fprintf('Perfect foresight solution found.\n\n') end oo_.deterministic_simulation.status = true; elseif options_.simul.homotopy_linearization_fallback && completed_share > 0 + if ~isempty(M_.perfect_foresight_controlled_paths) + error('Option homotopy_linearization_fallback not available with perfect_foresight_controlled_paths') + end oo_.deterministic_simulation.sim1.endo_simul = endo_simul; oo_.deterministic_simulation.sim1.exo_simul = exo_simul; oo_.deterministic_simulation.sim1.steady_state = steady_state; @@ -188,6 +194,9 @@ elseif options_.simul.homotopy_linearization_fallback && completed_share > 0 oo_.deterministic_simulation.status = true; oo_.deterministic_simulation.homotopy_linearization = true; elseif options_.simul.homotopy_marginal_linearization_fallback > 0 && completed_share > options_.simul.homotopy_marginal_linearization_fallback + if ~isempty(M_.perfect_foresight_controlled_paths) + error('Option homotopy_marginal_linearization_fallback not available with perfect_foresight_controlled_paths') + end oo_.deterministic_simulation.sim1.endo_simul = endo_simul; oo_.deterministic_simulation.sim1.exo_simul = exo_simul; oo_.deterministic_simulation.sim1.steady_state = steady_state; @@ -344,7 +353,7 @@ while step > options_.simul.homotopy_min_step_size end % Solve for the paths of the endogenous variables. - [endo_simul, success, maxerror, solver_iter, per_block_status] = perfect_foresight_solver_core(endo_simul, exo_simul, steady_state, exo_steady_state, M_, options_); + [endo_simul, success, maxerror, solver_iter, per_block_status, exo_simul] = perfect_foresight_solver_core(endo_simul, exo_simul, steady_state, exo_steady_state, M_, options_); else success = false; maxerror = NaN; @@ -362,6 +371,10 @@ while step > options_.simul.homotopy_min_step_size break end + if ~isempty(M_.perfect_foresight_controlled_paths) + error('Homotopy not available with perfect_foresight_controlled_paths') + end + if iteration == 1 && ~options_.noprint fprintf('\nEntering the homotopy method iterations...\n') iter_summary_table = { sprintf('\nIter. \t | Share \t | Status \t | Max. residual\t | Duration (sec)\n'), diff --git a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m index 05d20af691ee238380d258fadd288023224694f5..8d38504c3bbe71dad11095f20833e47d36be6906 100644 --- a/matlab/perfect-foresight-models/perfect_foresight_solver_core.m +++ b/matlab/perfect-foresight-models/perfect_foresight_solver_core.m @@ -1,4 +1,4 @@ -function [y, success, maxerror, iter, per_block_status] = perfect_foresight_solver_core(y, exo_simul, steady_state, exo_steady_state, M_, options_) +function [y, success, maxerror, iter, per_block_status, exo_simul] = perfect_foresight_solver_core(y, exo_simul, steady_state, exo_steady_state, M_, options_) % Core function calling solvers for perfect foresight model % @@ -17,7 +17,7 @@ function [y, success, maxerror, iter, per_block_status] = perfect_foresight_solv % - iter [integer] Number of iterations of the underlying nonlinear solver (empty for non-iterative methods) % - per_block_status [struct] In the case of block decomposition, provides per-block solver status information (empty if no block decomposition) -% Copyright © 2015-2024 Dynare Team +% Copyright © 2015-2025 Dynare Team % % This file is part of Dynare. % @@ -55,6 +55,25 @@ if options_.linear && (isequal(options_.stack_solve_algo, 0) || isequal(options_ options_.linear_approximation = true; end +if ~isempty(M_.perfect_foresight_controlled_paths) + assert(nargout >= 6); % Ensure modified exos are used + if options_.bytecode + error('perfect_foresight_controlled_paths is not compatible with bytecode option') + end + if options_.block + error('perfect_foresight_controlled_paths is not compatible with block option') + end + if options_.linear + error('perfect_foresight_controlled_paths is not compatible with linear option') + end + if ~ismember(options_.stack_solve_algo, [0 2 3]) + error('perfect_foresight_controlled_paths is only available with stack_solve_algo equal to 0, 2 or 3') + end + if M_.maximum_endo_lead == 0 || M_.maximum_endo_lag == 0 + error('perfect_foresight_controlled_paths is not compatible with purely backward, purely forward or static models') + end +end + maxerror = []; iter = []; per_block_status = []; @@ -103,7 +122,7 @@ else end [y, success, maxerror] = sim1_linear(y, exo_simul, steady_state, exo_steady_state, M_, options_); else - [y, success, maxerror, iter] = sim1(y, exo_simul, steady_state, M_, options_); + [y, success, maxerror, iter, exo_simul] = sim1(y, exo_simul, steady_state, M_, options_); end case {1 6} if options_.linear_approximation diff --git a/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m b/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m index 22a07c276f41b31940203543427607236ee170a6..c69029c9d90ffc4ff1949eb2fa88a3f21e8ec5ef 100644 --- a/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m +++ b/matlab/perfect-foresight-models/perfect_foresight_with_expectation_errors_solver.m @@ -7,7 +7,7 @@ function [oo_, ts] = perfect_foresight_with_expectation_errors_solver(M_, option % OUTPUTS % oo_ [structure] storing the results -% Copyright © 2021-2024 Dynare Team +% Copyright © 2021-2025 Dynare Team % % This file is part of Dynare. % @@ -28,6 +28,10 @@ if options_.pfwee.constant_simulation_length && ~isempty(options_.simul.last_sim error('Options constant_simulation_length and last_simulation_period cannot be used together') end +if ~isempty(M_.perfect_foresight_controlled_paths) + error('perfect_foresight_with_expectation_errors_solver not available with perfect_foresight_controlled_paths') +end + [periods, first_simulation_period] = get_simulation_periods(options_); % Retrieve initial paths built by pfwee_setup diff --git a/matlab/perfect-foresight-models/sim1.m b/matlab/perfect-foresight-models/sim1.m index 55f5f6f63fd658761d110a57119de0ae52dbe7bd..d9d7f5a39f84784b2308e2e10dc49a57dba055fc 100644 --- a/matlab/perfect-foresight-models/sim1.m +++ b/matlab/perfect-foresight-models/sim1.m @@ -1,4 +1,4 @@ -function [endogenousvariables, success, err, iter] = sim1(endogenousvariables, exogenousvariables, steadystate, M_, options_) +function [endogenousvariables, success, err, iter, exogenousvariables] = sim1(endogenousvariables, exogenousvariables, steadystate, M_, options_) % [endogenousvariables, success, err, iter] = sim1(endogenousvariables, exogenousvariables, steadystate, M_, options_) % Performs deterministic simulations with lead or lag of one period, using % a basic Newton solver on sparse matrices. @@ -6,7 +6,7 @@ function [endogenousvariables, success, err, iter] = sim1(endogenousvariables, e % % INPUTS % - endogenousvariables [double] N*(T+M_.maximum_lag+M_.maximum_lead) array, paths for the endogenous variables (initial condition + initial guess + terminal condition). -% - exogenousvariables [double] T*M array, paths for the exogenous variables. +% - exogenousvariables [double] (T+M_.maximum_lag+M_.maximum_lead)*M array, paths for the exogenous variables. % - steadystate [double] N*1 array, steady state for the endogenous variables. % - M_ [struct] contains a description of the model. % - options_ [struct] contains various options. @@ -15,8 +15,10 @@ function [endogenousvariables, success, err, iter] = sim1(endogenousvariables, e % - success [logical] Whether a solution was found % - err [double] ∞-norm of the residual % - iter [integer] Number of iterations +% - exogenousvariables [double] (T+M_.maximum_lag+M_.maximum_lead)*M array, paths for the exogenous variables +% (may be modified if perfect_foresight_controlled_paths present) -% Copyright © 1996-2024 Dynare Team +% Copyright © 1996-2025 Dynare Team % % This file is part of Dynare. % @@ -36,9 +38,26 @@ function [endogenousvariables, success, err, iter] = sim1(endogenousvariables, e verbose = options_.verbosity && ~options_.noprint; ny = M_.endo_nbr; -periods = get_simulation_periods(options_); +[periods, first_simulation_period] = get_simulation_periods(options_); vperiods = periods*ones(1,options_.simul.maxit); +if ~isempty(M_.perfect_foresight_controlled_paths) + for i=1:length(M_.perfect_foresight_controlled_paths) + prange = M_.perfect_foresight_controlled_paths(i).periods; + if isa(prange, 'dates') + prange = prange - first_simulation_period + 1; + end + endogenousvariables(M_.perfect_foresight_controlled_paths(i).exogenize_id, prange+M_.maximum_lag) = M_.perfect_foresight_controlled_paths(i).value; + end + + if options_.debug + error('Debugging not available with perfect_foresight_controlled_paths') + end + if options_.endogenous_terminal_period + error('The endogenous_terminal_period option not available with perfect_foresight_controlled_paths') + end +end + if M_.maximum_lag > 0 y0 = endogenousvariables(:, M_.maximum_lag); else @@ -70,6 +89,11 @@ for iter = 1:options_.simul.maxit [res, A] = perfect_foresight_problem(y, y0, yT, exogenousvariables, M_.params, steadystate, periods, M_, options_); % A is the stacked Jacobian with period x equations alongs the rows and % periods times variables (in declaration order) along the columns + + if ~isempty(M_.perfect_foresight_controlled_paths) + A = stacked_jacobian_substitute_controlled_paths(A, y, y0, yT, exogenousvariables, steadystate, M_, options_); + end + if options_.debug && iter==1 [row,col]=find(A); row=setdiff(1:periods*ny,row); @@ -136,6 +160,19 @@ for iter = 1:options_.simul.maxit end end y = y + dy; + + if ~isempty(M_.perfect_foresight_controlled_paths) + for i = 1:length(M_.perfect_foresight_controlled_paths) + prange = M_.perfect_foresight_controlled_paths(i).periods; + if isa(prange, 'dates') + prange = prange - first_simulation_period + 1; + end + endogenize_id = M_.perfect_foresight_controlled_paths(i).endogenize_id; + exogenize_id = M_.perfect_foresight_controlled_paths(i).exogenize_id; + y(exogenize_id+(prange-1)*M_.endo_nbr) = M_.perfect_foresight_controlled_paths(i).value; + exogenousvariables(prange+M_.maximum_lag,endogenize_id) = exogenousvariables(prange+M_.maximum_lag,endogenize_id) + dy(exogenize_id+(prange-1)*M_.endo_nbr); + end + end end endogenousvariables(:, M_.maximum_lag+(1:periods)) = reshape(y, ny, periods); @@ -278,7 +315,7 @@ end function display_critical_variables(dyy, M_, noprint) -if noprint +if noprint || ~isempty(M_.perfect_foresight_controlled_paths) return end diff --git a/matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m b/matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m new file mode 100644 index 0000000000000000000000000000000000000000..cf68ec6aa19d03da07b3b6e8703d5a6acbca560c --- /dev/null +++ b/matlab/perfect-foresight-models/stacked_jacobian_substitute_controlled_paths.m @@ -0,0 +1,71 @@ +function A = stacked_jacobian_substitute_controlled_paths(A, y, y0, yT, exogenousvariables, steadystate, M_, options_) + +% Copyright © 2025 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/>. + +[periods, first_simulation_period] = get_simulation_periods(options_); + +% Reorganize M_.perfect_foresight_controlled_paths period by period +flips_by_period = struct('exogenize_id', cell(periods, 1), 'endogenize_id', cell(periods, 1)); +for i=1:length(M_.perfect_foresight_controlled_paths) + prange = M_.perfect_foresight_controlled_paths(i).periods; + if isa(prange, 'dates') + prange = transpose(prange - first_simulation_period + 1); + end + for p = prange + flips_by_period(p).exogenize_id(end+1) = M_.perfect_foresight_controlled_paths(i).exogenize_id; + flips_by_period(p).endogenize_id(end+1) = M_.perfect_foresight_controlled_paths(i).endogenize_id; + end +end + +% TODO: move this check to perfect_foresight_solver_core and add check with shocks block +for p = 1:periods + exogenize_id = flips_by_period(p).exogenize_id; + [~, idx] = unique(exogenize_id, 'stable'); + duplicate_idx = setdiff(1:numel(exogenize_id), idx); + if ~isempty(duplicate_idx) + error('perfect_foresight_controlled_paths: variable %s is exogenized two times in period %d', M_.endo_names{exogenize_id(duplicate_idx(1))}, p); + end + + endogenize_id = flips_by_period(p).endogenize_id; + [~, idx] = unique(endogenize_id, 'stable'); + duplicate_idx = setdiff(1:numel(endogenize_id), idx); + if ~isempty(duplicate_idx) + error('perfect_foresight_controlled_paths: variable %s is endogenized two times in period %d', M_.exo_names{endogenize_id(duplicate_idx(1))}, p); + end +end + +dynamic_g1 = str2func([M_.fname '.sparse.dynamic_g1']); + +for p = 1:periods + nsubst = length(flips_by_period(p).exogenize_id); + if nsubst > 0 + if p == 1 + y3n = [ y0; y(1:2*M_.endo_nbr) ]; + elseif p == periods + y3n = [ y(end-2*M_.endo_nbr+1:end); yT ]; + else + y3n = y((p-2)*M_.endo_nbr+(1:3*M_.endo_nbr)); + end + g1 = dynamic_g1(y3n, exogenousvariables(p+M_.maximum_lag,:), M_.params, steadystate, M_.dynamic_g1_sparse_rowval, M_.dynamic_g1_sparse_colval, M_.dynamic_g1_sparse_colptr); + + subst = [ zeros((p-1)*M_.endo_nbr, nsubst); + g1(:,3*M_.endo_nbr+flips_by_period(p).endogenize_id); + zeros((periods-p)*M_.endo_nbr, nsubst) ]; + A(:,(p-1)*M_.endo_nbr+flips_by_period(p).exogenize_id) = subst; + end +end diff --git a/scripts/dynare.el b/scripts/dynare.el index e723bf6b78b673ffdc20983b6bdca6e7c9f620c8..50d6e1d252274c9895ea43418bdda4a2b3dabd47 100644 --- a/scripts/dynare.el +++ b/scripts/dynare.el @@ -1,7 +1,7 @@ ;;; dynare.el --- major mode for editing Dynare mod files ;; Copyright © 2010 Yannick Kalantzis -;; Copyright © 2019-2024 Dynare Team +;; Copyright © 2019-2025 Dynare Team ;; ;; This program is free software; you can redistribute it and/or modify ;; it under the terms of the GNU General Public License as published by @@ -88,7 +88,7 @@ '("stderr" "values" "periods" "scales" "restriction" "exclusion" "upper_cholesky" "lower_cholesky" "equation" "bind" "relax" "error_bind" "error_relax" "add" "multiply" "target" "auxname_target_nonstationary" - "component" "growth" "auxname" "kind" "weights") + "component" "growth" "auxname" "kind" "weights" "exogenize" "endogenize") "Dynare statements-like keywords.") ;; Those keywords that makes the lexer enter the DYNARE_BLOCK start condition @@ -106,7 +106,7 @@ "conditional_forecast_paths" "svar_identification" "moment_calibration" "irf_calibration" "ramsey_constraints" "generate_irfs" "matched_moments" "occbin_constraints" "model_replace" "pac_target_info" "matched_irfs" - "matched_irfs_weights" "verbatim") + "matched_irfs_weights" "perfect_foresight_controlled_paths" "verbatim") "Dynare block keywords.")) ;; Mathematical functions and operators used in model equations (see "hand_side" in Bison file) diff --git a/tests/deterministic_simulations/ramst_pf_controlled_paths.mod b/tests/deterministic_simulations/ramst_pf_controlled_paths.mod new file mode 100644 index 0000000000000000000000000000000000000000..e3bbdc6607d35d22e1545a7738e1144e9c228930 --- /dev/null +++ b/tests/deterministic_simulations/ramst_pf_controlled_paths.mod @@ -0,0 +1,58 @@ +var c k y; +varexo x e_y; + +parameters alph gam delt bet aa; +alph=0.5; +gam=0.5; +delt=0.02; +bet=0.05; +aa=0.5; + + +model; +c + k - aa*x*k(-1)^alph - (1-delt)*k(-1); +c^(-gam) - (1+bet)^(-1)*(aa*alph*x(+1)*k^(alph-1) + 1 - delt)*c(+1)^(-gam); +y = 0.98*y(-1)+e_y; +end; + +initval; +x = 1; +k = ((delt+bet)/(1.0*aa*alph))^(1/(alph-1)); +c = aa*k^alph-delt*k; +y = 0; +end; + +steady; + +check; + +perfect_foresight_controlled_paths; + exogenize c; + periods 2 4:5; + values 1.6 1.7; + endogenize x; + + exogenize k; + periods 7:9; + values 13; + endogenize x; + + exogenize y; + periods 2003Y:2005Y; + values 0.5; + endogenize e_y; +end; + +shocks; +var x; +periods 1; +values 1.2; +end; + + + +perfect_foresight_setup(periods=30, first_simulation_period = 2001Y); +perfect_foresight_solver; + +rplot c; +rplot k;