Make perfect_foresight_controlled_paths compatible with the “linear” option

doc/manual/source/the-model-file.rst

@@ -4672,8 +4672,8 @@ that have been left free.
     The ``perfect_foresight_controlled_paths`` block requires that the
     :opt:`stack_solve_algo <stack_solve_algo = INTEGER>` option be equal to
     either ``0``, ``1``, ``2``, ``3``, ``6`` or ``7``, and is incompatible with
-    the :opt:`block`, :opt:`linear` and :opt:`bytecode` options of the
-    :bck:`model` block and :comm:`model_options` command.
+    the :opt:`block` and :opt:`bytecode` options of the :bck:`model` block and
+    :comm:`model_options` command.
     *Options*

matlab/perfect-foresight-models/perfect_foresight_solver_core.m

@@ -90,9 +90,6 @@ if ~isempty(controlled_paths_by_period)
     if options_.block
         error('perfect_foresight_controlled_paths is not available with the block option')
     end
-    if options_.linear
-        error('perfect_foresight_controlled_paths is not available with the linear option')
-    end
     if M_.maximum_endo_lead == 0 || M_.maximum_endo_lag == 0
         error('perfect_foresight_controlled_paths is not available with purely backward, purely forward or static models')
     end
@@ -147,7 +144,7 @@ else
             switch options_.stack_solve_algo
               case {0 2 3}
                 if options_.linear_approximation
-                    [y, success, maxerror] = sim1_linear(y, exo_simul, steady_state, exo_steady_state, M_, options_);
+                    [y, success, maxerror, exo_simul] = sim1_linear(y, exo_simul, steady_state, exo_steady_state, controlled_paths_by_period, M_, options_);
                 else
                     [y, success, maxerror, iter, exo_simul] = sim1(y, exo_simul, steady_state, controlled_paths_by_period, M_, options_);
                 end

matlab/perfect-foresight-models/sim1_linear.m

-function [endogenousvariables, success, ERR] = sim1_linear(endogenousvariables, exogenousvariables, steadystate_y, steadystate_x, M_, options_)
+function [endogenousvariables, success, ERR, exogenousvariables] = sim1_linear(endogenousvariables, exogenousvariables, steadystate_y, steadystate_x, controlled_paths_by_period, M_, options_)
 % [endogenousvariables, success, ERR] = sim1_linear(endogenousvariables, exogenousvariables, steadystate_y, steadystate_x, M_, options_)
 % Solves a linear approximation of a perfect foresight model using sparse matrix.
 %
 % INPUTS
-% - endogenousvariables [double] N*T array, paths for the endogenous variables (initial guess).
+% - endogenousvariables [double] N*T array, paths for the endogenous variables (initial condition + initial guess + terminal condition).
 % - exogenousvariables  [double] T*M array, paths for the exogenous variables.
 % - steadystate_y       [double] N*1 array, steady state for the endogenous variables.
 % - steadystate_x       [double] M*1 array, steady state for the exogenous variables.
+% - controlled_paths_by_period [struct] data from perfect_foresight_controlled_paths block
 % - M_                  [struct] contains a description of the model.
 % - options_            [struct] contains various options.
 %
@@ -14,6 +15,8 @@ function [endogenousvariables, success, ERR] = sim1_linear(endogenousvariables,
 % - endogenousvariables [double] N*T array, paths for the endogenous variables (solution of the perfect foresight model).
 % - success             [logical] Whether a solution was found
 % - ERR                 [double] ∞-norm of the residual
+% - exogenousvariables  [double] T*M array, paths for the exogenous variables
+%                                (may be modified if perfect_foresight_controlled_paths present)
 %
 % NOTATIONS
 % - N is the number of endogenous variables.
@@ -39,7 +42,7 @@ function [endogenousvariables, success, ERR] = sim1_linear(endogenousvariables,
 % to center the variables around the deterministic steady state to solve the
 % perfect foresight model.
 
-% Copyright © 2015-2024 Dynare Team
+% Copyright © 2015-2025 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -67,6 +70,19 @@ periods = get_simulation_periods(options_);
 
 params = M_.params;
 
+if ~isempty(controlled_paths_by_period)
+    for p = 1:periods
+        if isempty(controlled_paths_by_period(p).exogenize_id)
+            continue
+        end
+        endogenousvariables(controlled_paths_by_period(p).exogenize_id, p+M_.maximum_lag) = controlled_paths_by_period(p).values;
+    end
+
+    if options_.debug
+        error('Debugging not available with perfect_foresight_controlled_paths')
+    end
+end
+
 if verbose
     skipline()
     printline(80)
@@ -125,6 +141,10 @@ h2 = clock;
 
 [res, A] = linear_perfect_foresight_problem(Y, jacobian, y0, yT, exogenousvariables, params, steadystate_y, maximum_lag, periods, ny);
 
+if ~isempty(controlled_paths_by_period)
+    A = controlled_paths_substitute_stacked_jacobian(A, Y, y0, yT, exogenousvariables, steadystate_y, controlled_paths_by_period, M_);
+end
+
 % Evaluation of the maximum residual at the initial guess (steady state for the endogenous variables).
 err = max(abs(res));
 
@@ -141,7 +161,8 @@ end
 
 % Try to update the vector of endogenous variables.
 try
-    Y =  Y - A\res;
+    dY = -A\res;
+    Y = Y + dY;
 catch
     % Normally, because the model is linear, the solution of the perfect foresight model should
     % be obtained in one Newton step. This is not the case if the model is singular.
@@ -154,6 +175,18 @@ catch
     return
 end
 
+if ~isempty(controlled_paths_by_period)
+    for p = 1:periods
+        endogenize_id = controlled_paths_by_period(p).endogenize_id;
+        exogenize_id = controlled_paths_by_period(p).exogenize_id;
+        if isempty(endogenize_id)
+            continue
+        end
+        Y(exogenize_id+(p-1)*M_.endo_nbr) = controlled_paths_by_period(p).values;
+        exogenousvariables(p+M_.maximum_lag,endogenize_id) = exogenousvariables(p+M_.maximum_lag,endogenize_id) + dY(exogenize_id+(p-1)*M_.endo_nbr)';
+    end
+end
+
 YY = reshape([y0; Y; yT], ny, periods+2);
 
 i_rows = (1:ny)';

meson.build

@@ -1385,6 +1385,7 @@ mod_and_m_tests = [
   { 'test' : [ 'deterministic_simulations/Irreversible_investment.mod',
                'deterministic_simulations/Irreversible_investment_leadlag.mod' ] },
   { 'test' : [ 'deterministic_simulations/linear_state_space_arma.mod' ] },
+  { 'test' : [ 'deterministic_simulations/linear_state_space_arma_condpaths.mod' ] },
   { 'test' : [ 'deterministic_simulations/purely_forward/ar1.mod' ] },
   { 'test' : [ 'deterministic_simulations/purely_forward/nk.mod' ] },
   { 'test' : [ 'deterministic_simulations/purely_backward/ar1.mod' ] },

tests/deterministic_simulations/linear_state_space_arma_condpaths.mod

+% Tests perfect simulation of linear model (sim1_linear.m) with controlled paths
+
+var x
+    y
+    z;
+
+varexo u
+       v;
+
+parameters a1 a2 a3 a4
+	   b1 b2 b3
+	   c1;
+
+a1 =  .50;
+a2 =  .00;
+a3 =  .70;
+a4 =  .40;
+b1 =  .90;
+b2 =  .00;
+b3 =  .80;
+c1 =  .95;
+
+
+model(linear);
+  y = a1*x(-1) + a2*x(+1) + a3*z + a4*y(-1);
+  z = b1*z(-1) + b2*z(+1) + b3*x + u;
+  x = c1*x(-1) + v +v(-1)+v(+1);
+end;
+
+initval;
+  y=-1;
+  x=-1;
+  z=-1;
+end;
+
+endval;
+  y=0;
+  x=0;
+  z=0;
+end;
+
+perfect_foresight_controlled_paths;
+  exogenize z;
+  periods 50;
+  values 0.5;
+  endogenize u;
+end;
+
+steady;
+
+perfect_foresight_setup(periods=100);
+perfect_foresight_solver;
+
+if ~oo_.deterministic_simulation.status
+   error('Perfect foresight simulation failed')
+end
+
+if ~(oo_.endo_simul(3, 51) == 0.5 && isequal(find(oo_.exo_simul ~= 0), 51))
+   error('Linear perfect foresight simulation with controlled paths failed')
+end