Commit e1f17fa9 authored by ferhat's avatar ferhat
Browse files

- adds a new option in simul command when block and bytecode options are used...

- adds a new option in simul command when block and bytecode options are used : minimal_solving_periods. With bytecode option, the linear solver is applied only on the first periods. As soon as the set of operations remains the same from one period to another, they are repeated for the remaining periods. In some cases the linear solver could encounter exact or near singularities. To avoid these effects, this new option allows the user to extend the minimal number of periods where the model has to be solved.
- corrects some mod files.
- In steady command, get check argument from bytecode dll

git-svn-id: https://www.dynare.org/svn/dynare/trunk@3051 ac1d8469-bf42-47a9-8791-bf33cf982152
parent d87129eb
......@@ -108,6 +108,8 @@ function global_initialization()
% Deterministic simulation
options_.stack_solve_algo = 0;
options_.markowitz = 0.5;
options_.minimal_solving_periods = 1;
% Solution
options_.order = 2;
......@@ -153,7 +155,6 @@ function global_initialization()
options_.load_mh_file = 0;
options_.logdata = 0;
options_.loglinear = 0;
options_.markowitz = 0.5;
options_.mh_conf_sig = 0.90;
options_.prior_interval = 0.90;
options_.mh_drop = 0.5;
......
......@@ -75,8 +75,7 @@ function steady_()
check1 = 1;
end
elseif options_.block && options_.bytecode
residuals = bytecode('evaluate','static');
check1 = 1;
[residuals, check1] = bytecode('evaluate','static');
else
check1 = 0;
check1 = max(abs(feval([M_.fname '_static'],...
......@@ -86,9 +85,7 @@ function steady_()
end
end
if check1
if ~options_.block && ~options_.bytecode
resid(1);
end
resid(1);
error(['The seadystate values returned by ' M_.fname ...
'_steadystate.m don''t solve the static model!' ])
end
......@@ -133,4 +130,4 @@ function steady_()
if check ~= 0
error('STEADY: convergence problems')
end
\ No newline at end of file
end
......@@ -28,7 +28,7 @@ Interpreter::Interpreter(double *params_arg, double *y_arg, double *ya_arg, doub
double *direction_arg, int y_size_arg,
int nb_row_x_arg, int nb_row_xd_arg, int periods_arg, int y_kmin_arg, int y_kmax_arg,
int maxit_arg_, double solve_tolf_arg, int size_of_direction_arg, double slowc_arg, int y_decal_arg, double markowitz_c_arg,
string &filename_arg)
string &filename_arg, int minimal_solving_periods_arg)
{
params=params_arg;
y=y_arg;
......@@ -53,6 +53,7 @@ Interpreter::Interpreter(double *params_arg, double *y_arg, double *ya_arg, doub
filename=filename_arg;
T=NULL;
error_not_printed = true;
minimal_solving_periods = minimal_solving_periods_arg;
}
double
......@@ -1359,10 +1360,10 @@ Interpreter::simulate_a_block(const int size, const int type, string file_name,
cvg = false;
else
cvg=(max_res<solve_tolf);
if(cvg)
continue;
/*if(cvg)
continue;*/
u_count=u_count_saved;
simulate_NG1(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, periods, true, cvg, iter);
simulate_NG1(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, periods, true, cvg, iter, minimal_solving_periods);
iter++;
}
if (!cvg)
......@@ -1395,7 +1396,7 @@ Interpreter::simulate_a_block(const int size, const int type, string file_name,
}
}
cvg = false;
simulate_NG1(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, periods, true, cvg, iter);
simulate_NG1(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, periods, true, cvg, iter, minimal_solving_periods);
}
mxFree(r);
mxFree(y_save);
......
......@@ -70,12 +70,13 @@ class Interpreter : SparseMatrix
map<pair<pair<int, int> ,int>, int> IM_i;
int equation, derivative_equation, derivative_variable;
string filename;
int minimal_solving_periods;
public :
Interpreter(double *params_arg, double *y_arg, double *ya_arg, double *x_arg, double *steady_y_arg, double *steady_x_arg,
double *direction_arg, int y_size_arg, int nb_row_x_arg,
int nb_row_xd_arg, int periods_arg, int y_kmin_arg, int y_kmax_arg, int maxit_arg_, double solve_tolf_arg, int size_o_direction_arg,
double slowc_arg, int y_decal_arg, double markowitz_c_arg, string &filename_arg);
double slowc_arg, int y_decal_arg, double markowitz_c_arg, string &filename_arg, int minimal_solving_periods_arg);
bool compute_blocks(string file_name, string bin_basename, bool steady_state, bool evaluate);
};
......
......@@ -1546,7 +1546,7 @@ SparseMatrix::Check_the_Solution(int periods, int y_kmin, int y_kmax, int Size,
}
int
SparseMatrix::simulate_NG1(int blck, int y_size, int it_, int y_kmin, int y_kmax, int Size, int periods, bool print_it, bool cvg, int &iter)
SparseMatrix::simulate_NG1(int blck, int y_size, int it_, int y_kmin, int y_kmax, int Size, int periods, bool print_it, bool cvg, int &iter, int minimal_solving_periods)
{
/*Triangularisation at each period of a block using a simple gaussian Elimination*/
t_save_op_s *save_op_s;
......@@ -1634,7 +1634,7 @@ SparseMatrix::simulate_NG1(int blck, int y_size, int it_, int y_kmin, int y_kmax
}
else
{
start_compare = y_kmin;
start_compare = max( y_kmin, minimal_solving_periods);
restart = 0;
}
res1a = res1;
......
......@@ -60,7 +60,7 @@ class SparseMatrix
{
public:
SparseMatrix();
int simulate_NG1(int blck, int y_size, int it_, int y_kmin, int y_kmax, int Size, int periods, bool print_it, bool cvg, int &iter);
int simulate_NG1(int blck, int y_size, int it_, int y_kmin, int y_kmax, int Size, int periods, bool print_it, bool cvg, int &iter, int minimal_solving_periods);
bool simulate_NG(int blck, int y_size, int it_, int y_kmin, int y_kmax, int Size, bool print_it, bool cvg, int &iter, bool steady_state);
void Direct_Simulate(int blck, int y_size, int it_, int y_kmin, int y_kmax, int Size, int periods, bool print_it, int iter);
void fixe_u(double **u, int u_count_int, int max_lag_plus_max_lead_plus_1);
......
......@@ -291,6 +291,7 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
int maxit_=int(floor(*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,"maxit_"))))));
double slowc=double(*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,"slowc")))));
double markowitz_c=double(*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,"markowitz")))));
int minimal_solving_periods=int(*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,"minimal_solving_periods")))));
double solve_tolf;
if(steady_state)
solve_tolf=*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,"solve_tolf"))));
......@@ -326,7 +327,7 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
for (int j = 0; j < y_size; j++)
mexPrintf(" variable %d at time %d and %d = %f\n", j+1, it_, it_+1, y[j+it_*y_size]);*/
clock_t t0= clock();
Interpreter interprete(params, y, ya, x, steady_yd, steady_xd, direction, y_size, nb_row_x, nb_row_xd, periods, y_kmin, y_kmax, maxit_, solve_tolf, size_of_direction, slowc, y_decal, markowitz_c, file_name);
Interpreter interprete(params, y, ya, x, steady_yd, steady_xd, direction, y_size, nb_row_x, nb_row_xd, periods, y_kmin, y_kmax, maxit_, solve_tolf, size_of_direction, slowc, y_decal, markowitz_c, file_name, minimal_solving_periods);
string f(fname);
bool result = interprete.compute_blocks(f, f, steady_state, evaluate);
clock_t t1= clock();
......
......@@ -106,7 +106,7 @@ class ParsingDriver;
%token KALMAN_ALGO KALMAN_TOL
%token LABELS LAPLACE LIK_ALGO LIK_INIT LINEAR LOAD_IDENT_FILES LOAD_MH_FILE LOAD_PARAMS_AND_STEADY_STATE LOGLINEAR
%token MARKOWITZ MARGINAL_DENSITY MAX
%token MFS MH_DROP MH_INIT_SCALE MH_JSCALE MH_MODE MH_NBLOCKS MH_REPLIC MH_RECOVER MIN
%token MFS MH_DROP MH_INIT_SCALE MH_JSCALE MH_MODE MH_NBLOCKS MH_REPLIC MH_RECOVER MIN MINIMAL_SOLVING_PERIODS
%token MODE_CHECK MODE_COMPUTE MODE_FILE MODEL MODEL_COMPARISON MODEL_INFO MSHOCKS
%token MODIFIEDHARMONICMEAN MOMENTS_VARENDO DIFFUSE_FILTER
%token <string_val> NAME
......@@ -141,16 +141,16 @@ class ParsingDriver;
%token KSSTAT_REDFORM ALPHA2_REDFORM NAMENDO NAMLAGENDO NAMEXO RMSE LIK_ONLY VAR_RMSE PFILT_RMSE ISTART_RMSE
%token ALPHA_RMSE ALPHA2_RMSE TRANS_IDENT
/* end of GSA analysis*/
%token FREQ INITIAL_YEAR INITIAL_SUBPERIOD FINAL_YEAR FINAL_SUBPERIOD DATA VLIST VARLIST LOG_VAR PERCENT_VAR
%token FREQ INITIAL_YEAR INITIAL_SUBPERIOD FINAL_YEAR FINAL_SUBPERIOD DATA VLIST VARLIST LOG_VAR PERCENT_VAR
%token VLISTLOG VLISTPER
%token RESTRICTION_FNAME NLAGS CROSS_RESTRICTIONS CONTEMP_REDUCED_FORM REAL_PSEUDO_FORECAST BAYESIAN_PRIOR
%token DUMMY_OBS NSTATES INDXSCALESSTATES
%token <string_val> ALPHA BETA ABAND NINV CMS NCMS CNUM
%token GSIG2_LMD GSIG2_LMDM Q_DIAG FLAT_PRIOR NCSK NSTD
%token DUMMY_OBS NSTATES INDXSCALESSTATES
%token <string_val> ALPHA BETA ABAND NINV CMS NCMS CNUM
%token GSIG2_LMD GSIG2_LMDM Q_DIAG FLAT_PRIOR NCSK NSTD
%token INDXPARR INDXOVR INDXAP APBAND INDXIMF IMFBAND INDXFORE FOREBAND INDXGFOREHAT INDXGIMFHAT
%token INDXESTIMA INDXGDLS EQ_MS
%token INDXESTIMA INDXGDLS EQ_MS
%token EQ_CMS TLINDX TLNUMBER BANACT CREATE_INITIALIZATION_FILE ESTIMATE_MSMODEL
%token COMPUTE_MDD COMPUTE_PROBABILITIES PRINT_DRAWS N_DRAWS THINNING_FACTOR PROPOSAL_DRAWS MARKOV_FILE
%token COMPUTE_MDD COMPUTE_PROBABILITIES PRINT_DRAWS N_DRAWS THINNING_FACTOR PROPOSAL_DRAWS MARKOV_FILE
%token MHM_FILE OUTPUT_FILE_TAG DRAWS_NBR_BURN_IN_1 DRAWS_NBR_BURN_IN_2 DRAWS_NBR_MEAN_VAR_ESTIMATE
%token DRAWS_NBR_MODIFIED_HARMONIC_MEAN DIRICHLET_SCALE
%token SBVAR MS_SBVAR
......@@ -677,6 +677,7 @@ simul_options : o_periods
| o_datafile
| o_stack_solve_algo
| o_markowitz
| o_minimal_solving_periods
;
stoch_simul : STOCH_SIMUL ';'
......@@ -1317,7 +1318,7 @@ sbvar_option : o_datafile
| o_tlindx
| o_tlnumber
| o_cnum
| o_forecast
| o_forecast
;
sbvar_options_list : sbvar_option COMMA sbvar_options_list
......@@ -1389,12 +1390,12 @@ ms_sbvar_option : o_datafile
| o_thinning_factor
| o_markov_file
| o_mhm_file
| o_proposal_draws
| o_proposal_draws
| o_draws_nbr_burn_in_1
| o_draws_nbr_burn_in_2
| o_draws_nbr_mean_var_estimate
| o_draws_nbr_modified_harmonic_mean
| o_dirichlet_scale
| o_dirichlet_scale
;
ms_sbvar_options_list : ms_sbvar_option COMMA ms_sbvar_options_list
......@@ -1562,6 +1563,7 @@ o_periods : PERIODS EQUAL INT_NUMBER
{ driver.option_num("periods", $3); driver.option_num("simul", "1"); };
o_cutoff : CUTOFF EQUAL number { driver.cutoff($3); }
o_markowitz : MARKOWITZ EQUAL number { driver.option_num("markowitz", $3); };
o_minimal_solving_periods : MINIMAL_SOLVING_PERIODS EQUAL number { driver.option_num("minimal_solving_periods", $3); };
o_mfs : MFS EQUAL INT_NUMBER { driver.mfs($3); };
o_simul : SIMUL { driver.option_num("simul", "1"); };
o_simul_seed : SIMUL_SEED EQUAL INT_NUMBER { driver.option_num("simul_seed", $3); } ;
......@@ -1803,13 +1805,13 @@ vec_value_1 : '[' value1
vec_value : vec_value_1 ']' { $1->append("]"); $$ = $1; };
symbol : NAME
| ALPHA
| BETA
| NINV
| ABAND
| CMS
| NCMS
| CNUM
| ALPHA
| BETA
| NINV
| ABAND
| CMS
| NCMS
| CNUM
;
%%
......
......@@ -210,6 +210,7 @@ int sigma_e = 0;
<DYNARE_STATEMENT>nocorr {return token::NOCORR;}
<DYNARE_STATEMENT>optim {return token::OPTIM;}
<DYNARE_STATEMENT>periods {return token::PERIODS;}
<DYNARE_STATEMENT>minimal_solving_periods {return token::MINIMAL_SOLVING_PERIODS;}
<DYNARE_STATEMENT>markowitz {return token::MARKOWITZ;}
<DYNARE_STATEMENT>marginal_density {return token::MARGINAL_DENSITY;}
<DYNARE_STATEMENT>laplace {return token::LAPLACE;}
......
......@@ -17,8 +17,9 @@
// Note that there is an initial minus sign missing in equation (A1), p. S63.
//
// Michel Juillard, February 2004
var m P c e W R k d n l gy_obs gp_obs y dA vv ww;
@#define bytecode = 1
@#define block = 1
var m P c e W R k d n l gy_obs gp_obs y dA;
varexo e_a e_m;
parameters alp bet gam mst rho psi del;
......@@ -30,11 +31,16 @@ mst = 1.011;
rho = 0.7;
psi = 0.787;
del = 0.02;
toto = [2 3];
model(block, bytecode);
//model(sparse);
//model;
@#if block == 1
@#if bytecode == 1
model(block, bytecode);
@#else
model(block);
@#endif
@#else
model;
@#endif
/*0*/ exp(gam+e_a) = dA ;
/*1*/ log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m;
/*2*/ -P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0;
......@@ -49,19 +55,6 @@ model(block, bytecode);
/*11*/ k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a)) = y ;
/*12*/ gy_obs = dA*y/y(-1);
/*13*/ gp_obs = (P/P(-1))*m(-1)/dA;
/*14*/ vv = 0.2*ww+0.5*vv(-1)+1+c(-1)+e_a;
/*15*/ ww = 0.1*vv+0.5*ww(-1)+2;
/* A lt=
0.5*vv-0.2*ww = 1
-0.1*vv+0.5*ww = 2
[ 0.5 -0.2][vv] [1]
=
[-0.1 0.5][ww] [2]
det = 0.25-0.02 = 0.23
[vv] [0.5 0.2] [1] [0.9] [3.91304]
= 1/0.23* = 1/0.23* =
[ww] [0.1 0.5] [2] [1.1] [4.7826]
*/
end;
initval;
......@@ -81,8 +74,6 @@ gp_obs = exp(-gam);
dA = exp(gam);
e_a=0;
e_m=0;
vv = 0;
ww = 0;
end;
shocks;
......@@ -91,39 +82,37 @@ var e_m; stderr 0.005;
end;
//options_.solve_tolf=1e-10;
options_.maxit_=10;
steady;//(block_mfs);
@#if bytecode == 1
steady(solve_algo = 5);
@#else
steady(solve_algo = 4);
@#endif
model_info;
//check;
@#if block == 0
check;
@#endif
shocks;
var e_a;
periods 1;
values 0.16;
end;
disp(toto(1,2));
simul(periods=20, method=lu);
//stoch_simul(periods=200,order=1);
@#if block == 1
@#if bytecode == 1
simul(periods=200, stack_solve_algo = 5);
@#else
simul(periods=200, stack_solve_algo = 1);
@#endif
@#else
simul(periods=200, stack_solve_algo = 0);
@#endif
rplot y;
rplot k;
rplot c;
/*estimated_params;
alp, beta_pdf, 0.356, 0.02;
bet, beta_pdf, 0.993, 0.002;
gam, normal_pdf, 0.0085, 0.003;
mst, normal_pdf, 1.0002, 0.007;
rho, beta_pdf, 0.129, 0.223;
psi, beta_pdf, 0.65, 0.05;
del, beta_pdf, 0.01, 0.005;
stderr e_a, inv_gamma_pdf, 0.035449, inf;
stderr e_m, inv_gamma_pdf, 0.008862, inf;
end;
varobs gp_obs gy_obs;
estimation(datafile=fsdat,nobs=192,loglinear,mh_replic=2000,mh_nblocks=5,mh_jscale=0.8);
*/
options_.maxit_ = 100;
var c k /*s*/;
var c k;
varexo x;
parameters alph gam delt bet aa;
......@@ -10,12 +10,8 @@ bet=0.05;
aa=0.5;
//model(block, mfs=3);
model(block, bytecode, mfs=3);
//model;
//s = aa*x*k(-1)^alph - c;
c + k - aa*x*k(-1)^alph - (1-delt)*k(-1);// + 0.00000001*s;
//k = s - (1-delt)*k(-1);
model(block, bytecode);
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);
end;
......@@ -25,7 +21,7 @@ k = ((delt+bet)/(1.0*aa*alph))^(1/(alph-1));
c = aa*k^alph-delt*k;
end;
steady;
steady(solve_algo = 5);
//check;
model_info;
......@@ -35,7 +31,7 @@ periods 1;
values 1.02;
end;
simul(periods=200, METHOD=LU);
simul(periods = 200, stack_solve_algo = 5, markowitz = 2);
rplot c;
rplot k;
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment