From 5d43070041a42def472e688c8b6c192de4e79de4 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?St=C3=A9phane=20Adjemian=20=28Argos=29?= <stepan@adjemian.eu> Date: Fri, 14 Mar 2025 11:45:54 +0100 Subject: [PATCH] Remove unused routine. --- matlab/ep/extended_path_homotopy.m | 4 +- ...solve_stochastic_perfect_foresight_model.m | 254 ------------------ 2 files changed, 2 insertions(+), 256 deletions(-) delete mode 100644 matlab/ep/solve_stochastic_perfect_foresight_model.m diff --git a/matlab/ep/extended_path_homotopy.m b/matlab/ep/extended_path_homotopy.m index ce5c207a47..c4df1c1592 100644 --- a/matlab/ep/extended_path_homotopy.m +++ b/matlab/ep/extended_path_homotopy.m @@ -40,7 +40,7 @@ if ismember(method, [1, 2]) switch(algo) case 0 [flag, endo_simul_new] = ... - solve_stochastic_perfect_foresight_model(endo_simul, exo_simul, pfm, ep.stochastic.quadrature.nodes, ep.stochastic.order); + solve_stochastic_perfect_foresight_model_0(endo_simul, exo_simul, pfm, ep.stochastic.quadrature.nodes, ep.stochastic.order); case 1 [flag, endo_simul_new] = ... solve_stochastic_perfect_foresight_model_1(endo_simul, exo_simul, options_, pfm, ep.stochastic.order); @@ -107,7 +107,7 @@ if isequal(method, 3) || (isequal(method, 2) && noconvergence) switch(algo) case 0 [flag, endo_simul_new] = ... - solve_stochastic_perfect_foresight_model(endo_simul, exo_simul, pfm, ep.stochastic.quadrature.nodes, ep.stochastic.order); + solve_stochastic_perfect_foresight_model_0(endo_simul, exo_simul, pfm, ep.stochastic.quadrature.nodes, ep.stochastic.order); case 1 [flag, endo_simul_new] = ... solve_stochastic_perfect_foresight_model_1(endo_simul, exo_simul, options_, pfm, ep.stochastic.order); diff --git a/matlab/ep/solve_stochastic_perfect_foresight_model.m b/matlab/ep/solve_stochastic_perfect_foresight_model.m deleted file mode 100644 index 8263f5e706..0000000000 --- a/matlab/ep/solve_stochastic_perfect_foresight_model.m +++ /dev/null @@ -1,254 +0,0 @@ -function [flag,endo_simul,err] = solve_stochastic_perfect_foresight_model(endo_simul,exo_simul,pfm,nnodes,order) - -% Copyright © 2012-2017 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/>. - -flag = 0; -err = 0; -stop = 0; - -params = pfm.params; -steady_state = pfm.steady_state; -ny = pfm.ny; -periods = pfm.periods; -dynamic_model = pfm.dynamic_model; -lead_lag_incidence = pfm.lead_lag_incidence; -lead_lag_incidence_t = transpose(lead_lag_incidence); -nyp = pfm.nyp; -nyf = pfm.nyf; -i_cols_1 = pfm.i_cols_1; -i_cols_j = pfm.i_cols_j; -i_cols_T = nonzeros(lead_lag_incidence(1:2,:)'); -hybrid_order = pfm.hybrid_order; -if hybrid_order > 0 - if hybrid_order == 2 - h_correction = 0.5*pfm.dr.ghs2(pfm.dr.inv_order_var); - end -else - h_correction = 0; -end -%h_correction = pfm.h_correction; - - -maxit = pfm.maxit_; -tolerance = pfm.tolerance; -verbose = pfm.verbose; - -number_of_shocks = size(exo_simul,2); - -[nodes,weights] = gauss_hermite_weights_and_nodes(nnodes); - -if number_of_shocks>1 - nodes = repmat(nodes,1,number_of_shocks)*chol(pfm.Sigma); - % to be fixed for Sigma ~= I - for i=number_of_shocks:-1:1 - rr(i) = {nodes(:,i)}; - ww(i) = {weights}; - end - nodes = cartesian_product_of_sets(rr{:}); - weights = prod(cartesian_product_of_sets(ww{:}),2); - nnodes = nnodes^number_of_shocks; -else - nodes = nodes*sqrt(pfm.Sigma); -end - -if verbose - disp (' -----------------------------------------------------'); - disp ('MODEL SIMULATION :'); - fprintf('\n'); -end - -z = endo_simul(lead_lag_incidence_t(:)>0); -[~, jacobian] = dynamic_model(z, exo_simul, params,steady_state, 2); - -% Each column of Y represents a different world -% The upper right cells are unused -% The first row block is ny x 1 -% The second row block is ny x nnodes -% The third row block is ny x nnodes^2 -% and so on until size ny x nnodes^order -world_nbr = nnodes^order; -Y = repmat(endo_simul(:),1,world_nbr); - -% The columns of A map the elements of Y such that -% each block of Y with ny rows are unfolded column wise -dimension = ny*(sum(nnodes.^(0:order-1),2)+(periods-order)*world_nbr); -if order == 0 - i_upd_r = (1:ny*periods); - i_upd_y = i_upd_r + ny; -else - i_upd_r = zeros(dimension,1); - i_upd_y = i_upd_r; - i_upd_r(1:ny) = (1:ny); - i_upd_y(1:ny) = ny+(1:ny); - i1 = ny+1; - i2 = 2*ny; - n1 = ny+1; - n2 = 2*ny; - for i=2:periods - for j=1:nnodes^min(i-1,order) - i_upd_r(i1:i2) = (n1:n2)+(j-1)*ny*periods; - i_upd_y(i1:i2) = (n1:n2)+ny+(j-1)*ny*(periods+2); - i1 = i2+1; - i2 = i2+ny; - end - n1 = n2+1; - n2 = n2+ny; - end -end -if rows(lead_lag_incidence)>2 - icA = [find(lead_lag_incidence(1,:)) find(lead_lag_incidence(2,:))+world_nbr*ny ... - find(lead_lag_incidence(3,:))+2*world_nbr*ny]'; -else - if nyf - icA = [find(lead_lag_incidence(2,:))+world_nbr*ny find(lead_lag_incidence(3,:))+2*world_nbr*ny ]'; - else - icA = [find(lead_lag_incidence(1,:)) find(lead_lag_incidence(2,:))+world_nbr*ny ]'; - end -end -h1 = clock; -for iter = 1:maxit - A1 = sparse([],[],[],ny*(sum(nnodes.^(0:order-1),2)+1),dimension,(order+1)*world_nbr*nnz(jacobian)); - res = zeros(ny,periods,world_nbr); - i_rows = 1:ny; - i_cols = find(lead_lag_incidence'); - i_cols_p = i_cols(1:nyp); - i_cols_s = i_cols(nyp+(1:ny)); - i_cols_f = i_cols(nyp+ny+(1:nyf)); - i_cols_Ap = i_cols_p; - i_cols_As = i_cols_s; - i_cols_Af = i_cols_f - ny; - i_hc = i_cols_f - 2*ny; - for i = 1:order+1 - i_w_p = 1; - for j = 1:nnodes^(i-1) - innovation = exo_simul; - if i > 1 - innovation(i+1,:) = nodes(mod(j-1,nnodes)+1,:); - end - if i <= order - for k=1:nnodes - y = [Y(i_cols_p,i_w_p); - Y(i_cols_s,j); - Y(i_cols_f,(j-1)*nnodes+k)]; - [d1,jacobian] = dynamic_model(y,innovation,params,steady_state,i+1); - if i == 1 - % in first period we don't keep track of - % predetermined variables - i_cols_A = [i_cols_As - ny; i_cols_Af]; - A1(i_rows,i_cols_A) = A1(i_rows,i_cols_A) + weights(k)*jacobian(:,i_cols_1); - else - i_cols_A = [i_cols_Ap; i_cols_As; i_cols_Af]; - A1(i_rows,i_cols_A) = A1(i_rows,i_cols_A) + weights(k)*jacobian(:,i_cols_j); - end - res(:,i,j) = res(:,i,j)+weights(k)*d1; - i_cols_Af = i_cols_Af + ny; - end - else - % i==order+1 - if hybrid_order==2 - y = [Y(i_cols_p,i_w_p); - Y(i_cols_s,j); - Y(i_cols_f,j)+h_correction(i_hc)]; - else - y = [Y(i_cols_p,i_w_p); - Y(i_cols_s,j); - Y(i_cols_f,j)]; - end - [d1,jacobian] = dynamic_model(y,innovation,params,steady_state,i+1); - if i == 1 - % in first period we don't keep track of - % predetermined variables - i_cols_A = [i_cols_As - ny; i_cols_Af]; - A1(i_rows,i_cols_A) = jacobian(:,i_cols_1); - else - i_cols_A = [i_cols_Ap; i_cols_As; i_cols_Af]; - A1(i_rows,i_cols_A) = jacobian(:,i_cols_j); - end - res(:,i,j) = d1; - i_cols_Af = i_cols_Af + ny; - end - i_rows = i_rows + ny; - if mod(j,nnodes) == 0 - i_w_p = i_w_p + 1; - end - if i > 1 - if mod(j,nnodes) == 0 - i_cols_Ap = i_cols_Ap + ny; - end - i_cols_As = i_cols_As + ny; - end - end - i_cols_p = i_cols_p + ny; - i_cols_s = i_cols_s + ny; - i_cols_f = i_cols_f + ny; - end - nzA = cell(periods,world_nbr); - for j=1:world_nbr - i_rows_y = find(lead_lag_incidence')+(order+1)*ny; - offset_c = ny*(sum(nnodes.^(0:order-1),2)+j-1); - offset_r = (j-1)*ny; - for i=order+2:periods - [d1,jacobian] = dynamic_model(Y(i_rows_y,j), ... - exo_simul,params, ... - steady_state,i+1); - if i == periods - [ir,ic,v] = find(jacobian(:,i_cols_T)); - else - [ir,ic,v] = find(jacobian(:,i_cols_j)); - end - nzA{i,j} = [offset_r+ir,offset_c+icA(ic), v]'; - res(:,i,j) = d1; - i_rows_y = i_rows_y + ny; - offset_c = offset_c + world_nbr*ny; - offset_r = offset_r + world_nbr*ny; - end - end - err = max(abs(res(i_upd_r))); - if err < tolerance - stop = 1; - if verbose - fprintf('\n') ; - disp([' Total time of simulation :' num2str(etime(clock,h1))]) ; - fprintf('\n') ; - disp(' Convergency obtained.') ; - fprintf('\n') ; - end - flag = 0;% Convergency obtained. - endo_simul = reshape(Y(:,1),ny,periods+2); - break - end - A2 = [nzA{:}]'; - A = [A1; sparse(A2(:,1),A2(:,2),A2(:,3),ny*(periods-order-1)*world_nbr,dimension)]; - dy = -A\res(i_upd_r); - Y(i_upd_y) = Y(i_upd_y) + dy; -end - -if ~stop - if verbose - fprintf('\n') ; - disp([' Total time of simulation :' num2str(etime(clock,h1))]) ; - fprintf('\n') ; - disp('WARNING : maximum number of iterations is reached (modify options_.simul.maxit).') ; - fprintf('\n') ; - end - flag = 1;% more iterations are needed. - endo_simul = 1; -end -if verbose - disp ('-----------------------------------------------------') ; -end -- GitLab