solve_perfect_foresight_model

0001 function [flag,endo_simul,err] = solve_perfect_foresight_model(endo_simul,exo_simul,pfm)
0002 
0003     flag = 0;
0004     err = 0;
0005     stop = 0;
0006     nan_flag = 0;
0007 
0008     model_dynamic = pfm.dynamic_model;
0009 
0010     Y = endo_simul(:);
0011 
0012     if pfm.verbose
0013         disp (['-----------------------------------------------------']) ;
0014         disp (['MODEL SIMULATION :']) ;
0015         fprintf('\n') ;
0016     end
0017 
0018     z = Y(find(pfm.lead_lag_incidence'));
0019     [d1,jacobian] = model_dynamic(z,exo_simul,pfm.params,pfm.steady_state,2);
0020 
0021     A = sparse([],[],[],pfm.periods*pfm.ny,pfm.periods*pfm.ny,pfm.periods*nnz(jacobian));
0022     res = zeros(pfm.periods*pfm.ny,1);
0023 
0024     h1 = clock;
0025     for iter = 1:pfm.maxit_
0026         h2 = clock;
0027         i_rows = 1:pfm.ny;
0028         i_cols = find(pfm.lead_lag_incidence');
0029         i_cols_A = i_cols;
0030         for it = 2:(pfm.periods+1)
0031             [d1,jacobian] = model_dynamic(Y(i_cols),exo_simul,pfm.params,pfm.steady_state,it);
0032             if it == 2
0033                 A(i_rows,pfm.i_cols_A1) = jacobian(:,pfm.i_cols_1);
0034             elseif it == pfm.periods+1
0035                 A(i_rows,i_cols_A(pfm.i_cols_T)) = jacobian(:,pfm.i_cols_T);
0036             else
0037                 A(i_rows,i_cols_A) = jacobian(:,pfm.i_cols_j);
0038             end
0039             res(i_rows) = d1;
0040             i_rows = i_rows + pfm.ny;
0041             i_cols = i_cols + pfm.ny;
0042             if it > 2
0043                 i_cols_A = i_cols_A + pfm.ny;
0044             end
0045         end
0046         err = max(abs(res));
0047         if err < pfm.tolerance
0048             stop = 1 ;
0049             if pfm.verbose
0050                 fprintf('\n') ;
0051                 disp([' Total time of simulation        :' num2str(etime(clock,h1))]) ;
0052                 fprintf('\n') ;
0053                 disp([' Convergency obtained.']) ;
0054                 fprintf('\n') ;
0055             end
0056             flag = 0;% Convergency obtained.
0057             endo_simul = reshape(Y,pfm.ny,pfm.periods+2);
0058             break
0059         end
0060         dy = -A\res;
0061         if any(isnan(dy))
0062             nan_flag = 1;
0063             break
0064         end
0065         Y(pfm.i_upd) =   Y(pfm.i_upd) + dy;
0066     end
0067 
0068     if ~stop
0069         if pfm.verbose
0070             fprintf('\n') ;
0071             disp(['     Total time of simulation        :' num2str(etime(clock,h1))]) ;
0072             fprintf('\n') ;
0073             disp(['WARNING : maximum number of iterations is reached (modify options_.maxit_).']) ;
0074             fprintf('\n') ;
0075         end
0076         flag = 1;% more iterations are needed.
0077         endo_simul = 1;
0078     end
0079     if nan_flag
0080         if pfm.verbose
0081             fprintf('\n') ;
0082             disp(['     Total time of simulation        :' num2str(etime(clock,h1))]) ;
0083             fprintf('\n') ;
0084             disp(['WARNING : NaNs!']) ;
0085             fprintf('\n') ;
0086         end
0087         flag = 1;
0088         endo_simul = 1;
0089     end
0090     if pfm.verbose
0091         disp (['-----------------------------------------------------']) ;
0092     end