New release of deterministic simulation.

git-svn-id: https://www.dynare.org/svn/dynare/dynare_v4@1292 ac1d8469-bf42-47a9-8791-bf33cf982152
time-shift
ferhat 2007-06-01 11:43:49 +00:00
parent 56da813fec
commit 9af2e461c3
34 changed files with 18667 additions and 1853 deletions

Binary file not shown.

View File

@ -2,9 +2,9 @@ function erased_compiled_function(func)
% erase compiled function with name 'func'
if exist([func '.dll'])
clear [func '.dll']
delete [func '.dll']
clear([func '.dll'])
delete([func '.dll'])
elseif exist ([func '.mexw32'])
clear [func '.dll']
delete [func '.mexw32']
clear([func '.dll'])
delete([func '.mexw32'])
end

64
matlab/formdata.m Normal file
View File

@ -0,0 +1,64 @@
function formdata(fname,date)
% function formdata(fname,date)
% store endogenous and exogenous variables in a "FRM" TROLL text format file
% INPUT
% fname: name of the FRM file
% date: the date of first observation (i.e. 2007A for an annual dataset)
% OUTPUT
% none
% ALGORITHM
% none
% SPECIAL REQUIREMENT
% none
%
% part of DYNARE, copyright (2007)
% Gnu Public License.
global M_ oo_
fid = fopen([fname '_endo.frm'],'w');
n=size(oo_.endo_simul,1);
t=size(oo_.endo_simul,2);
SN=upper(cellstr(M_.endo_names));
for i=1:n
str=strvcat(SN(i));
fprintf(fid,'USER: x x DATAFILE: x %s\n',str);
fprintf(fid,'PER: 1 YEAR: %s FRAC: 1 NOBS: %d CLINES: 0 DLINES: ???\n',date,t);
fprintf(fid,'%10.5f %10.5f %10.5f %10.5f\n',reshape(oo_.endo_simul(i,1:floor(t/4)*4),floor(t/4),4));
if(floor(t/4)*4<t)
switch(t-floor(t/4)*4)
case 1
fprintf(fid,'%10.5f\n',oo_.endo_simul(i,floor(t/4)*4+1:t));
case 2
fprintf(fid,'%10.5f %10.5f\n',oo_.endo_simul(i,floor(t/4)*4+1:t));
case 3
fprintf(fid,'%10.5f %10.5f %10.5f\n',oo_.endo_simul(i,floor(t/4)*4+1:t));
end;
%else
% fprintf(fid,'\n');
end;
end;
fclose(fid);
fid = fopen([fname '_exo.frm'],'w');
n=size(oo_.exo_simul,2);
t=size(oo_.exo_simul,1);
SN=upper(cellstr(M_.exo_names));
for i=1:n
str=strvcat(SN(i));
fprintf(fid,'USER: x x DATAFILE: x %s\n',str);
fprintf(fid,'PER: 1 YEAR: %s FRAC: 1 NOBS: %d CLINES: 0 DLINES: ???\n',date,t);
fprintf(fid,'%10.5f %10.5f %10.5f %10.5f\n',reshape(oo_.exo_simul(1:floor(t/4)*4,i),floor(t/4),4));
if(floor(t/4)*4<t)
switch(t-floor(t/4)*4)
case 1
fprintf(fid,'%10.5f\n',oo_.exo_simul(floor(t/4)*4+1:t,i)');
case 2
fprintf(fid,'%10.5f %10.5f\n',oo_.exo_simul(floor(t/4)*4+1:t,i)');
case 3
fprintf(fid,'%10.5f %10.5f %10.5f\n',oo_.exo_simul(floor(t/4)*4+1:t,i)');
end;
%else
% fprintf(fid,'\n');
end;
end;
fclose(fid);
return;

48
matlab/read_data_.m Normal file
View File

@ -0,0 +1,48 @@
function read_data_
% function read_data_
% read endogenous and exogenous variables from a text file
% Used by datafile option in simulate
% INPUT
% none
% OUTPUT
% none
% ALGORITHM
% none
% SPECIAL REQUIREMENT
% none
%
% part of DYNARE, copyright (2007)
% Gnu Public License.
global options_ M_ oo_;
dname= options_.datafile;
fid = fopen([dname '_endo.dat'],'r');
names_line = fgetl(fid);
allVariables = '';
positions = ones(0);
while (any(names_line))
[chopped,names_line] = strtok(names_line);
allVariables = strvcat(allVariables, chopped);
positions = [positions ; strmatch(chopped,M_.endo_names,'exact')];
end
Values=fscanf(fid,'%f',inf);
Values=reshape(Values,M_.endo_nbr,size(Values,1)/M_.endo_nbr);
oo_.endo_simul=Values(positions,:);
fclose(fid);
fid = fopen([dname '_exo.dat'],'r');
names_line = fgetl(fid);
allVariables = '';
positions = ones(0);
while (any(names_line))
[chopped,names_line] = strtok(names_line);
allVariables = strvcat(allVariables, chopped);
positions = [positions ; strmatch(chopped,M_.exo_names,'exact')];
end
Values=fscanf(fid,'%f',inf);
Values=reshape(Values,M_.exo_nbr,size(Values,1)/M_.exo_nbr);
oo_.exo_simul=(Values(positions,:))';
fclose(fid);
%disp([allVariables M_.endo_names]);
%disp(positions);
end

32
matlab/writedata.m Normal file
View File

@ -0,0 +1,32 @@
function writedata(fname)
% function writedata(fname)
% store endogenous and exogenous variables in a XLS spreadsheet file
% INPUT
% fname: name of the XLS file
% OUTPUT
% none
% ALGORITHM
% none
% SPECIAL REQUIREMENT
% none
%
% part of DYNARE, copyright (2007)
% Gnu Public License.
global M_ oo_
S=[fname '_endo.xls'];
n=size(oo_.endo_simul,2);
delete(S);
S=upper(cellstr(M_.endo_names));
S1=cellstr([num2str((1:n)') char(65*ones(1,n))']);
xlswrite([fname '_endo'], S', 'endogenous', 'B1');
xlswrite([fname '_endo'], S1, 'endogenous', 'A2');
xlswrite([fname '_endo'], oo_.endo_simul', 'endogenous', 'B2');
S=[fname '_exo.xls'];
n=size(oo_.exo_simul,1);
delete(S);
S=upper(cellstr(M_.exo_names));
S1=cellstr([num2str((1:n)') char(65*ones(1,n))']);
xlswrite([fname '_exo'], S','exogenous', 'B1');
xlswrite([fname '_exo'], S1, 'exogenous', 'A2');
xlswrite([fname '_exo'], oo_.exo_simul,'exogenous', 'B2');
return;

32
matlab/writedata_text.m Normal file
View File

@ -0,0 +1,32 @@
function writedata(fname)
% function writedata(fname)
% store endogenous and exogenous variables in a XLS spreadsheet file
% INPUT
% fname: name of the XLS file
% OUTPUT
% none
% ALGORITHM
% none
% SPECIAL REQUIREMENT
% none
%
% part of DYNARE, copyright (2007)
% Gnu Public License.
global M_ oo_
S=[fname '_endo.xls'];
n=size(oo_.endo_simul,2);
delete(S);
S=upper(cellstr(M_.endo_names));
S1=cellstr([num2str((1:n)') char(65*ones(1,n))']);
xlswrite([fname '_endo'], S', 'endogenous', 'B1');
xlswrite([fname '_endo'], S1, 'endogenous', 'A2');
xlswrite([fname '_endo'], oo_.endo_simul', 'endogenous', 'B2');
S=[fname '_exo.xls'];
n=size(oo_.exo_simul,1);
delete(S);
S=upper(cellstr(M_.exo_names));
S1=cellstr([num2str((1:n)') char(65*ones(1,n))']);
xlswrite([fname '_exo'], S','exogenous', 'B1');
xlswrite([fname '_exo'], S1, 'exogenous', 'A2');
xlswrite([fname '_exo'], oo_.exo_simul,'exogenous', 'B2');
return;

View File

@ -13,13 +13,15 @@
#include <time.h>
#include <malloc.h>
#include <string.h>
#include <math.h>
using namespace std;
//------------------------------------------------------------------------------
#include "BlockTriangular.hh"
//------------------------------------------------------------------------------
BlockTriangular::BlockTriangular(const SymbolTable &symbol_table_arg) :
symbol_table(symbol_table_arg)
symbol_table(symbol_table_arg),
normalization(symbol_table_arg)
{
bt_verbose = 0;
ModelBlock = NULL;
@ -699,11 +701,162 @@ BlockTriangular::Free_Block(Model_Block* ModelBlock)
free(ModelBlock);
}
string
BlockTriangular::getnamebyID(Type type, int id)
{
return symbol_table.getNameByID(type,id);
}
//------------------------------------------------------------------------------
// Normalize each equation of the model (endgenous_i = f_i(endogenous_1, ..., endogenous_n) - in order to apply strong connex components search algorithm -
// and find the optimal blocks triangular decomposition
bool
BlockTriangular::Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock, int n, int* prologue, int* epilogue, simple* Index_Var_IM, simple* Index_Equ_IM, bool Do_Normalization, bool mixing, bool* IM_0 )
BlockTriangular::Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock, int n, int* prologue, int* epilogue, simple* Index_Var_IM, simple* Index_Equ_IM, bool Do_Normalization, bool mixing, bool* IM_0, jacob_map j_m )
{
int i, j, Nb_TotalBlocks, Nb_RecursBlocks;
int count_Block, count_Equ;
block_result_t* res;
List_IM * p_First_IM, *p_Cur_IM, *Cur_IM;
Equation_set* Equation_gr = (Equation_set*) malloc(sizeof(Equation_set));
bool* SIM0, *SIM00;
p_First_IM = (List_IM*)malloc(sizeof(*p_First_IM));
p_Cur_IM = p_First_IM;
Cur_IM = First_IM;
i = endo_nbr * endo_nbr;
while(Cur_IM)
{
p_Cur_IM->lead_lag = Cur_IM->lead_lag;
p_Cur_IM->IM = (bool*)malloc(i * sizeof(bool));
memcpy ( p_Cur_IM->IM, Cur_IM->IM, i);
Cur_IM = Cur_IM->pNext;
if(Cur_IM)
{
p_Cur_IM->pNext = (List_IM*)malloc(sizeof(*p_Cur_IM));
p_Cur_IM = p_Cur_IM->pNext;
}
else
p_Cur_IM->pNext = NULL;
}
Prologue_Epilogue(IM, prologue, epilogue, n, Index_Var_IM, Index_Equ_IM);
if(bt_verbose)
{
cout << "prologue : " << *prologue << " epilogue : " << *epilogue << "\n";
Print_SIM(IM, n);
for(i = 0;i < n;i++)
cout << "Index_Var_IM[" << i << "]=" << Index_Var_IM[i].index << " Index_Equ_IM[" << i << "]=" << Index_Equ_IM[i].index << "\n";
}
if(Do_Normalization)
{
cout << "Normalizing the model ...\n";
if(mixing)
{
double* max_val=(double*)malloc(n*sizeof(double));
memset(max_val,0,n*sizeof(double));
for( map< pair< int, int >, double >::iterator iter = j_m.begin(); iter != j_m.end(); iter++ )
{
if(fabs(iter->second)>max_val[iter->first.first])
max_val[iter->first.first]=fabs(iter->second);
}
for( map< pair< int, int >, double >::iterator iter = j_m.begin(); iter != j_m.end(); iter++ )
iter->second/=max_val[iter->first.first];
free(max_val);
bool OK=false;
double bi=0.99999999;
int suppressed=0;
while(!OK && bi>1e-14)
{
int suppress=0;
SIM0 = (bool*)malloc(n * n * sizeof(bool));
memset(SIM0,0,n*n*sizeof(bool));
SIM00 = (bool*)malloc(n * n * sizeof(bool));
memset(SIM00,0,n*n*sizeof(bool));
for( map< pair< int, int >, double >::iterator iter = j_m.begin(); iter != j_m.end(); iter++ )
{
if(fabs(iter->second)>bi)
{
SIM0[iter->first.first*n+iter->first.second]=1;
if(!IM_0[iter->first.first*n+iter->first.second])
{
cout << "error nothing at IM_0[" << iter->first.first << ", " << iter->first.second << "]=" << IM_0[iter->first.first*n+iter->first.second] << "\n";
}
}
else
suppress++;
}
for(i = 0;i < n;i++)
for(j = 0;j < n;j++)
SIM00[i*n + j] = SIM0[Index_Equ_IM[i].index * n + Index_Var_IM[j].index];
free(SIM0);
if(suppress!=suppressed)
OK=normalization.Normalize(n, *prologue, *epilogue, SIM00, Index_Equ_IM, Equation_gr, 1, IM);
suppressed=suppress;
if(!OK)
bi/=1.05;
if(bi>1e-14)
free(SIM00);
}
if(!OK)
{
normalization.Set_fp_verbose(true);
OK=normalization.Normalize(n, *prologue, *epilogue, SIM00, Index_Equ_IM, Equation_gr, 1, IM);
cout << "Error\n";
exit(-1);
}
}
else
normalization.Normalize(n, *prologue, *epilogue, IM, Index_Equ_IM, Equation_gr, 0, 0);
}
else
normalization.Gr_to_IM_basic(n, *prologue, *epilogue, IM, Equation_gr, false);
cout << "Finding the optimal block decomposition of the model ...\n";
if(bt_verbose)
blocks.Print_Equation_gr(Equation_gr);
res = blocks.sc(Equation_gr);
if(bt_verbose)
blocks.block_result_print(res);
if ((*prologue) || (*epilogue))
j = 1;
else
j = 0;
for(i = 0;i < res->n_sets;i++)
if ((res->sets_f[i] - res->sets_s[i] + 1) > j)
j = res->sets_f[i] - res->sets_s[i] + 1;
Nb_RecursBlocks = *prologue + *epilogue;
Nb_TotalBlocks = res->n_sets + Nb_RecursBlocks;
cout << Nb_TotalBlocks << " block(s) found:\n";
cout << " " << Nb_RecursBlocks << " recursive block(s) and " << res->n_sets << " simultaneous block(s). \n";
cout << " the largest simultaneous block has " << j << " equation(s). \n";
ModelBlock->Size = Nb_TotalBlocks;
ModelBlock->Periods = periods;
ModelBlock->in_Block_Equ = (int*)malloc(n * sizeof(int));
ModelBlock->in_Block_Var = (int*)malloc(n * sizeof(int));
ModelBlock->in_Equ_of_Block = (int*)malloc(n * sizeof(int));
ModelBlock->in_Var_of_Block = (int*)malloc(n * sizeof(int));
ModelBlock->Block_List = (Block*)malloc(sizeof(ModelBlock->Block_List[0]) * Nb_TotalBlocks);
blocks.block_result_to_IM(res, IM, *prologue, endo_nbr, Index_Equ_IM, Index_Var_IM);
Free_IM(p_First_IM);
count_Equ = count_Block = 0;
if (*prologue)
Allocate_Block(*prologue, &count_Equ, &count_Block, PROLOGUE, ModelBlock, Table, TableSize);
for(j = 0;j < res->n_sets;j++)
{
if(res->sets_f[res->ordered[j]] == res->sets_s[res->ordered[j]])
Allocate_Block(res->sets_f[res->ordered[j]] - res->sets_s[res->ordered[j]] + 1, &count_Equ, &count_Block, PROLOGUE, ModelBlock, Table, TableSize);
else
Allocate_Block(res->sets_f[res->ordered[j]] - res->sets_s[res->ordered[j]] + 1, &count_Equ, &count_Block, SIMULTANS, ModelBlock, Table, TableSize);
}
if (*epilogue)
Allocate_Block(*epilogue, &count_Equ, &count_Block, EPILOGUE, ModelBlock, Table, TableSize);
return 0;
}
//------------------------------------------------------------------------------
// Normalize each equation of the model (endgenous_i = f_i(endogenous_1, ..., endogenous_n) - in order to apply strong connex components search algorithm -
// and find the optimal blocks triangular decomposition
bool
BlockTriangular::Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock, int n, int* prologue, int* epilogue, simple* Index_Var_IM, simple* Index_Equ_IM, bool Do_Normalization, bool mixing, bool* IM_0)
{
int i, j, Nb_TotalBlocks, Nb_RecursBlocks;
int count_Block, count_Equ;
@ -743,6 +896,7 @@ BlockTriangular::Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock,
{
bool* SIM0;
SIM0 = (bool*)malloc(n * n * sizeof(bool));
for(i = 0;i < n*n;i++)
SIM0[i] = IM_0[i];
for(i = 0;i < n;i++)
@ -800,7 +954,6 @@ BlockTriangular::Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock,
//------------------------------------------------------------------------------
// For the contemparenous simultaneities
// normalize each equation of the model
@ -944,7 +1097,7 @@ BlockTriangular::Normalize_and_BlockDecompose_Static_Model()
// normalize each equation of the dynamic model
// and find the optimal block triangular decomposition of the static model
void
BlockTriangular::Normalize_and_BlockDecompose_Static_0_Model()
BlockTriangular::Normalize_and_BlockDecompose_Static_0_Model(const jacob_map &j_m)
{
bool* SIM, *SIM_0;
List_IM* Cur_IM;
@ -985,7 +1138,7 @@ BlockTriangular::Normalize_and_BlockDecompose_Static_0_Model()
SIM_0 = (bool*)malloc(endo_nbr * endo_nbr * sizeof(*SIM_0));
for(i = 0;i < endo_nbr*endo_nbr;i++)
SIM_0[i] = Cur_IM->IM[i];
Normalize_and_BlockDecompose(SIM, ModelBlock, endo_nbr, &prologue, &epilogue, Index_Var_IM, Index_Equ_IM, 1, 1, SIM_0);
Normalize_and_BlockDecompose(SIM, ModelBlock, endo_nbr, &prologue, &epilogue, Index_Var_IM, Index_Equ_IM, 1, 1, SIM_0, j_m);
if(bt_verbose)
for(i = 0;i < endo_nbr;i++)
cout << "Block=" << Index_Equ_IM[i].block << " Equ=" << Index_Equ_IM[i].index << " Var= " << Index_Var_IM[i].index << " " << symbol_table.getNameByID(eEndogenous, Index_Var_IM[i].index) << "\n";

View File

@ -78,14 +78,20 @@ SimulSparseStatement::writeOutput(ostream &output, const string &basename) const
{
options_list.writeOutput(output);
output << "if ~ options_.initval_file\n";
output << " make_y_;\n";
output << " make_ex_;\n";
output << " if ~isfield(options_,'datafile')\n";
output << " make_y_;\n";
output << " make_ex_;\n";
output << " else\n";
output << " read_data_;\n";
output << " end\n";
output << "end\n";
output << "disp('compiling...');\n";
output << "t0=clock;\n";
if (compiler == 0)
output << "mex " << basename << "_dynamic.c;\n";
else
output << "mex " << basename << "_dynamic.cc;\n";
output << "disp(['compiling time: ' num2str(etime(clock,t0))]);\n";
output << "oo_.endo_simul=" << basename << "_dynamic;\n";
}
@ -236,7 +242,7 @@ PeriodsStatement::writeOutput(ostream &output, const string &basename) const
output << "options_.simul = 1;" << endl;
}
CutoffStatement::CutoffStatement(int cutoff_arg) : cutoff(cutoff_arg)
CutoffStatement::CutoffStatement(double cutoff_arg) : cutoff(cutoff_arg)
{
}

File diff suppressed because it is too large Load Diff

View File

@ -36,7 +36,7 @@ class ParsingDriver;
%token AR AUTOCORR
%token BAYESIAN_IRF BETA_PDF
%token CALIB CALIB_VAR CHECK CONF_SIG CONSTANT CORR COVAR CUTOFF
%token DATAFILE DR_ALGO DROP DSAMPLE DYNASAVE DYNATYPE
%token DATAFILE DR_ALGO DROP DSAMPLE DYNASAVE DYNATYPE
%token END ENDVAL EQUAL ESTIMATION ESTIMATED_PARAMS ESTIMATED_PARAMS_BOUNDS ESTIMATED_PARAMS_INIT
%token FILENAME FILTER_STEP_AHEAD FILTERED_VARS FIRST_OBS
%token <string_val> FLOAT_NUMBER
@ -50,14 +50,14 @@ class ParsingDriver;
%token LAPLACE LCC_COMPILER LIK_ALGO LIK_INIT LINEAR LOAD_MH_FILE LOGLINEAR
%token MH_DROP MH_INIT_SCALE MH_JSCALE MH_MODE MH_NBLOCKS MH_REPLIC MH_RECOVER
%token MODE_CHECK MODE_COMPUTE MODE_FILE MODEL MODEL_COMPARISON MSHOCKS
%token MODEL_COMPARISON_APPROXIMATION MODIFIEDHARMONICMEAN MOMENTS_VARENDO
%token MODEL_COMPARISON_APPROXIMATION MODIFIEDHARMONICMEAN MOMENTS_VARENDO
%token <string_val> NAME
%token NOBS NOCONSTANT NOCORR NODIAGNOSTIC NOFUNCTIONS NOGRAPH NOMOMENTS NOPRINT NORMAL_PDF
%token OBSERVATION_TRENDS OLR OLR_INST OLR_BETA OPTIM OPTIM_WEIGHTS ORDER OSR OSR_PARAMS
%token OBSERVATION_TRENDS OLR OLR_INST OLR_BETA OPTIM OPTIM_WEIGHTS ORDER OSR OSR_PARAMS
%token PARAMETERS PERIODS PLANNER_OBJECTIVE PREFILTER PRESAMPLE PRINT PRIOR_TRUNC PRIOR_ANALYSIS POSTERIOR_ANALYSIS
%token QZ_CRITERIUM
%token RELATIVE_IRF REPLIC RPLOT
%token SHOCKS SIGMA_E SIMUL SIMUL_ALGO SIMUL_SEED SMOOTHER SOLVE_ALGO SPARSE_DLL STDERR STEADY STOCH_SIMUL
%token SHOCKS SIGMA_E SIMUL SIMUL_ALGO SIMUL_SEED SMOOTHER SOLVE_ALGO SPARSE_DLL STDERR STEADY STOCH_SIMUL
%token TEX RAMSEY_POLICY PLANNER_DISCOUNT
%token <string_val> TEX_NAME
%token UNIFORM_PDF UNIT_ROOT_VARS USE_DLL
@ -67,7 +67,7 @@ class ParsingDriver;
%left PLUS MINUS
%left TIMES DIVIDE
%left UMINUS
%nonassoc POWER
%nonassoc POWER
%token EXP LOG LOG10 SIN COS TAN ASIN ACOS ATAN SINH COSH TANH ASINH ACOSH ATANH SQRT
%type <node_val> expression
@ -126,7 +126,7 @@ class ParsingDriver;
| ramsey_policy
;
declaration
: parameters
| var
@ -134,44 +134,44 @@ class ParsingDriver;
| varexo_det
;
dsample : DSAMPLE INT_NUMBER ';' { driver.dsample($2);}
| DSAMPLE INT_NUMBER INT_NUMBER ';' {driver.dsample($2, $3);}
;
;
rplot : RPLOT tmp_var_list ';' {driver.rplot();}
;
var
: VAR var_list ';'
;
var
: VAR var_list ';'
;
varexo
varexo
: VAREXO varexo_list ';'
;
varexo_det
: VAREXO_DET varexo_det_list ';'
;
parameters
: PARAMETERS parameter_list ';'
;
var_list
: var_list NAME
: var_list NAME
{ driver.declare_endogenous($2); }
| var_list COMMA NAME
| var_list COMMA NAME
{ driver.declare_endogenous($3); }
| NAME
{ driver.declare_endogenous($1); }
| var_list NAME TEX_NAME
| var_list NAME TEX_NAME
{ driver.declare_endogenous($2, $3); }
| var_list COMMA NAME TEX_NAME
{ driver.declare_endogenous($3, $4); }
| NAME TEX_NAME
{ driver.declare_endogenous($1, $2); }
;
;
varexo_list
: varexo_list NAME
{ driver.declare_exogenous($2); }
@ -185,7 +185,7 @@ class ParsingDriver;
{ driver.declare_exogenous($3, $4); }
| NAME TEX_NAME
{ driver.declare_exogenous($1, $2); }
;
;
varexo_det_list
: varexo_det_list NAME
@ -200,7 +200,7 @@ class ParsingDriver;
{ driver.declare_exogenous_det($3, $4); }
| NAME TEX_NAME
{ driver.declare_exogenous_det($1, $2); }
;
;
parameter_list
: parameter_list NAME
@ -215,9 +215,9 @@ class ParsingDriver;
{ driver.declare_parameter($3, $4); }
| NAME TEX_NAME
{ driver.declare_parameter($1, $2); }
;
;
periods
periods
: PERIODS INT_NUMBER ';'
{
driver.periods($2);
@ -238,12 +238,12 @@ cutoff
driver.cutoff($3);
}
;
init_param
: NAME EQUAL expression ';'
{driver.init_param($1, $3);}
{driver.init_param($1, $3);}
;
expression
: '(' expression ')'
{ $$ = $2;}
@ -253,15 +253,15 @@ cutoff
{$$ = driver.add_constant($1);}
| INT_NUMBER
{$$ = driver.add_constant($1);}
| expression PLUS expression
| expression PLUS expression
{$$ = driver.add_plus($1, $3);}
| expression MINUS expression
{$$ = driver.add_minus($1, $3);}
| expression DIVIDE expression
| expression DIVIDE expression
{$$ = driver.add_divide($1, $3);}
| expression TIMES expression
| expression TIMES expression
{$$ = driver.add_times($1, $3);}
| expression POWER expression
| expression POWER expression
{$$ = driver.add_power($1, $3);}
| MINUS expression %prec UMINUS
{$$ = driver.add_uminus($2);}
@ -287,9 +287,9 @@ cutoff
{$$ = driver.add_atan($3);}
| SQRT '(' expression ')'
{$$ = driver.add_sqrt($3);}
| NAME '(' comma_expression ')'
| NAME '(' comma_expression ')'
{$$ = driver.add_unknown_function($1);}
;
;
comma_expression :
expression
@ -307,22 +307,22 @@ cutoff
initval_option
: FILENAME EQUAL NAME {driver.init_val_filename($3);}
;
endval
: ENDVAL ';' initval_list END
{driver.end_endval();}
;
;
initval_list
: initval_list initval_elem
| initval_elem
;
initval_elem
initval_elem
: NAME EQUAL expression ';'
{driver.init_val($1, $3);}
;
;
histval
: HISTVAL ';' histval_list END
{ driver.end_histval(); }
@ -332,7 +332,7 @@ cutoff
: histval_list histval_elem
| histval_elem
;
histval_elem
: NAME '(' signed_integer ')' EQUAL expression ';'
{driver.hist_val($1, $3, $6);}
@ -347,10 +347,10 @@ cutoff
| GCC_COMPILER { driver.init_compiler(1); }
| o_cutoff
;
model
: MODEL ';' { driver.begin_model(); } equation_list END { driver.reset_data_tree(); }
| MODEL '(' o_linear ')' ';' { driver.begin_model(); }
| MODEL '(' o_linear ')' ';' { driver.begin_model(); }
equation_list END { driver.reset_data_tree(); }
| MODEL '(' USE_DLL ')' ';' { driver.begin_model(); driver.use_dll(); }
equation_list END { driver.reset_data_tree(); }
@ -361,19 +361,19 @@ cutoff
;
equation_list
: equation_list equation
: equation_list equation
| equation_list pound_expression
| equation
| pound_expression
;
equation
: hand_side EQUAL hand_side ';'
{$$ = driver.add_model_equal($1, $3);}
| hand_side ';'
{$$ = driver.add_model_equal_with_zero_rhs($1);}
;
hand_side
: '(' hand_side ')' {$$ = $2;}
| model_var
@ -381,16 +381,16 @@ cutoff
{$$ = driver.add_constant($1);}
| INT_NUMBER
{$$ = driver.add_constant($1);}
| hand_side PLUS hand_side
| hand_side PLUS hand_side
{$$ = driver.add_plus($1, $3);}
| hand_side MINUS hand_side
{$$ = driver.add_minus($1, $3);}
| hand_side DIVIDE hand_side
| hand_side DIVIDE hand_side
{$$ = driver.add_divide($1, $3);}
| hand_side TIMES hand_side
| hand_side TIMES hand_side
{$$ = driver.add_times($1, $3);}
| hand_side POWER hand_side
{$$ = driver.add_power($1, $3);}
| hand_side POWER hand_side
{$$ = driver.add_power($1, $3);}
| MINUS hand_side %prec UMINUS
{ $$ = driver.add_uminus($2);}
| PLUS hand_side
@ -416,17 +416,17 @@ cutoff
| SQRT '(' hand_side ')'
{$$ = driver.add_sqrt($3);}
;
pound_expression: '#' NAME EQUAL hand_side ';'
{driver.declare_and_init_model_local_variable($2, $4);}
model_var
: NAME
: NAME
{$$ = driver.add_model_variable($1);}
| NAME '(' signed_integer ')'
{$$ = driver.add_model_variable($1, $3);}
;
shocks
: SHOCKS ';' shock_list END {driver.end_shocks();}
;
@ -435,18 +435,18 @@ cutoff
: MSHOCKS ';' shock_list END {driver.end_mshocks();}
;
shock_list
shock_list
: shock_list shock_elem
| shock_elem
;
shock_elem
shock_elem
: VAR NAME ';' PERIODS period_list ';' VALUES value_list ';'
{driver.add_det_shock($2);}
| VAR NAME ';' STDERR expression ';'
{driver.add_stderr_shock($2, $5);}
| VAR NAME EQUAL expression ';'
{driver.add_var_shock($2, $4);}
{driver.add_var_shock($2, $4);}
| VAR NAME COMMA NAME EQUAL expression ';'
{driver.add_covar_shock($2, $4, $6);}
| CORR NAME COMMA NAME EQUAL expression ';'
@ -456,11 +456,11 @@ cutoff
period_list
: period_list INT_NUMBER
{driver.add_period($2);}
| period_list INT_NUMBER ':' INT_NUMBER
| period_list INT_NUMBER ':' INT_NUMBER
{driver.add_period($2,$4);}
| period_list COMMA INT_NUMBER
{driver.add_period($3);}
| period_list COMMA INT_NUMBER ':' INT_NUMBER
| period_list COMMA INT_NUMBER ':' INT_NUMBER
{driver.add_period($3, $5);}
| INT_NUMBER ':' INT_NUMBER
{driver.add_period($1, $3);}
@ -470,7 +470,7 @@ cutoff
value_list
: value_list signed_float
: value_list signed_float
{driver.add_value_const($2);}
| value_list signed_integer
{driver.add_value_const($2);}
@ -487,49 +487,49 @@ cutoff
| '(' expression ')'
{driver.add_value($2);}
;
sigma_e
: SIGMA_E EQUAL '[' triangular_matrix ']' ';'
sigma_e
: SIGMA_E EQUAL '[' triangular_matrix ']' ';'
{driver.do_sigma_e();}
;
triangular_matrix
: triangular_matrix ';' triangular_row
triangular_matrix
: triangular_matrix ';' triangular_row
{driver.end_of_row();}
| triangular_row
| triangular_row
{driver.end_of_row();}
;
triangular_row
: triangular_row COMMA '(' expression ')'
triangular_row
: triangular_row COMMA '(' expression ')'
{driver.add_to_row($4);}
| triangular_row COMMA FLOAT_NUMBER
| triangular_row COMMA FLOAT_NUMBER
{driver.add_to_row_const($3);}
| triangular_row COMMA INT_NUMBER
| triangular_row COMMA INT_NUMBER
{driver.add_to_row_const($3);}
| triangular_row '(' expression ')'
{driver.add_to_row($3);}
| triangular_row FLOAT_NUMBER
| triangular_row FLOAT_NUMBER
{driver.add_to_row_const($2);}
| triangular_row INT_NUMBER
| triangular_row INT_NUMBER
{driver.add_to_row_const($2);}
| '(' expression ')'
{driver.add_to_row($2);}
| FLOAT_NUMBER
| FLOAT_NUMBER
{driver.add_to_row_const($1);}
| INT_NUMBER
| INT_NUMBER
{driver.add_to_row_const($1);}
;
steady
: STEADY ';'
steady
: STEADY ';'
{
driver.steady();
}
| STEADY '(' steady_options_list ')' ';'
{driver.steady();}
;
steady_options_list : steady_options_list COMMA steady_options
| steady_options
;
@ -537,11 +537,11 @@ cutoff
steady_options: o_solve_algo
;
check
: CHECK ';'
check
: CHECK ';'
{driver.check();}
| CHECK '(' check_options_list ')' ';'
{driver.check();}
| CHECK '(' check_options_list ')' ';'
{driver.check();}
;
check_options_list : check_options_list COMMA check_options
@ -551,8 +551,8 @@ cutoff
check_options : o_solve_algo
;
simul
: SIMUL ';'
simul
: SIMUL ';'
{driver.simulate();}
| SIMUL '(' simul_options_list ')' ';'
{driver.simulate();}
@ -563,19 +563,20 @@ cutoff
;
simul_options: o_periods
| o_datafile
;
stoch_simul
: STOCH_SIMUL ';'
{driver.stoch_simul();}
| STOCH_SIMUL '(' stoch_simul_options_list ')' ';'
| STOCH_SIMUL '(' stoch_simul_options_list ')' ';'
{driver.stoch_simul();}
| STOCH_SIMUL tmp_var_list ';'
{driver.stoch_simul();}
| STOCH_SIMUL '(' stoch_simul_options_list ')' tmp_var_list ';'
| STOCH_SIMUL '(' stoch_simul_options_list ')' tmp_var_list ';'
{driver.stoch_simul();}
;
stoch_simul_options_list: stoch_simul_options_list COMMA stoch_simul_options
| stoch_simul_options
;
@ -584,7 +585,7 @@ cutoff
| o_solve_algo
| o_simul_algo
| o_linear
| o_order
| o_order
| o_replic
| o_drop
| o_ar
@ -625,7 +626,7 @@ cutoff
| INT_NUMBER
{$$ = $1;}
;
signed_float
: PLUS FLOAT_NUMBER
{$$ = $2;}
@ -635,24 +636,24 @@ cutoff
{$$ = $1;}
;
estimated_params
estimated_params
: ESTIMATED_PARAMS ';' estimated_list END
{ driver.estimated_params(); }
;
estimated_list
: estimated_list estimated_elem
estimated_list
: estimated_list estimated_elem
{driver.add_estimated_params_element();}
| estimated_elem
| estimated_elem
{driver.add_estimated_params_element();}
;
estimated_elem
estimated_elem
: estimated_elem1 COMMA estimated_elem2 ';'
;
estimated_elem1
: STDERR NAME
estimated_elem1
: STDERR NAME
{driver.estim_params.type = 1;
driver.estim_params.name = *$2;
delete $2;
@ -671,19 +672,19 @@ cutoff
}
;
estimated_elem2
: prior COMMA estimated_elem3
estimated_elem2
: prior COMMA estimated_elem3
{
driver.estim_params.prior=*$1;
delete $1;
}
| value COMMA prior COMMA estimated_elem3
| value COMMA prior COMMA estimated_elem3
{driver.estim_params.init_val=*$1;
driver.estim_params.prior=*$3;
delete $1;
delete $3;
}
| value COMMA value COMMA value COMMA prior COMMA estimated_elem3
| value COMMA value COMMA value COMMA prior COMMA estimated_elem3
{driver.estim_params.init_val=*$1;
driver.estim_params.low_bound=*$3;
driver.estim_params.up_bound=*$5;
@ -693,12 +694,12 @@ cutoff
delete $5;
delete $7;
}
| value
| value
{
driver.estim_params.init_val=*$1;
delete $1;
}
| value COMMA value COMMA value
| value COMMA value COMMA value
{driver.estim_params.init_val=*$1;
driver.estim_params.low_bound=*$3;
driver.estim_params.up_bound=*$5;
@ -708,8 +709,8 @@ cutoff
}
;
estimated_elem3
: value COMMA value
estimated_elem3
: value COMMA value
{driver.estim_params.mean=*$1;
driver.estim_params.std=*$3;
delete $1;
@ -722,8 +723,8 @@ cutoff
delete $1;
delete $3;
delete $5;
}
| value COMMA value COMMA value COMMA value
}
| value COMMA value COMMA value COMMA value
{driver.estim_params.mean=*$1;
driver.estim_params.std=*$3;
driver.estim_params.p3=*$5;
@ -733,7 +734,7 @@ cutoff
delete $5;
delete $7;
}
| value COMMA value COMMA value COMMA value COMMA value
| value COMMA value COMMA value COMMA value COMMA value
{driver.estim_params.mean=*$1;
driver.estim_params.std=*$3;
driver.estim_params.p3=*$5;
@ -824,20 +825,20 @@ cutoff
;
prior
: BETA_PDF
: BETA_PDF
{$$ = new string("1");}
| GAMMA_PDF
| GAMMA_PDF
{$$ = new string("2");}
| NORMAL_PDF
| NORMAL_PDF
{$$ = new string("3");}
| INV_GAMMA_PDF
| INV_GAMMA_PDF
{$$ = new string("4");}
| UNIFORM_PDF
| UNIFORM_PDF
{$$ = new string("5");}
;
value
: {$$ = new string("NaN");}
value
: {$$ = new string("NaN");}
| INT_NUMBER
| FLOAT_NUMBER
| NAME
@ -847,20 +848,20 @@ cutoff
{$2->insert(0, "-"); $$ = $2;}
;
estimation
estimation
: ESTIMATION ';'
{driver.run_estimation();}
| ESTIMATION '(' estimation_options_list ')' ';'
| ESTIMATION '(' estimation_options_list ')' ';'
{driver.run_estimation();}
| ESTIMATION tmp_var_list ';'
{driver.run_estimation();}
| ESTIMATION '(' estimation_options_list ')' tmp_var_list ';'
| ESTIMATION '(' estimation_options_list ')' tmp_var_list ';'
{driver.run_estimation();}
;
estimation_options_list
estimation_options_list
: estimation_options_list COMMA estimation_options
| estimation_options
;
@ -870,22 +871,22 @@ cutoff
| o_first_obs
| o_prefilter
| o_presample
| o_lik_algo
| o_lik_init
| o_lik_algo
| o_lik_init
| o_nograph
| o_conf_sig
| o_conf_sig
| o_mh_replic
| o_mh_drop
| o_mh_jscale
| o_optim
| o_mh_init_scale
| o_mode_file
| o_mode_compute
| o_mh_init_scale
| o_mode_file
| o_mode_compute
| o_mode_check
| o_prior_trunc
| o_mh_mode
| o_mh_nblcks
| o_load_mh_file
| o_prior_trunc
| o_mh_mode
| o_mh_nblcks
| o_load_mh_file
| o_loglinear
| o_nodiagnostic
| o_bayesian_irf
@ -905,23 +906,23 @@ cutoff
| o_noconstant
| o_mh_recover
;
prior_analysis
: PRIOR_ANALYSIS '(' prior_posterior_options_list ')' ';'
prior_analysis
: PRIOR_ANALYSIS '(' prior_posterior_options_list ')' ';'
{driver.run_prior_analysis();}
| PRIOR_ANALYSIS '(' prior_posterior_options_list ')' tmp_var_list ';'
| PRIOR_ANALYSIS '(' prior_posterior_options_list ')' tmp_var_list ';'
{driver.run_prior_analysis();}
;
prior_posterior_options_list
prior_posterior_options_list
: prior_posterior_options_list COMMA prior_posterior_options
| prior_posterior_options
;
prior_posterior_options
prior_posterior_options
: o_nograph
| o_conf_sig
| o_prior_trunc
| o_conf_sig
| o_prior_trunc
| o_bayesian_irf
| o_irf
| o_tex
@ -933,11 +934,11 @@ cutoff
| o_xls_range
| o_filter_step_ahead
;
posterior_analysis
: POSTERIOR_ANALYSIS '(' prior_posterior_options_list ')' ';'
posterior_analysis
: POSTERIOR_ANALYSIS '(' prior_posterior_options_list ')' ';'
{driver.run_posterior_analysis();}
| POSTERIOR_ANALYSIS '(' prior_posterior_options_list ')' tmp_var_list ';'
| POSTERIOR_ANALYSIS '(' prior_posterior_options_list ')' tmp_var_list ';'
{driver.run_posterior_analysis();}
;
@ -950,9 +951,9 @@ cutoff
: list_optim_option
| optim_options COMMA list_optim_option;
;
varobs
: VAROBS tmp_var_list ';'
varobs
: VAROBS tmp_var_list ';'
{driver.set_varobs();}
;
@ -961,12 +962,12 @@ cutoff
{ driver.set_trends(); }
;
trend_list
trend_list
: trend_list trend_element
| trend_element
;
trend_element : NAME '(' expression ')' ';'
trend_element : NAME '(' expression ')' ';'
{driver.set_trend_element($1, $3);}
;
@ -977,7 +978,7 @@ cutoff
{ driver.optim_weights(); }
;
optim_weights_list : optim_weights_list NAME expression ';'
optim_weights_list : optim_weights_list NAME expression ';'
{driver.set_optim_weights($2, $3);}
| optim_weights_list NAME COMMA NAME expression ';'
{driver.set_optim_weights($2, $4, $5);}
@ -995,17 +996,17 @@ cutoff
| OSR tmp_var_list ';' {driver.run_osr();}
| OSR '(' olr_options ')' tmp_var_list ';' {driver.run_osr();}
;
olr : OLR ';' {driver.run_olr();}
| OLR '(' olr_options ')' ';' {driver.run_olr();}
| OLR tmp_var_list ';' {driver.run_olr();}
| OLR '(' olr_options ')' tmp_var_list ';' {driver.run_olr();}
;
olr_option : o_olr_beta
| stoch_simul_options
;
olr_options : olr_option
| olr_options COMMA olr_option
;
@ -1081,22 +1082,22 @@ cutoff
planner_objective : PLANNER_OBJECTIVE { driver.begin_planner_objective(); } hand_side { driver.end_planner_objective($3); } ';'
;
ramsey_policy : RAMSEY_POLICY ';'
{driver.ramsey_policy();}
| RAMSEY_POLICY '(' ramsey_policy_options_list ')' ';'
| RAMSEY_POLICY '(' ramsey_policy_options_list ')' ';'
{driver.ramsey_policy();}
| RAMSEY_POLICY tmp_var_list ';'
{driver.ramsey_policy();}
| RAMSEY_POLICY '(' ramsey_policy_options_list ')' tmp_var_list ';'
| RAMSEY_POLICY '(' ramsey_policy_options_list ')' tmp_var_list ';'
{driver.ramsey_policy();}
;
ramsey_policy_options_list :
ramsey_policy_options_list :
ramsey_policy_options_list COMMA ramsey_policy_options
| ramsey_policy_options
;
ramsey_policy_options : stoch_simul_options
| o_planner_discount
;
@ -1125,18 +1126,18 @@ cutoff
o_datafile: DATAFILE EQUAL NAME {driver.option_str("datafile", $3);};
o_nobs: NOBS EQUAL vec_int {driver.option_num("nobs", $3);}
| NOBS EQUAL INT_NUMBER {driver.option_num("nobs", $3);}
;
;
o_first_obs: FIRST_OBS EQUAL INT_NUMBER {driver.option_num("first_obs", $3);};
o_prefilter: PREFILTER EQUAL INT_NUMBER {driver.option_num("prefilter", $3);};
o_presample: PRESAMPLE EQUAL INT_NUMBER {driver.option_num("presample", $3);};
o_lik_algo: LIK_ALGO EQUAL INT_NUMBER {driver.option_num("lik_algo", $3);};
o_lik_init: LIK_INIT EQUAL INT_NUMBER {driver.option_num("lik_init", $3);};
o_nograph: NOGRAPH {driver.option_num("nograph","1");};
| GRAPH {driver.option_num("nograph", "0");};
o_conf_sig: CONF_SIG EQUAL FLOAT_NUMBER {driver.option_num("conf_sig", $3);};
o_mh_replic: MH_REPLIC EQUAL INT_NUMBER {driver.option_num("mh_replic", $3);};
o_mh_drop: MH_DROP EQUAL FLOAT_NUMBER {driver.option_num("mh_drop", $3);};
o_mh_jscale: MH_JSCALE EQUAL FLOAT_NUMBER {driver.option_num("mh_jscale", $3);};
o_lik_algo: LIK_ALGO EQUAL INT_NUMBER {driver.option_num("lik_algo", $3);};
o_lik_init: LIK_INIT EQUAL INT_NUMBER {driver.option_num("lik_init", $3);};
o_nograph: NOGRAPH {driver.option_num("nograph","1");};
| GRAPH {driver.option_num("nograph", "0");};
o_conf_sig: CONF_SIG EQUAL FLOAT_NUMBER {driver.option_num("conf_sig", $3);};
o_mh_replic: MH_REPLIC EQUAL INT_NUMBER {driver.option_num("mh_replic", $3);};
o_mh_drop: MH_DROP EQUAL FLOAT_NUMBER {driver.option_num("mh_drop", $3);};
o_mh_jscale: MH_JSCALE EQUAL FLOAT_NUMBER {driver.option_num("mh_jscale", $3);};
o_optim: OPTIM EQUAL '(' optim_options ')';
o_mh_init_scale: MH_INIT_SCALE EQUAL FLOAT_NUMBER {driver.option_num("mh_init_scale", $3);};
| MH_INIT_SCALE EQUAL INT_NUMBER {driver.option_num("mh_init_scale", $3);};
@ -1167,8 +1168,8 @@ cutoff
;
o_print : PRINT {driver.option_num("noprint", "0");};
o_noprint : NOPRINT {driver.option_num("noprint", "1");};
o_xls_sheet : XLS_SHEET EQUAL NAME {driver.option_str("xls_sheet", $3);}
o_xls_range : XLS_RANGE EQUAL range {driver.option_str("xls_range", $3);}
o_xls_sheet : XLS_SHEET EQUAL NAME {driver.option_str("xls_sheet", $3);}
o_xls_range : XLS_RANGE EQUAL range {driver.option_str("xls_range", $3);}
o_filter_step_ahead : FILTER_STEP_AHEAD EQUAL vec_int {driver.option_num("filter_step_ahead", $3);}
o_constant : CONSTANT {driver.option_num("noconstant", "0");}
o_noconstant : NOCONSTANT {driver.option_num("noconstant", "1");}

View File

@ -15,11 +15,11 @@ typedef yy::parser::token token;
#define yyterminate() return yy::parser::token_type (0);
int comment_caller;
/* Particular value : when sigma_e command is found
/* Particular value : when sigma_e command is found
this flag is set to 1, when command finished it is set to 0
*/
int sigma_e = 0;
%}
%}
%option case-insensitive noyywrap nounput batch debug never-interactive yylineno
@ -44,10 +44,10 @@ int sigma_e = 0;
/* Comments */
<INITIAL,DYNARE_STATEMENT,DYNARE_BLOCK>["%"].*
<INITIAL,DYNARE_STATEMENT,DYNARE_BLOCK>["/"]["/"].*
<INITIAL,DYNARE_STATEMENT,DYNARE_BLOCK>["/"]["/"].*
<INITIAL,DYNARE_STATEMENT,DYNARE_BLOCK>"/*" {comment_caller = YY_START; BEGIN COMMENT;}
<COMMENT>[^*\n]*
<COMMENT>[^*\n]*
<COMMENT>"*"+[^/\n]
<COMMENT>"*"+"/" {BEGIN comment_caller;}
@ -106,24 +106,24 @@ int sigma_e = 0;
<INITIAL>calib_var {BEGIN DYNARE_BLOCK; return token::CALIB_VAR;}
/* End of a Dynare block */
<DYNARE_BLOCK>end[ \t\n]*; {BEGIN INITIAL; return token::END;}
<DYNARE_BLOCK>end[ \t\n]*; {BEGIN INITIAL; return token::END;}
/* Inside of a Dynare statement */
<DYNARE_STATEMENT>datafile {return token::DATAFILE;}
<DYNARE_STATEMENT>nobs {return token::NOBS;}
<DYNARE_STATEMENT>first_obs {return token::FIRST_OBS;}
<DYNARE_STATEMENT>prefilter {return token::PREFILTER;}
<DYNARE_STATEMENT>presample {return token::PRESAMPLE;}
<DYNARE_STATEMENT>lik_algo {return token::LIK_ALGO;}
<DYNARE_STATEMENT>lik_init {return token::LIK_INIT;}
<DYNARE_STATEMENT>graph {return token::GRAPH;}
<DYNARE_STATEMENT>nograph {return token::NOGRAPH;}
<DYNARE_STATEMENT>print {return token::PRINT;}
<DYNARE_STATEMENT>noprint {return token::NOPRINT;}
<DYNARE_STATEMENT>conf_sig {return token::CONF_SIG;}
<DYNARE_STATEMENT>mh_replic {return token::MH_REPLIC;}
<DYNARE_STATEMENT>mh_drop {return token::MH_DROP;}
<DYNARE_STATEMENT>mh_jscale {return token::MH_JSCALE;}
<DYNARE_STATEMENT>prefilter {return token::PREFILTER;}
<DYNARE_STATEMENT>presample {return token::PRESAMPLE;}
<DYNARE_STATEMENT>lik_algo {return token::LIK_ALGO;}
<DYNARE_STATEMENT>lik_init {return token::LIK_INIT;}
<DYNARE_STATEMENT>graph {return token::GRAPH;}
<DYNARE_STATEMENT>nograph {return token::NOGRAPH;}
<DYNARE_STATEMENT>print {return token::PRINT;}
<DYNARE_STATEMENT>noprint {return token::NOPRINT;}
<DYNARE_STATEMENT>conf_sig {return token::CONF_SIG;}
<DYNARE_STATEMENT>mh_replic {return token::MH_REPLIC;}
<DYNARE_STATEMENT>mh_drop {return token::MH_DROP;}
<DYNARE_STATEMENT>mh_jscale {return token::MH_JSCALE;}
<DYNARE_STATEMENT>mh_init_scale {return token::MH_INIT_SCALE;}
<DYNARE_STATEMENT>mode_file {return token::MODE_FILE;}
<DYNARE_STATEMENT>mode_compute {return token::MODE_COMPUTE;}
@ -159,7 +159,7 @@ int sigma_e = 0;
<DYNARE_STATEMENT>[\$][^$]*[\$] {
strtok(yytext+1, "$");
yylval->string_val = new string(yytext + 1);
yylval->string_val = new string(yytext + 1);
return token::TEX_NAME;
}
@ -169,6 +169,7 @@ int sigma_e = 0;
<DYNARE_BLOCK>values {return token::VALUES;}
<DYNARE_BLOCK>corr {return token::CORR;}
<DYNARE_BLOCK>periods {return token::PERIODS;}
<DYNARE_BLOCK>cutoff {return token::CUTOFF;}
<DYNARE_BLOCK>filename {return token::FILENAME;}
<DYNARE_BLOCK>gamma_pdf {return token::GAMMA_PDF;}
<DYNARE_BLOCK>beta_pdf {return token::BETA_PDF;}
@ -216,7 +217,7 @@ int sigma_e = 0;
<DYNARE_STATEMENT,DYNARE_BLOCK>linear {return token::LINEAR;}
<DYNARE_STATEMENT,DYNARE_BLOCK>[,] {return token::COMMA;}
<DYNARE_STATEMENT,DYNARE_BLOCK>[:] {return yy::parser::token_type (yytext[0]);}
<DYNARE_STATEMENT,DYNARE_BLOCK>[\(\)] {return yy::parser::token_type (yytext[0]);}
<DYNARE_STATEMENT,DYNARE_BLOCK>[\(\)] {return yy::parser::token_type (yytext[0]);}
<DYNARE_STATEMENT,DYNARE_BLOCK>[\[] {return yy::parser::token_type (yytext[0]);}
<DYNARE_STATEMENT,DYNARE_BLOCK>[\]] {
if (sigma_e)
@ -261,12 +262,12 @@ int sigma_e = 0;
yylval->string_val = new string(yytext);
return token::INT_NUMBER;
}
/* an instruction starting with a recognized symbol (which is not a modfile local variable)
is passed as NAME,
otherwise it is a native statement until the end of the line
*/
<INITIAL>[A-Za-z_][A-Za-z0-9_]* {
<INITIAL>[A-Za-z_][A-Za-z0-9_]* {
if (driver.symbol_exists_and_is_not_modfile_local_variable(yytext))
{
BEGIN DYNARE_STATEMENT;

View File

@ -335,13 +335,17 @@ VariableNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
double
VariableNode::eval(const eval_context_type &eval_context) const throw (EvalException)
{
if (lag != 0)
/*if (lag != 0)
throw EvalException();*/
if(&eval_context==NULL)
throw EvalException();
eval_context_type::const_iterator it = eval_context.find(make_pair(symb_id, type));
if (it == eval_context.end())
throw EvalException();
{
cout << "Error: the variable or parameter (" << datatree.symbol_table.getNameByID( type, symb_id) << ") has not been initialized (in derivatives evaluation)\n";
cout.flush();
throw EvalException();
}
return it->second;
}
@ -517,7 +521,7 @@ UnaryOpNode::cost(const temporary_terms_type &temporary_terms, bool is_matlab) c
}
cerr << "Impossible case!" << endl;
exit(-1);
}
}
void
UnaryOpNode::computeTemporaryTerms(map<NodeID, int> &reference_count,
@ -967,7 +971,7 @@ BinaryOpNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
// Treat special case of power operator in C
if (op_code == oPower && (!OFFSET(output_type)))
{
output << "pow(";
output << "pow1(";
arg1->writeOutput(output, output_type, temporary_terms);
output << ",";
arg2->writeOutput(output, output_type, temporary_terms);
@ -1110,5 +1114,7 @@ UnknownFunctionNode::collectEndogenous(NodeID &Id)
double
UnknownFunctionNode::eval(const eval_context_type &eval_context) const throw (EvalException)
{
cout << "Unknown function\n";
cout.flush();
throw EvalException();
}

View File

@ -10,9 +10,11 @@
using namespace std;
Normalization::Normalization()
Normalization::Normalization(const SymbolTable &symbol_table_arg) :
symbol_table(symbol_table_arg)
{
//Empty
fp_verbose=false;
};
Normalization::~Normalization()
@ -488,10 +490,18 @@ Normalization::ErrorHandling(int n, bool* IM, simple* Index_Equ_IM)
}
}
void
Normalization::Set_fp_verbose(bool ok)
{
fp_verbose=ok;
}
bool
Normalization::Normalize(int n, int prologue, int epilogue, bool* IM, simple* Index_Equ_IM, Equation_set* Equation, bool mixing, bool* IM_s)
{
int matchingSize, effective_n;
int save_fp_verbose=fp_verbose;
fp_verbose = 0;
Variable_set* Variable = (Variable_set*) malloc(sizeof(Variable_set));
#ifdef DEBUG
@ -502,23 +512,40 @@ Normalization::Normalize(int n, int prologue, int epilogue, bool* IM, simple* In
MaximumMatching(Equation, Variable);
matchingSize = MeasureMatching(Equation);
effective_n = n - prologue - epilogue;
if(matchingSize < effective_n)
fp_verbose=save_fp_verbose;
if(matchingSize < effective_n && fp_verbose)
{
cout << "Error: dynare could not normalize the model\n";
ErrorHandling(n, IM, Index_Equ_IM);
system("PAUSE");
cout << "Error: dynare could not normalize the model.\n The following equations:\n - ";
int i;
for(i = 0; i < Equation->size; i++)
if(Equation->Number[i].matched == -1)
cout << i << " ";
cout << "\n and the following variables:\n - ";
for(i = 0; i < Variable->size; i++)
if(Variable->Number[i].matched == -1)
cout << symbol_table.getNameByID(eEndogenous, Index_Equ_IM[i].index) << " ";
cout << "\n could not be normalized\n";
//ErrorHandling(n, IM, Index_Equ_IM);
//system("PAUSE");
exit( -1);
}
Gr_to_IM(n, prologue, epilogue, IM, Index_Equ_IM, Equation, mixing, IM_s);
if(fp_verbose)
if(matchingSize >= effective_n )
{
OutputMatching(Equation);
for(int i = 0;i < n;i++)
cout << "Index_Equ_IM[" << i << "]=" << Index_Equ_IM[i].index /*<< " == " "Index_Var_IM[" << i << "]=" << Index_Var_IM[i].index*/ << "\n";
Gr_to_IM(n, prologue, epilogue, IM, Index_Equ_IM, Equation, mixing, IM_s);
if(fp_verbose)
{
OutputMatching(Equation);
for(int i = 0;i < n;i++)
cout << "Index_Equ_IM[" << i << "]=" << Index_Equ_IM[i].index /*<< " == " "Index_Var_IM[" << i << "]=" << Index_Var_IM[i].index*/ << "\n";
}
}
Free_Other(Variable);
#ifdef DEBUG
cout << "end of Normalize\n";
#endif
if(matchingSize < effective_n )
return(0);
else
return(1);
}

View File

@ -8,14 +8,14 @@
#include "Interface.hh"
#include "Model_Graph.hh"
#include "SymbolGaussElim.hh"
ModelTree::ModelTree(SymbolTable &symbol_table_arg,
NumericalConstants &num_constants_arg) :
DataTree(symbol_table_arg, num_constants_arg),
mode(eStandardMode),
compiler(LCC_COMPILE),
cutoff(1e-6),
cutoff(1e-12),
new_SGE(true),
computeJacobian(false),
computeJacobianExo(false),
computeHessian(false),
@ -67,7 +67,7 @@ ModelTree::derive(int order)
int eq = it->first.first;
int var1 = it->first.second;
NodeID d1 = it->second;
// Store only second derivatives with var2 <= var1
for(int var2 = 0; var2 <= var1; var2++)
{
@ -220,6 +220,7 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
eq_node = equations[ModelBlock->Block_List[j].Equation[0]];
lhs = eq_node->arg1;
rhs = eq_node->arg2;
tmp_s.str("");
tmp_output.str("");
lhs->writeOutput(tmp_output, oCDynamicModelSparseDLL, temporary_terms);
tmp_s << "y[Per_y_+" << ModelBlock->Block_List[j].Variable[0] << "]";
@ -230,6 +231,18 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
else if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_FOREWARD_SIMPLE)
ModelBlock->Block_List[j].Simulation_Type=EVALUATE_FOREWARD;
}
else
{
tmp_output.str("");
rhs->writeOutput(tmp_output, oCDynamicModelSparseDLL, temporary_terms);
if (tmp_output.str()==tmp_s.str())
{
if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_BACKWARD_SIMPLE)
ModelBlock->Block_List[j].Simulation_Type=EVALUATE_BACKWARD_R;
else if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_FOREWARD_SIMPLE)
ModelBlock->Block_List[j].Simulation_Type=EVALUATE_FOREWARD_R;
}
}
}
for(i = 0;i < ModelBlock->Block_List[j].Size;i++)
{
@ -237,7 +250,9 @@ ModelTree::computeTemporaryTermsOrdered(int order, Model_Block *ModelBlock)
eq_node->computeTemporaryTerms(reference_count, temporary_terms, first_occurence, j, ModelBlock);
}
if (ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FOREWARD)
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FOREWARD
&&ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD_R
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FOREWARD_R)
{
if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_COMPLETE ||
ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_SIMPLE)
@ -295,7 +310,7 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
int prev_Simulation_Type=-1;
temporary_terms_type::const_iterator it_temp=temporary_terms.begin();
//----------------------------------------------------------------------
//Temporary variables décalaration
//Temporary variables declaration
OK=true;
for(temporary_terms_type::const_iterator it = temporary_terms.begin();
it != temporary_terms.end(); it++)
@ -330,7 +345,9 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
lhs_rhs_done=false;
if (prev_Simulation_Type==ModelBlock->Block_List[j].Simulation_Type
&& (ModelBlock->Block_List[j].Simulation_Type==EVALUATE_BACKWARD
||ModelBlock->Block_List[j].Simulation_Type==EVALUATE_FOREWARD ))
||ModelBlock->Block_List[j].Simulation_Type==EVALUATE_FOREWARD
||ModelBlock->Block_List[j].Simulation_Type==EVALUATE_BACKWARD_R
||ModelBlock->Block_List[j].Simulation_Type==EVALUATE_FOREWARD_R ))
skip_the_head=true;
else
skip_the_head=false;
@ -346,6 +363,10 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
" // Simulation type ";
output << BlockTriangular::BlockSim(ModelBlock->Block_List[j].Simulation_Type) << " //\n" <<
" ////////////////////////////////////////////////////////////////////////\n";
#ifdef CONDITION
if(ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_COMPLETE)
output << " longd condition[" << ModelBlock->Block_List[j].Size << "]; /*to improve condition*/\n";
#endif
}
//The Temporary terms
temporary_terms_type tt2;
@ -389,6 +410,13 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
rhs->writeOutput(output, oCDynamicModelSparseDLL, temporary_terms);
output << ";\n";
break;
case EVALUATE_BACKWARD_R:
case EVALUATE_FOREWARD_R:
rhs->writeOutput(output, oCDynamicModelSparseDLL, temporary_terms);
output << " = ";
lhs->writeOutput(output, oCDynamicModelSparseDLL, temporary_terms);
output << ";\n";
break;
case SOLVE_BACKWARD_COMPLETE:
case SOLVE_FOREWARD_COMPLETE:
Uf[ModelBlock->Block_List[j].Equation[i]] << " u[" << i << "] = residual[" << i << "]";
@ -403,11 +431,17 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
output << ") - (";
rhs->writeOutput(output, oCDynamicModelSparseDLL, temporary_terms);
output << ");\n";
#ifdef CONDITION
if(ModelBlock->Block_List[j].Simulation_Type==SOLVE_TWO_BOUNDARIES_COMPLETE)
output << " condition[" << i << "]=0;\n";
#endif
}
}
// The Jacobian if we have to solve the block
if (ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FOREWARD)
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FOREWARD
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_BACKWARD_R
&& ModelBlock->Block_List[j].Simulation_Type!=EVALUATE_FOREWARD_R)
{
output << " /* Jacobian */\n";
switch(ModelBlock->Block_List[j].Simulation_Type)
@ -427,11 +461,11 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
{
int eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_Index[i];
int var=ModelBlock->Block_List[j].IM_lead_lag[m].Var_Index[i];
int u=ModelBlock->Block_List[j].IM_lead_lag[m].u[i];
int u=ModelBlock->Block_List[j].IM_lead_lag[m].us[i];
int eqr=ModelBlock->Block_List[j].IM_lead_lag[m].Equ[i];
int varr=ModelBlock->Block_List[j].IM_lead_lag[m].Var[i];
Uf[ModelBlock->Block_List[j].Equation[eqr]] << "-u[" << u << "]*y[Per_y_+" << var << "]";
output << " u[" << u << "] = g1[" << eqr << "*" << ModelBlock->Block_List[j].Size << "+" << varr << "] = ";
output << " u[" << u << "] = "/*g1[" << eqr << "*" << ModelBlock->Block_List[j].Size << "+" << varr << "] = "*/;
writeDerivative(output, eq, var, 0, oCDynamicModelSparseDLL, temporary_terms);
output << "; // variable=" << symbol_table.getNameByID(eEndogenous, var)
<<"(" << variable_table.getLag(variable_table.getSymbolID(var))<< ") " << var
@ -461,10 +495,36 @@ ModelTree::writeModelEquationsOrdered(ostream &output, Model_Block *ModelBlock)
output << "; // variable=" << symbol_table.getNameByID(eEndogenous, var)
<<"(" << k << ") " << var
<< ", equation=" << eq << "\n";
#ifdef CONDITION
output << " if(fabs(condition[" << eqr << "])<fabs(u[" << u << "+Per_u_]))\n";
output << " condition[" << eqr << "]=u[" << u << "+Per_u_];\n";
#endif
}
}
for(i = 0;i < ModelBlock->Block_List[j].Size;i++)
output << Uf[ModelBlock->Block_List[j].Equation[i]].str() << ";\n";
{
output << Uf[ModelBlock->Block_List[j].Equation[i]].str() << ";\n";
#ifdef CONDITION
output << " if(fabs(condition[" << i << "])<fabs(u[" << i << "+Per_u_]))\n";
output << " condition[" << i << "]=u[" << i << "+Per_u_];\n";
#endif
}
#ifdef CONDITION
for(m=0;m<=ModelBlock->Block_List[j].Max_Lead+ModelBlock->Block_List[j].Max_Lag;m++)
{
k=m-ModelBlock->Block_List[j].Max_Lag;
for(i=0;i<ModelBlock->Block_List[j].IM_lead_lag[m].size;i++)
{
int eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_Index[i];
int var=ModelBlock->Block_List[j].IM_lead_lag[m].Var_Index[i];
int u=ModelBlock->Block_List[j].IM_lead_lag[m].u[i];
int eqr=ModelBlock->Block_List[j].IM_lead_lag[m].Equ[i];
output << " u[" << u << "+Per_u_] /= condition[" << eqr << "];\n";
}
}
for(i = 0;i < ModelBlock->Block_List[j].Size;i++)
output << " u[" << i << "+Per_u_] /= condition[" << i << "];\n";
#endif
break;
}
}
@ -686,11 +746,10 @@ ModelTree::writeDynamicCFile(const string &dynamic_basename) const
<< " }" << endl
<< " params = mxGetPr(mxGetFieldByNumber(M_, 0, mxGetFieldNumber(M_,\"params\")));" << endl
<< " /* Gets it_ from global workspace of Matlab */" << endl
<< " it_ = (int) floor(mxGetScalar(mexGetVariable(\"global\", \"it_\")))-1;" << endl
<< " //it_ = (int) floor(mxGetScalar(mexGetVariable(\"global\", \"it_\")))-1;" << endl
<< " /* Call the C subroutines. */" << endl
<< " Dynamic(y, x, residual, g1, g2);" << endl
<< "}" << endl;
mDynamicModelFile.close();
}
@ -723,7 +782,7 @@ ModelTree::writeStaticModel(ostream &StaticOutput) const
ostringstream g1;
g1 << " g1" << LPAR(output_type) << eq + 1 << ", "
<< variable_table.getSymbolID(var) + 1 << RPAR(output_type);
jacobian_output << g1.str() << "=" << g1.str() << "+";
d1->writeOutput(jacobian_output, output_type, temporary_terms);
jacobian_output << ";" << endl;
@ -947,13 +1006,81 @@ ModelTree::writeSparseDLLDynamicHFile(const string &dynamic_basename) const
mDynamicModelFile << "#endif\n";
mDynamicModelFile.close();
}
void
ModelTree::Write_Inf_To_Bin_File(const string &dynamic_basename, const string &bin_basename, const int &num,
int &u_count_int, bool &file_open) const
{
int j;
std::ofstream SaveCode;
/*cout << "bin_basename=" << bin_basename << "\n";
system("pause");*/
if(file_open)
SaveCode.open((bin_basename + ".bin").c_str(), ios::out | ios::in | ios::binary | ios ::ate );
else
SaveCode.open((bin_basename + ".bin").c_str(), ios::out | ios::binary);
if(!SaveCode.is_open())
{
cout << "Error : Can't open file \"" << bin_basename << ".bin\" for writing\n";
exit( -1);
}
u_count_int=0;
for(int m=0;m<=block_triangular.ModelBlock->Block_List[num].Max_Lead+block_triangular.ModelBlock->Block_List[num].Max_Lag;m++)
{
int k1=m-block_triangular.ModelBlock->Block_List[num].Max_Lag;
//mDynamicModelFile << " Lead_Lag.push_back(k1);\n";
for(j=0;j<block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].size;j++)
{
//int eq=block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].Equ_Index[j];
//int var=block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].Var_Index[j];
//int eqr=block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].Equ[j];
int varr=block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].Var[j]+k1*block_triangular.ModelBlock->Block_List[num].Size;
int u=block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].u[j];
int eqr1=block_triangular.ModelBlock->Block_List[num].IM_lead_lag[m].Equ[j];
/*cout << " ! IM_i[std::make_pair(std::make_pair(" << eqr1 << ", " << varr+k1*block_triangular.ModelBlock->Block_List[num].Size << "), " << k1 << ")] = " << u << ";\n";
cout << " ? IM_i[std::make_pair(std::make_pair(" << eqr1 << ", " << varr << "), " << k1 << ")] = " << u << ";\n";*/
SaveCode.write(reinterpret_cast<char *>(&eqr1), sizeof(eqr1));
SaveCode.write(reinterpret_cast<char *>(&varr), sizeof(varr));
SaveCode.write(reinterpret_cast<char *>(&k1), sizeof(k1));
SaveCode.write(reinterpret_cast<char *>(&u), sizeof(u));
u_count_int++;
}
}
for(j=0;j<block_triangular.ModelBlock->Block_List[num].Size;j++)
{
//int eq=block_triangular.ModelBlock->Block_List[i].Equation[j];
int eqr1=j;
int varr=block_triangular.ModelBlock->Block_List[num].Size*(block_triangular.periods
+/*block_triangular.ModelBlock->Block_List[num].Max_Lead*/block_triangular.Model_Max_Lead);
int k1=0;
//mDynamicModelFile << " var_in_equ_and_lag[std::make_pair(std::make_pair(" << j << ", 0)] = -1;\n";
//mDynamicModelFile << " /*periods=" << block_triangular.periods << " Size=" << block_triangular.ModelBlock->Block_List[i].Size << "*/\n";
//mDynamicModelFile << " var_in_equ_and_lag.insert(std::make_pair(std::make_pair(" << eqr1 << ", 0), " << block_triangular.ModelBlock->Block_List[i].Size*block_triangular.periods << "));\n";
//mDynamicModelFile << " equ_in_var_and_lag[std::make_pair(std::make_pair(-1, " << k1 << ")] = " << equ << ";\n";
SaveCode.write(reinterpret_cast<char *>(&eqr1), sizeof(eqr1));
SaveCode.write(reinterpret_cast<char *>(&varr), sizeof(varr));
SaveCode.write(reinterpret_cast<char *>(&k1), sizeof(k1));
SaveCode.write(reinterpret_cast<char *>(&eqr1), sizeof(eqr1));
//cout << " IM_i[std::make_pair(std::make_pair(" << eqr1 << ", " << varr << "), " << k1 << ")] = " << eqr1 << ";\n";
u_count_int++;
}
for(j=0;j<block_triangular.ModelBlock->Block_List[num].Size;j++)
{
//mDynamicModelFile << " index_var[" << j << "]=" << block_triangular.ModelBlock->Block_List[i].Variable[j] << ";\n";
int varr=block_triangular.ModelBlock->Block_List[num].Variable[j];
SaveCode.write(reinterpret_cast<char *>(&varr), sizeof(varr));
}
SaveCode.close();
}
void
ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename, const string &bin_basename) const
{
string filename;
ofstream mDynamicModelFile;
SymbolicGaussElimination SGE;
if (compiler == LCC_COMPILE)
filename = dynamic_basename + ".c";
else
@ -995,23 +1122,22 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
int i, j, k, Nb_SGE=0;
bool printed = false, skip_head, open_par=false;
SymbolicGaussElimination SGE;
if (computeJacobian || computeJacobianExo || computeHessian)
{
mDynamicModelFile << "void Dynamic_Init(tModel_Block *Model_Block)\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " int i;\n";
//mDynamicModelFile << "void Dynamic_Init(tModel_Block *Model_Block)\n";
mDynamicModelFile << "void Dynamic_Init()\n";
mDynamicModelFile << " {\n";
int prev_Simulation_Type=-1;
for(i = 0;i < block_triangular.ModelBlock->Size;i++)
{
k = block_triangular.ModelBlock->Block_List[i].Simulation_Type;
if (prev_Simulation_Type==k &&
(k==EVALUATE_FOREWARD || k==EVALUATE_BACKWARD))
(k==EVALUATE_FOREWARD || k==EVALUATE_BACKWARD || k==EVALUATE_FOREWARD_R || k==EVALUATE_BACKWARD_R))
skip_head=true;
else
skip_head=false;
if ((k == EVALUATE_FOREWARD) && (block_triangular.ModelBlock->Block_List[i].Size))
if ((k == EVALUATE_FOREWARD || k == EVALUATE_FOREWARD_R) && (block_triangular.ModelBlock->Block_List[i].Size))
{
if (!skip_head)
{
@ -1027,10 +1153,10 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << "#ifdef DEBUG\n";
}
for(j = 0;j < block_triangular.ModelBlock->Block_List[i].Size;j++)
mDynamicModelFile << " mexPrintf(\"y[%d, %d]=%f \\n\",it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << ",y[it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << "]);\n";
mDynamicModelFile << " mexPrintf(\"y[%d, %d]=%f \\n\",it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << ",double(y[it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << "]));\n";
open_par=true;
}
else if ((k == EVALUATE_BACKWARD) && (block_triangular.ModelBlock->Block_List[i].Size))
else if ((k == EVALUATE_BACKWARD || k == EVALUATE_BACKWARD_R) && (block_triangular.ModelBlock->Block_List[i].Size))
{
if (!skip_head)
{
@ -1046,7 +1172,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << "#ifdef DEBUG\n";
}
for(j = 0;j < block_triangular.ModelBlock->Block_List[i].Size;j++)
mDynamicModelFile << " mexPrintf(\"y[%d, %d]=%f \\n\",it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << ",y[it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << "]);\n";
mDynamicModelFile << " mexPrintf(\"y[%d, %d]=%f \\n\",it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << ",double(y[it_," << block_triangular.ModelBlock->Block_List[i].Variable[j] << "]));\n";
open_par=true;
}
else if ((k == SOLVE_FOREWARD_SIMPLE) && (block_triangular.ModelBlock->Block_List[i].Size))
@ -1129,7 +1255,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
{
printed = true;
}
SGE.SGE_compute(block_triangular.ModelBlock, i, true, bin_basename, /*mod_param.endo_nbr*/symbol_table.endo_nbr);
SGE.SGE_compute(block_triangular.ModelBlock, i, true, bin_basename, symbol_table.endo_nbr);
Nb_SGE++;
#ifdef PRINT_OUT
cout << "end of Gaussian elimination\n";
@ -1163,7 +1289,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " }\n";
mDynamicModelFile << " iter++;\n";
mDynamicModelFile << " cvg=(max_res<solve_tolf);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", periods, y_kmin, y_kmax);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", periods, true);\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " if (!cvg)\n";
mDynamicModelFile << " {\n";
@ -1186,9 +1312,9 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " mexPrintf(\"\\n\");\n";
mDynamicModelFile << "#endif\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", periods, true);\n";
}
mDynamicModelFile << "#ifdef DEBUG\n";
/*mDynamicModelFile << "#ifdef DEBUG\n";
mDynamicModelFile << " for(it_=y_kmin;it_<periods+y_kmin;it_++)\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " for(i=0;i<Model_Block->List[" << i << "].Size;i++)\n";
@ -1198,7 +1324,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " }";
mDynamicModelFile << " mexPrintf(\" \\n \");\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << "#endif\n";
mDynamicModelFile << "#endif\n";*/
mDynamicModelFile << " mxFree(g1);\n";
mDynamicModelFile << " mxFree(r);\n";
mDynamicModelFile << " mxFree(u);\n";
@ -1218,7 +1344,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
}
SGE.SGE_compute(block_triangular.ModelBlock, i, false, bin_basename, /*mod_param.endo_nbr*/symbol_table.endo_nbr);
Nb_SGE++;
mDynamicModelFile << " Read_file(\"" << reform(bin_basename) << "\", periods, 0, 0, 0, 0);\n";
mDynamicModelFile << " Read_file(\"" << reform(bin_basename) << "\", periods, 0, " << symbol_table.endo_nbr << ", " << block_triangular.ModelBlock->Block_List[i].Max_Lag << ", " << block_triangular.ModelBlock->Block_List[i].Max_Lead << " );\n";
mDynamicModelFile << " g1=(double*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size*block_triangular.ModelBlock->Block_List[i].Size << "*sizeof(double));\n";
mDynamicModelFile << " r=(double*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size << "*sizeof(double));\n";
mDynamicModelFile << " for(it_=y_kmin;it_<periods+y_kmin;it_++)\n";
@ -1231,7 +1357,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " while(!((cvg)||(iter>maxit_)))\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " Dynamic" << i + 1 << "(y, x, r, g1, g2);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", 0, false);\n";
mDynamicModelFile << " res2=0;\n";
mDynamicModelFile << " res1=0;\n";
mDynamicModelFile << " max_res=0;\n";
@ -1254,13 +1380,13 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
else
{
mDynamicModelFile << " Dynamic" << i + 1 << "(y, x, r, g1, g2);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", 0, false);\n";
}
mDynamicModelFile << "#ifdef DEBUG\n";
/*mDynamicModelFile << "#ifdef DEBUG\n";
mDynamicModelFile << " for(i=0;i<Model_Block->List[" << i << "].Size;i++)\n";
mDynamicModelFile << " mexPrintf(\" y[%d, %d]=%f \",it_,Model_Block->List[" << i << "].Variable[i],y[it_*" << /*mod_param.endo_nbr*/symbol_table.endo_nbr << "+Model_Block->List[" << i << "].Variable[i]]);\n";
mDynamicModelFile << " mexPrintf(\" y[%d, %d]=%f \",it_,Model_Block->List[" << i << "].Variable[i],y[it_*" << symbol_table.endo_nbr << "+Model_Block->List[" << i << "].Variable[i]]);\n";
mDynamicModelFile << " mexPrintf(\" \\n \");\n";
mDynamicModelFile << "#endif\n";
mDynamicModelFile << "#endif\n";*/
mDynamicModelFile << " }\n";
mDynamicModelFile << " mxFree(g1);\n";
mDynamicModelFile << " mxFree(r);\n";
@ -1276,7 +1402,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
open_par=false;
SGE.SGE_compute(block_triangular.ModelBlock, i, false, bin_basename, /*mod_param.endo_nbr*/symbol_table.endo_nbr);
Nb_SGE++;
mDynamicModelFile << " Read_file(\"" << reform(bin_basename) << "\", periods, 0, 0, 0, 0);\n";
mDynamicModelFile << " Read_file(\"" << reform(bin_basename) << "\", periods, 0, " << symbol_table.endo_nbr << ", " << block_triangular.ModelBlock->Block_List[i].Max_Lag << ", " << block_triangular.ModelBlock->Block_List[i].Max_Lead << " );\n";
mDynamicModelFile << " g1=(double*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size*block_triangular.ModelBlock->Block_List[i].Size << "*sizeof(double));\n";
mDynamicModelFile << " r=(double*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size << "*sizeof(double));\n";
mDynamicModelFile << " for(it_=periods+y_kmin;it_>y_kmin;it_--)\n";
@ -1289,7 +1415,7 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " while(!((cvg)||(iter>maxit_)))\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " Dynamic" << i + 1 << "(y, x, r, g1, g2);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", 0, false);\n";
mDynamicModelFile << " res2=0;\n";
mDynamicModelFile << " for(i=0;i<" << block_triangular.ModelBlock->Block_List[i].Size << ";i++)\n";
mDynamicModelFile << " res2+=r[i]*r[i];\n";
@ -1305,13 +1431,13 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
else
{
mDynamicModelFile << " Dynamic" << i + 1 << "(y, x, r, g1, g2);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", 0, false);\n";
}
mDynamicModelFile << "#ifdef DEBUG\n";
/*mDynamicModelFile << "#ifdef DEBUG\n";
mDynamicModelFile << " for(i=0;i<Model_Block->List[" << i << "].Size;i++)\n";
mDynamicModelFile << " mexPrintf(\" y[%d, %d]=%f \",it_,Model_Block->List[" << i << "].Variable[i],y[it_*" << /*mod_param.endo_nbr*/symbol_table.endo_nbr << "+Model_Block->List[" << i << "].Variable[i]]);\n";
mDynamicModelFile << " mexPrintf(\" y[%d, %d]=%f \",it_,Model_Block->List[" << i << "].Variable[i],y[it_*" << symbol_table.endo_nbr << "+Model_Block->List[" << i << "].Variable[i]]);\n";
mDynamicModelFile << " mexPrintf(\" \\n \");\n";
mDynamicModelFile << "#endif\n";
mDynamicModelFile << "#endif\n";*/
mDynamicModelFile << " }\n";
mDynamicModelFile << " mxFree(g1);\n";
mDynamicModelFile << " mxFree(r);\n";
@ -1329,11 +1455,33 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
{
printed = true;
}
SGE.SGE_compute(block_triangular.ModelBlock, i, true, bin_basename, /*mod_param.endo_nbr*/symbol_table.endo_nbr);
Nb_SGE++;
mDynamicModelFile << " Read_file(\"" << reform(bin_basename) << "\",periods," <<
block_triangular.ModelBlock->Block_List[i].IM_lead_lag[block_triangular.ModelBlock->Block_List[i].Max_Lag + block_triangular.ModelBlock->Block_List[i].Max_Lead].u_finish + 1 << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr <<
", " << block_triangular.ModelBlock->Block_List[i].Max_Lag << ", " << block_triangular.ModelBlock->Block_List[i].Max_Lead << ");\n";
//cout << "new_SGE=" << new_SGE << "\n";
if(new_SGE)
{
int u_count_int=0;
Write_Inf_To_Bin_File(dynamic_basename, bin_basename, i, u_count_int,SGE.file_open);
SGE.file_is_open();
mDynamicModelFile << " u_count=" << u_count_int << "*periods;\n";
mDynamicModelFile << " u_count_alloc = 2*u_count;\n";
mDynamicModelFile << " u=(longd*)mxMalloc(u_count_alloc*sizeof(longd));\n";
mDynamicModelFile << " memset(u, 0, u_count_alloc*sizeof(longd));\n";
mDynamicModelFile << " u_count_init=" <<
block_triangular.ModelBlock->Block_List[i].IM_lead_lag[block_triangular.ModelBlock->Block_List[i].Max_Lag +
block_triangular.ModelBlock->Block_List[i].Max_Lead].u_finish + 1 << ";\n";
//mDynamicModelFile << " index_var=(int*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size<< "*sizeof(*index_var));\n";
mDynamicModelFile << " Read_SparseMatrix(\"" << reform(bin_basename) << "\","
<< block_triangular.ModelBlock->Block_List[i].Size << ", periods, y_kmin, y_kmax"
<< ");\n";
mDynamicModelFile << " u_count=" << u_count_int << "*(periods+y_kmax);\n";
}
else
{
SGE.SGE_compute(block_triangular.ModelBlock, i, true, bin_basename, /*mod_param.endo_nbr*/symbol_table.endo_nbr);
mDynamicModelFile << " Read_file(\"" << reform(bin_basename) << "\",periods," <<
block_triangular.ModelBlock->Block_List[i].IM_lead_lag[block_triangular.ModelBlock->Block_List[i].Max_Lag + block_triangular.ModelBlock->Block_List[i].Max_Lead].u_finish + 1 << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr <<
", " << block_triangular.ModelBlock->Block_List[i].Max_Lag << ", " << block_triangular.ModelBlock->Block_List[i].Max_Lead << ");\n";
}
mDynamicModelFile << " g1=(double*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size*block_triangular.ModelBlock->Block_List[i].Size << "*sizeof(double));\n";
mDynamicModelFile << " r=(double*)mxMalloc(" << block_triangular.ModelBlock->Block_List[i].Size << "*sizeof(double));\n";
if (!block_triangular.ModelBlock->Block_List[i].is_linear)
@ -1350,6 +1498,8 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " Per_u_=(it_-y_kmin)*" << block_triangular.ModelBlock->Block_List[i].IM_lead_lag[block_triangular.ModelBlock->Block_List[i].Max_Lag + block_triangular.ModelBlock->Block_List[i].Max_Lead].u_finish + 1 << ";\n";
mDynamicModelFile << " Per_y_=it_*y_size;\n";
mDynamicModelFile << " Dynamic" << i + 1 << "(y, x, r, g1, g2);\n";
mDynamicModelFile << " if(isnan(res1)||isinf(res1))\n";
mDynamicModelFile << " break;\n";
mDynamicModelFile << " for(i=0;i<" << block_triangular.ModelBlock->Block_List[i].Size << ";i++)\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " if (max_res<fabs(r[i]))\n";
@ -1359,7 +1509,11 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " }\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " cvg=(max_res<solve_tolf);\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", periods, y_kmin, y_kmax);\n";
mDynamicModelFile << " if(!cvg)\n";
if(new_SGE)
mDynamicModelFile << " simulate_NG1(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", periods, true);\n";
else
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", periods, true);\n";
mDynamicModelFile << " iter++;\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " if (!cvg)\n";
@ -1383,22 +1537,27 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " mexPrintf(\"\\n\");\n";
mDynamicModelFile << "#endif\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax);\n";
if(new_SGE)
mDynamicModelFile << " simulate_NG1(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", periods, true);\n";
else
mDynamicModelFile << " simulate(" << i << ", " << /*mod_param.endo_nbr*/symbol_table.endo_nbr << ", it_, y_kmin, y_kmax," << block_triangular.ModelBlock->Block_List[i].Size << ", periods, true);\n";
}
mDynamicModelFile << "#ifdef DEBUG\n";
/*mDynamicModelFile << "#ifdef DEBUG\n";
mDynamicModelFile << " for(it_=y_kmin;it_<periods+y_kmin;it_++)\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " for(i=0;i<Model_Block->List[" << i << "].Size;i++)\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " Per_y_=it_*y_size;\n";
mDynamicModelFile << " mexPrintf(\" y[%d, %d]=%f \",it_,Model_Block->List[" << i << "].Variable[i],y[it_*" << /*mod_param.endo_nbr*/symbol_table.endo_nbr << "+Model_Block->List[" << i << "].Variable[i]]);\n";
mDynamicModelFile << " mexPrintf(\" y[%d, %d]=%f \",it_,Model_Block->List[" << i << "].Variable[i],y[it_*" << symbol_table.endo_nbr << "+Model_Block->List[" << i << "].Variable[i]]);\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " mexPrintf(\" \\n \");\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << "#endif\n";
mDynamicModelFile << "#endif\n";*/
mDynamicModelFile << " mxFree(g1);\n";
mDynamicModelFile << " mxFree(r);\n";
mDynamicModelFile << " mxFree(u);\n";
mDynamicModelFile << " mxFree(index_vara);\n";
mDynamicModelFile << " memset(direction,0,size_of_direction);\n";
mDynamicModelFile << " //mexErrMsgTxt(\"Exit from Dynare\");\n";
}
prev_Simulation_Type=k;
@ -1409,19 +1568,18 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " }\n";
}
mDynamicModelFile << " }\n";
// Writing the gateway routine
mDynamicModelFile << " int max(int a, int b)\n";
/*mDynamicModelFile << " int max(int a, int b)\n";
mDynamicModelFile << " {\n";
mDynamicModelFile << " if (a>b) return(a); else return(b);\n";
mDynamicModelFile << " }\n\n\n";
mDynamicModelFile << " }\n\n\n";*/
mDynamicModelFile << "/* The gateway routine */\n";
mDynamicModelFile << "void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])\n";
mDynamicModelFile << "{\n";
mDynamicModelFile << " tModel_Block *Model_Block;\n";
/*mDynamicModelFile << " tModel_Block *Model_Block;\n";*/
mDynamicModelFile << " mxArray *M_, *oo_, *options_;\n";
mDynamicModelFile << " int i, j, row_y, col_y, row_x, col_x, x_FieldNumber;\n";
mDynamicModelFile << " mxArray *x_FieldByNumber;\n";
mDynamicModelFile << " int i, row_y, col_y, row_x, col_x;\n";
mDynamicModelFile << " double * pind ;\n";
mDynamicModelFile << "\n";
mDynamicModelFile << " /* Gets model parameters from global workspace of Matlab */\n";
@ -1445,9 +1603,10 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " mexErrMsgTxt(\"options_ \\n\");\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " params = mxGetPr(mxGetFieldByNumber(M_, 0, mxGetFieldNumber(M_,\"params\")));\n";
mDynamicModelFile << " y= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"endo_simul\")));\n";
mDynamicModelFile << " double *yd, *xd;\n";
mDynamicModelFile << " yd= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"endo_simul\")));\n";
mDynamicModelFile << " row_y=mxGetM(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"endo_simul\")));\n";
mDynamicModelFile << " x= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"exo_simul\")));\n";
mDynamicModelFile << " xd= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"exo_simul\")));\n";
mDynamicModelFile << " row_x=mxGetM(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"exo_simul\")));\n";
mDynamicModelFile << " col_x=mxGetN(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,\"exo_simul\")));\n";
if (compiler==GCC_COMPILE)
@ -1474,8 +1633,21 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " col_y=row_x;\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " solve_tolf=*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,\"dynatol\"))));\n";
mDynamicModelFile << " size_of_direction=col_y*row_y*sizeof(longd);\n";
mDynamicModelFile << " y=(longd*)mxMalloc(size_of_direction);\n";
mDynamicModelFile << " ya=(longd*)mxMalloc(size_of_direction);\n";
mDynamicModelFile << " direction=(longd*)mxMalloc(size_of_direction);\n";
mDynamicModelFile << " memset(direction,0,size_of_direction);\n";
mDynamicModelFile << " x=(longd*)mxMalloc(col_x*row_x*sizeof(longd));\n";
mDynamicModelFile << " for(i=0;i<row_x*col_x;i++)\n";
mDynamicModelFile << " x[i]=longd(xd[i]);\n";
mDynamicModelFile << " for(i=0;i<row_y*col_y;i++)\n";
mDynamicModelFile << " y[i]=longd(yd[i]);\n";
mDynamicModelFile << " \n";
mDynamicModelFile << " y_size=row_y;\n";
mDynamicModelFile << " x_size=periods+y_kmin+y_kmax;\n";
mDynamicModelFile << " x_size=row_x;\n";
mDynamicModelFile << " nb_row_x=row_x;\n";
mDynamicModelFile << "#ifdef DEBUG\n";
mDynamicModelFile << " for(j=0;j<periods+y_kmin+y_kmax;j++)\n";
mDynamicModelFile << " {\n";
@ -1496,10 +1668,10 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << "#endif\n";
mDynamicModelFile << " /* Gets it_ from global workspace of Matlab */\n";
mDynamicModelFile << " it_ = (int) floor(mxGetScalar(mexGetVariable(\"global\", \"it_\")))-1;\n";
mDynamicModelFile << " //it_ = (int) floor(mxGetScalar(mexGetVariable(\"global\", \"it_\")))-1;\n";
mDynamicModelFile << " /* Call the C subroutines. */\n";
mDynamicModelFile << " t0= pctimer();\n";
mDynamicModelFile << " Dynamic_Init(Model_Block);\n";
mDynamicModelFile << " Dynamic_Init();\n";
mDynamicModelFile << " t1= pctimer();\n";
mDynamicModelFile << " mexPrintf(\"Simulation Time=%f milliseconds\\n\",1000*(t1-t0));\n";
if (compiler==GCC_COMPILE )
@ -1519,6 +1691,10 @@ ModelTree::writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename,
mDynamicModelFile << " for(i=0;i<row_y*col_y;i++)\n";
mDynamicModelFile << " pind[i]=y[i];\n";
mDynamicModelFile << " }\n";
mDynamicModelFile << " mxFree(x);\n";
mDynamicModelFile << " mxFree(y);\n";
mDynamicModelFile << " mxFree(ya);\n";
mDynamicModelFile << " mxFree(direction);\n";
mDynamicModelFile << "}\n";
}
mDynamicModelFile.close();
@ -1810,9 +1986,10 @@ ModelTree::checkPass() const
}
void
ModelTree::evaluateJacobian(const eval_context_type &eval_context)
ModelTree::evaluateJacobian(const eval_context_type &eval_context, jacob_map *j_m)
{
int i=0;
int j=0;
bool *IM;
int a_variable_lag=-9999;
for(first_derivatives_type::iterator it = first_derivatives.begin();
@ -1830,6 +2007,11 @@ ModelTree::evaluateJacobian(const eval_context_type &eval_context)
IM=block_triangular.bGet_IM(k1);
a_variable_lag=k1;
}
if(k1==0)
{
j++;
(*j_m)[make_pair(eq,var)]=val;
}
if (IM[eq*symbol_table.endo_nbr+var] && (fabs(val) < cutoff))
{
//cout << "the coefficient related to variable " << var << " with lag " << k1 << " in equation " << eq << " is equal to " << interprete_.u1 << " and is set to 0 in the incidence matrix\n";
@ -1839,33 +2021,43 @@ ModelTree::evaluateJacobian(const eval_context_type &eval_context)
}
}
if (i>0)
cout << i << " elements in the incidence matrices are below the cutoff (" << cutoff << ") and are discarded\n";
{
cout << i << " elements among " << first_derivatives.size() << " in the incidence matrices are below the cutoff (" << cutoff << ") and are discarded\n";
cout << "the contemporaneous incidence matrix has " << j << " elements\n";
}
}
void
ModelTree::BlockLinear(Model_Block *ModelBlock)
{
int i,j,l,m;
int i,j,l,m,ll;
for(j = 0;j < ModelBlock->Size;j++)
{
if (ModelBlock->Block_List[j].Simulation_Type==SOLVE_BACKWARD_COMPLETE ||
ModelBlock->Block_List[j].Simulation_Type==SOLVE_FOREWARD_COMPLETE)
{
for(i=0;i<ModelBlock->Block_List[j].IM_lead_lag[0].size;i++)
ll=ModelBlock->Block_List[j].Max_Lag;
for(i=0;i<ModelBlock->Block_List[j].IM_lead_lag[ll].size;i++)
{
int eq=ModelBlock->Block_List[j].IM_lead_lag[0].Equ_Index[i];
int var=ModelBlock->Block_List[j].IM_lead_lag[0].Var_Index[i];
int eq=ModelBlock->Block_List[j].IM_lead_lag[ll].Equ_Index[i];
int var=ModelBlock->Block_List[j].IM_lead_lag[ll].Var_Index[i];
first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,variable_table.getmVariableSelector(var,0)));
NodeID Id = it->second;
Id->collectEndogenous(Id);
if (Id->present_endogenous_size()>0)
if(it!= first_derivatives.end())
{
for(l=0;l<ModelBlock->Block_List[j].Size;l++)
NodeID Id = it->second;
/*cout << "i=" << i << " j=" << j << "\n";
cout << "eq=" << eq << " var=" << var << "\n";
cout << "Id=" << Id << "\n";*/
Id->collectEndogenous(Id);
if (Id->present_endogenous_size()>0)
{
if (Id->present_endogenous_find(ModelBlock->Block_List[j].Variable[l],0))
for(l=0;l<ModelBlock->Block_List[j].Size;l++)
{
ModelBlock->Block_List[j].is_linear=false;
goto follow;
if (Id->present_endogenous_find(ModelBlock->Block_List[j].Variable[l],0))
{
ModelBlock->Block_List[j].is_linear=false;
goto follow;
}
}
}
}
@ -1882,15 +2074,18 @@ ModelTree::BlockLinear(Model_Block *ModelBlock)
int var=ModelBlock->Block_List[j].IM_lead_lag[m].Var_Index[i];
first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,variable_table.getmVariableSelector(var,k1)));
NodeID Id = it->second;
Id->collectEndogenous(Id);
if (Id->present_endogenous_size()>0)
if(it!= first_derivatives.end())
{
for(l=0;l<ModelBlock->Block_List[j].Size;l++)
Id->collectEndogenous(Id);
if (Id->present_endogenous_size()>0)
{
if (Id->present_endogenous_find(ModelBlock->Block_List[j].Variable[l],k1))
for(l=0;l<ModelBlock->Block_List[j].Size;l++)
{
ModelBlock->Block_List[j].is_linear=false;
goto follow;
if (Id->present_endogenous_find(ModelBlock->Block_List[j].Variable[l],k1))
{
ModelBlock->Block_List[j].is_linear=false;
goto follow;
}
}
}
}
@ -1924,11 +2119,13 @@ ModelTree::computingPass(const eval_context_type &eval_context)
if (mode == eSparseDLLMode)
{
jacob_map j_m;
int Size;
int HSize;
int *Table=variable_table.GetVariableTable(&Size,&HSize);
evaluateJacobian(eval_context);
evaluateJacobian(eval_context, &j_m);
if (block_triangular.bt_verbose)
{
@ -1936,7 +2133,7 @@ ModelTree::computingPass(const eval_context_type &eval_context)
block_triangular.Print_IM( symbol_table.endo_nbr);
}
block_triangular.SetVariableTable(Table, Size, HSize);
block_triangular.Normalize_and_BlockDecompose_Static_0_Model();
block_triangular.Normalize_and_BlockDecompose_Static_0_Model(j_m);
BlockLinear(block_triangular.ModelBlock);
computeTemporaryTermsOrdered(order, block_triangular.ModelBlock);

View File

@ -203,10 +203,15 @@ ModelBlock_Graph(Model_Block *ModelBlock, int Blck_num, bool dynamic, t_model_gr
}
else
{
int sup;
Lead = ModelBlock->Block_List[Blck_num].Max_Lead;
Lag = ModelBlock->Block_List[Blck_num].Max_Lag;
int sup = Lead + Lag +3;
*periods = Lead + Lag + sup;
cout << "---> *periods=" << *periods << "\n";
if(*periods>3)
{
sup = Lead + Lag +3;
*periods = Lead + Lag + sup;
}
#ifdef PRINT_OUT
cout << "Lag=" << Lag << " Lead=" << Lead << "\n";
cout << "periods=Lead+2*Lag+2= " << *periods << "\n";

View File

@ -187,7 +187,7 @@ ParsingDriver::periods(string *periods)
void
ParsingDriver::cutoff(string *cutoff)
{
int cutoff_val = atoi(cutoff->c_str());
double cutoff_val = atof(cutoff->c_str());
mod_file->addStatement(new CutoffStatement(cutoff_val));
delete cutoff;
}
@ -225,7 +225,7 @@ ParsingDriver::init_param(string *name, NodeID rhs)
double val = rhs->eval(mod_file->global_eval_context);
int symb_id = mod_file->symbol_table.getID(*name);
mod_file->global_eval_context[make_pair(symb_id, eParameter)] = val;
}
}
catch(ExprNode::EvalException &e)
{
}
@ -427,7 +427,7 @@ ParsingDriver::add_covar_shock(string *var1, string *var2, NodeID value)
|| corr_shocks.find(key_inv) != corr_shocks.end())
error("shocks: covariance or correlation shock on variable pair (" + *var1 + ", "
+ *var2 + ") declared twice");
covar_shocks[key] = value;
delete var1;
@ -448,7 +448,7 @@ ParsingDriver::add_correl_shock(string *var1, string *var2, NodeID value)
|| corr_shocks.find(key_inv) != corr_shocks.end())
error("shocks: covariance or correlation shock on variable pair (" + *var1 + ", "
+ *var2 + ") declared twice");
corr_shocks[key] = value;
delete var1;

View File

@ -11,6 +11,7 @@
//#include "pctimer_h.hh"
#include "Model_Graph.hh"
#include "SymbolGaussElim.hh"
//#define SIMPLIFY
using namespace std;
int max_nb_table_y, max_nb_in_degree_edges=0, max_nb_out_degree_edges=0;
@ -297,10 +298,11 @@ SymbolicGaussElimination::write_to_file_table_u(t_table_u *save_table_u, t_table
}
void
SymbolicGaussElimination::write_to_file_table_u_b(t_table_u *save_table_u, t_table_u *last_table_u, int *nb_save_table_u)
SymbolicGaussElimination::write_to_file_table_u_b(t_table_u *save_table_u, t_table_u *last_table_u, int *nb_save_table_u, bool chk)
{
t_table_u *table_u;
while(save_table_u!=last_table_u)
bool OK=true;
while(/*save_table_u!=last_table_u*/OK)
{
#ifdef PRINT_OUT
cout << "**save_table_u->type=" << int(save_table_u->type) << "\n";
@ -349,9 +351,23 @@ SymbolicGaussElimination::write_to_file_table_u_b(t_table_u *save_table_u, t_tab
#endif /**PRINT_OUT**/
break;
}
table_u = save_table_u->pNext;
free(save_table_u);
save_table_u = table_u;
if(chk)
{
OK=(save_table_u!=last_table_u);
if(OK)
{
table_u = save_table_u->pNext;
free(save_table_u);
save_table_u = table_u;
}
}
else
{
table_u = save_table_u->pNext;
free(save_table_u);
save_table_u = table_u;
OK=(save_table_u!=last_table_u);
}
}
}
@ -525,7 +541,7 @@ SymbolicGaussElimination::interpolation(t_model_graph* model_graph, t_table_y* t
t_table_u *tmp_table_u, *old_table_u;
t_table_u *c_first_table_u, *c_second_table_u, *c_third_table_u;
int c_first_y_blck, c_second_y_blck;
int i, j, k, up_to, op_count, k1, k2;
int i, j, k, up_to=0, op_count, k1, k2;
bool OK, modify_u_count;
int cur_pos, nb_table_u=0;
#ifdef PRINT_OUT
@ -701,6 +717,7 @@ SymbolicGaussElimination::interpolation(t_model_graph* model_graph, t_table_y* t
}
op_count++;
}
SaveCode.flush();
k=SaveCode.tellp();
SaveCode.seekp(cur_pos);
SaveCode.write(reinterpret_cast<char *>(&nb_table_u), sizeof(nb_table_u));
@ -974,6 +991,7 @@ SymbolicGaussElimination::Vertex_Elimination(t_model_graph* model_graph, int pos
cout << "--------------------------------------------------------------- \n";
cout << "Elimination of vertex " << lvertex[vertex_to_eliminate].index << " interpolate=" << *interpolate << "\n";
cout << "min_edge=" << min_edge << " length_markowitz=" << length_markowitz << "\n";
print_Graph(model_graph);
#endif /**PRINT_OUT**/
#ifdef DIRECT_COMPUTE
table_y[vertex_count].index = lvertex[vertex_to_eliminate].index;
@ -1087,7 +1105,7 @@ SymbolicGaussElimination::Vertex_Elimination(t_model_graph* model_graph, int pos
#ifdef DEBUGR
cout << " creation of edge between vertices " << lvertex[i1].index << " and " << lvertex[j1].index << "\n";
#endif /**DEBUGR**/
#ifdef SYMPLIFY
#ifdef SYMPLIFY /*no reuse of free numbers, to be sure that there is no cyclcial behaviour in numbering*/
curr_u_count = get_free_u_list(dynamic);
#else
curr_u_count = u_count;
@ -1441,13 +1459,13 @@ SymbolicGaussElimination::Gaussian_Elimination(t_model_graph* model_graph
cout << "going to open file file_open=" << file_open << " file_name={" << file_name << "}\n";
#endif
if(file_open)
SaveCode.open((file_name + ".bin").c_str(), ios::app | ios::binary);
SaveCode.open((file_name + ".bin").c_str(), ios::out | ios::in | ios::binary | ios ::ate );
else
SaveCode.open((file_name + ".bin").c_str(), ios::out | ios::binary);
file_open = true;
if(!SaveCode.is_open())
{
cout << "Error : Can't open file \"CurrentModel.bin\" for writing\n";
cout << "Error : Can't open file \"" << file_name << ".bin\" for writing\n";
exit( -1);
}
#ifdef PRINT_OUT
@ -1728,14 +1746,14 @@ SymbolicGaussElimination::Gaussian_Elimination(t_model_graph* model_graph
prologue_save_table_y[i-first_nb_prologue_save_table_y]=table_y[i];
k=SaveCode.tellp();
SaveCode.seekp(cur_pos);
i=nb_second_u_blck-prologue_nb_table_u;
i=nb_second_u_blck-prologue_nb_table_u+1;
//cout << "nb_prologue_save_table_u=" << i << "\n";
SaveCode.write(reinterpret_cast<char *>(&i), sizeof(i));
SaveCode.seekp(k);
write_to_file_table_u_b(first_u_blck,second_u_blck, &nb_prologue_save_table_u);
write_to_file_table_u_b(first_u_blck,second_u_blck, &nb_prologue_save_table_u,true);
//cout << "nb_prologue_save_table_u(1)=" << nb_prologue_save_table_u << "\n";
nb_first_u_blck=nb_second_u_blck;
nb_second_u_blck=nb_third_u_blck;
nb_second_u_blck=nb_third_u_blck+1;
nb_first_y_blck=nb_second_y_blck;
nb_second_y_blck=nb_third_y_blck;
first_u_blck = second_u_blck;
@ -1780,7 +1798,7 @@ SymbolicGaussElimination::Gaussian_Elimination(t_model_graph* model_graph
else
{
SaveCode.write(reinterpret_cast<char *>(&nb_table_u), sizeof(nb_table_u));
write_to_file_table_u_b(First_table_u->pNext, table_u->pNext, &nb_last_save_table_u );
write_to_file_table_u_b(First_table_u->pNext, table_u->pNext, &nb_last_save_table_u, false );
nb_last_save_table_y = vertex_count;
last_save_table_y=(t_table_y*)malloc(nb_last_save_table_y*sizeof(*last_save_table_y));
for(i=0;i<nb_last_save_table_y;i++)
@ -1834,6 +1852,18 @@ SymbolicGaussElimination::Gaussian_Elimination(t_model_graph* model_graph
free(s_j2);
}
void
SymbolicGaussElimination::file_is_open1()
{
file_open=true;
}
void
SymbolicGaussElimination::file_is_open()
{
file_open=true;
//file_is_open1();
}
void
SymbolicGaussElimination::SGE_compute(Model_Block *ModelBlock, int blck, bool dynamic, string file_name, int endo_nbr)
@ -1843,14 +1873,12 @@ SymbolicGaussElimination::SGE_compute(Model_Block *ModelBlock, int blck, bool dy
// pctimer_t t1, t2;
int i;
int mean_var_in_equation;
#ifdef PRINT_OUT
cout << "sizeof(i)=" << sizeof(i) << " sizeof(ttt)=" << sizeof(ttt) << " sizeof(short int)=" << sizeof(tti) << "\n" ;
system("pause");
#endif
init_glb();
model_graph = (t_model_graph*)malloc(sizeof(*model_graph));
nstacked = dynamic;
#ifdef PRINT_OUT
periods = ModelBlock->Periods;
cout << "nstacked=" << nstacked << "\n";
// t1 = pctimer();
cout << "periods=" << periods << "\n";
@ -1872,6 +1900,9 @@ SymbolicGaussElimination::SGE_compute(Model_Block *ModelBlock, int blck, bool dy
y_kmin = ModelBlock->Block_List[blck].Max_Lag;
y_kmax = ModelBlock->Block_List[blck].Max_Lead;
periods = ModelBlock->Periods;
if(periods<=3)
nstacked=false;
//cout << "periods=" << periods << " y_kmin=" << y_kmin << " y_kmax=" << y_kmax << "\n";
u_count_init = u_count;
#ifdef PRINT_OUT
cout << "size : " << size << "\n";

View File

@ -37,15 +37,17 @@ typedef struct vari
};
typedef map<pair<int ,int >,double> jacob_map;
class BlockTriangular
{
public:
Normalization normalization;
BlockTriangular(const SymbolTable &symbol_table_arg);
BlockTriangular(const SymbolTable &symbol_table_arg);
~BlockTriangular();
/*! The incidence matrix for each lead and lags */
Blocks blocks;
SymbolTable symbol_table;
Normalization normalization;
const SymbolTable &symbol_table;
List_IM* Build_IM(int lead_lag);
List_IM* Get_IM(int lead_lag);
bool* bGet_IM(int lead_lag);
@ -57,8 +59,9 @@ public:
void Free_IM(List_IM* First_IM);
void Print_SIM(bool* IM, int n);
void Normalize_and_BlockDecompose_Static_Model();
void Normalize_and_BlockDecompose_Static_0_Model();
void Normalize_and_BlockDecompose_Static_0_Model(const jacob_map &j_m);
bool Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock, int n, int* prologue, int* epilogue, simple* Index_Var_IM, simple* Index_Equ_IM, bool Do_Normalization, bool mixing, bool* IM_0 );
bool Normalize_and_BlockDecompose(bool* IM, Model_Block* ModelBlock, int n, int* prologue, int* epilogue, simple* Index_Var_IM, simple* Index_Equ_IM, bool Do_Normalization, bool mixing, bool* IM_0 , jacob_map j_m);
void Normalize_and_BlockDecompose_0();
void Normalize_and_BlockDecompose_Inside_Earth();
void Prologue_Epilogue(bool* IM, int* prologue, int* epilogue, int n, simple* Index_Var_IM, simple* Index_Equ_IM);
@ -69,6 +72,7 @@ public:
void Allocate_Block(int size, int *count_Equ, int *count_Block, int type, Model_Block * ModelBlock, int* Table, int TableSize);
void Free_Block(Model_Block* ModelBlock);
void SetVariableTable(int *Table,int Size,int HSize);
string getnamebyID(Type type, int id);
List_IM *First_IM ;
List_IM *Last_IM ;
simple *Index_Equ_IM;
@ -104,28 +108,30 @@ public:
{
switch (type)
{
case 0:
case EVALUATE_FOREWARD:
case EVALUATE_FOREWARD_R:
return ("EVALUATE FOREWARD ");
break;
case 1:
case EVALUATE_BACKWARD:
case EVALUATE_BACKWARD_R:
return ("EVALUATE BACKWARD ");
break;
case 2:
case SOLVE_FOREWARD_SIMPLE:
return ("SOLVE FOREWARD SIMPLE ");
break;
case 3:
case SOLVE_BACKWARD_SIMPLE:
return ("SOLVE BACKWARD SIMPLE ");
break;
case 4:
case SOLVE_TWO_BOUNDARIES_SIMPLE:
return ("SOLVE TWO BOUNDARIES SIMPLE ");
break;
case 5:
case SOLVE_FOREWARD_COMPLETE:
return ("SOLVE FOREWARD COMPLETE ");
break;
case 6:
case SOLVE_BACKWARD_COMPLETE:
return ("SOLVE BACKWARD COMPLETE ");
break;
case 7:
case SOLVE_TWO_BOUNDARIES_COMPLETE:
return ("SOLVE TWO BOUNDARIES COMPLETE");
break;
default:
@ -137,28 +143,30 @@ public:
{
switch (type)
{
case 0:
case EVALUATE_FOREWARD:
case EVALUATE_FOREWARD_R:
return ("EVALUATE_FOREWARD ");
break;
case 1:
case EVALUATE_BACKWARD:
case EVALUATE_BACKWARD_R:
return ("EVALUATE_BACKWARD ");
break;
case 2:
case SOLVE_FOREWARD_SIMPLE:
return ("SOLVE_FOREWARD_SIMPLE ");
break;
case 3:
case SOLVE_BACKWARD_SIMPLE:
return ("SOLVE_BACKWARD_SIMPLE ");
break;
case 4:
case SOLVE_TWO_BOUNDARIES_SIMPLE:
return ("SOLVE_TWO_BOUNDARIES_SIMPLE ");
break;
case 5:
case SOLVE_FOREWARD_COMPLETE:
return ("SOLVE_FOREWARD_COMPLETE ");
break;
case 6:
case SOLVE_BACKWARD_COMPLETE:
return ("SOLVE_BACKWARD_COMPLETE ");
break;
case 7:
case SOLVE_TWO_BOUNDARIES_COMPLETE:
return ("SOLVE_TWO_BOUNDARIES_COMPLETE");
break;
default:

View File

@ -104,9 +104,9 @@ public:
class CutoffStatement : public Statement
{
private:
const int cutoff;
const double cutoff;
public:
CutoffStatement(int cutoff_arg);
CutoffStatement(double cutoff_arg);
virtual void writeOutput(ostream &output, const string &basename) const;
};

View File

@ -25,7 +25,7 @@ protected:
SymbolTable &symbol_table;
//! Reference to numerical constants table
NumericalConstants &num_constants;
typedef list<NodeID> node_list_type;
//! The list of nodes
node_list_type node_list;
@ -128,7 +128,7 @@ DataTree::AddPossiblyNegativeConstant(double v)
}
ostringstream ost;
ost << v;
NodeID cnode = AddNumConstant(ost.str());
if (neg)
@ -152,7 +152,8 @@ DataTree::AddUnaryOp(UnaryOpcode op_code, NodeID arg)
{
try
{
double argval = arg->eval(eval_context_type());
eval_context_type *evc=NULL;
double argval = arg->eval(/*eval_context_type()*/*evc);
double val = UnaryOpNode::eval_opcode(op_code, argval);
return AddPossiblyNegativeConstant(val);
}
@ -160,7 +161,6 @@ DataTree::AddUnaryOp(UnaryOpcode op_code, NodeID arg)
{
}
}
return new UnaryOpNode(*this, op_code, arg);
}
@ -174,15 +174,15 @@ DataTree::AddBinaryOp(NodeID arg1, BinaryOpcode op_code, NodeID arg2)
// Try to reduce to a constant
try
{
double argval1 = arg1->eval(eval_context_type());
double argval2 = arg2->eval(eval_context_type());
eval_context_type *evc=NULL;
double argval1 = arg1->eval(/*eval_context_type()*/*evc);
double argval2 = arg2->eval(/*eval_context_type()*/*evc);
double val = BinaryOpNode::eval_opcode(argval1, op_code, argval2);
return AddPossiblyNegativeConstant(val);
}
catch(ExprNode::EvalException &e)
{
}
return new BinaryOpNode(*this, arg1, op_code, arg2);
}

View File

@ -1,6 +1,7 @@
#ifndef MODELNORMALIZATION
#define MODELNORMALIZATION
#include "SymbolTableTypes.hh"
#include "SymbolTable.hh"
const int SIMULTANS=0;
const int PROLOGUE=1;
const int EPILOGUE=2;
@ -14,6 +15,8 @@ const int SOLVE_TWO_BOUNDARIES_SIMPLE=4;
const int SOLVE_FOREWARD_COMPLETE=5;
const int SOLVE_BACKWARD_COMPLETE=6;
const int SOLVE_TWO_BOUNDARIES_COMPLETE=7;
const int EVALUATE_FOREWARD_R=8;
const int EVALUATE_BACKWARD_R=9;
typedef struct Edge
{
@ -62,11 +65,13 @@ private:
int v; /* current matched of u */
};
public:
Normalization(/*t_getNameByID gdi*/);
Normalization(const SymbolTable &symbol_table_arg);
~Normalization();
void Normalize(int n, int prologue, int epilogue, bool* IM, simple* Index_Var_IM, Equation_set* Equation,bool mixing, bool* IM_s);
bool Normalize(int n, int prologue, int epilogue, bool* IM, simple* Index_Var_IM, Equation_set* Equation,bool mixing, bool* IM_s);
void Gr_to_IM_basic(int n0, int prologue, int epilogue, bool* IM, Equation_set *Equation,bool transpose);
t_getNameByID getnamebyID;
const SymbolTable &symbol_table;
void Set_fp_verbose(bool ok);
private:
void IM_to_Gr(int n0, int prologue, int epilogue, bool* IM, Equation_set *Equation, Variable_set *Variable );
void Inits(Equation_set *Equation);

View File

@ -12,9 +12,11 @@ using namespace std;
#include "NumericalConstants.hh"
#include "DataTree.hh"
#include "BlockTriangular.hh"
#include "SymbolGaussElim.hh"
#define LCC_COMPILE 0
#define GCC_COMPILE 1
//#define CONDITION
//! The three in which ModelTree can work
enum ModelTreeMode
@ -94,7 +96,7 @@ private:
void writeSparseDLLDynamicHFile(const string &dynamic_basename) const;
//! Writes dynamic model file when SparseDLL option is on
void writeSparseDLLDynamicCFileAndBinFile(const string &dynamic_basename, const string &bin_basename) const;
void evaluateJacobian(const eval_context_type &eval_context);
void evaluateJacobian(const eval_context_type &eval_context, jacob_map *j_m);
void BlockLinear(Model_Block *ModelBlock);
string reform(string name) const;
@ -106,7 +108,9 @@ public:
int compiler;
//! Absolute value under which a number is considered to be zero
double cutoff;
//! Use a graphical and symbolic version of the symbolic gaussian elimination new_SGE = false
//! or use direct gaussian elimination new_SGE = true
bool new_SGE;
//! Declare a node as an equation of the model
void addEquation(NodeID eq);
//! Do some checking
@ -132,6 +136,10 @@ public:
void writeDynamicFile(const string &basename) const;
//! Complete set to block decompose the model
BlockTriangular block_triangular;
//! Adds informations for simulation in a binary file
void Write_Inf_To_Bin_File(const string &dynamic_basename, const string &bin_basename, const int &num,
int &u_count_int, bool &file_open) const;
int equation_number() const;
};

View File

@ -6,7 +6,8 @@
#define SIMPLIFYS
#define SAVE
#define COMPUTE
#define PRINT_OUT_OUT
//#define PRINT_OUT_OUT
//#define PRINT_OUT_1
#define DIRECT_SAVE
#include "ModelTree.hh"
#include "BlockTriangular.hh"

View File

@ -75,7 +75,7 @@ public:
void init_glb();
void write_to_file_table_y( t_table_y *save_table_y, t_table_y *save_i_table_y, int nb_save_table_y, int nb_save_i_table_y);
void write_to_file_table_u(t_table_u *save_table_u,t_table_u *save_i_table_u, int nb_save_table_u);
void write_to_file_table_u_b(t_table_u *save_table_u, t_table_u *last_table_u, int *nb_save_table_u);
void write_to_file_table_u_b(t_table_u *save_table_u, t_table_u *last_table_u, int *nb_save_table_u, bool chk);
bool Check_Regularity(t_table_u *first_u_blck, t_table_u *second_u_blck, t_table_u *third_u_blck);
t_table_u* interpolation(t_model_graph* model_graph,t_table_y* table_y, int to_add, bool middle,int per);
bool Loop_Elimination(t_model_graph* model_graph);
@ -87,6 +87,9 @@ public:
, bool dynamic);
int SGE_all(int endo,int Time, List_IM *First_IM);
void SGE_compute(Model_Block *ModelBlock, int blck, bool dynamic, string file_name, int endo_nbr);
void file_is_open();
void file_is_open1();
};
//------------------------------------------------------------------------------
#endif

8024
tests/ferhat/Doquier.mod Normal file

File diff suppressed because it is too large Load Diff

131
tests/ferhat/MARK3_endo.dat Normal file

File diff suppressed because one or more lines are too long

131
tests/ferhat/MARK3_exo.dat Normal file

File diff suppressed because one or more lines are too long

View File

@ -13,14 +13,14 @@ m=1;
n=1;
o=1;
model(SPARSE_DLL,gcc_compiler);
model(SPARSE_DLL,gcc_compiler,cutoff=1e-12);
/*0*/ k=(1-h)*k(-1)+i; /*k:0*/
/*1*/ y=l^j*k^m; /*l:1*/
/*2*/ c=y*a+b+0.3*c(-1)+0.1*c(+1)+0.*g_bar(-10); /*c:2*/
/*3*/ infl=0.02*y+0.5*r; /*infl:3*/
/*4*/ i=d*(y-y(-1))+e/**r*/; /*i4*/
/*5*/ g=f*g_bar; /*g:5*/
/*6*/ y=0.6*(c+i+g)+0.1*y(-2)+0.1*y(+2)+0.1*y(-1)+0.1*y(+1); /*y:6*/
/*6*/ y=0.6*(c+i+g)+/*0.1*y(-2)+0.1*y(+2)+*/0.1*y(-1)+0.1*y(+1); /*y:6*/
/*7*/ r=y-1+infl-0.02; /*r7*/
/*8*/ p1=i+0.5*q1;
/*9*/ q1=0.5*p1+c;
@ -50,16 +50,16 @@ end;
steady(solve_algo=2);
//check;
/*shocks;
shocks;
var g_bar;
periods 1;
values 0.16;
end;*/
end;
options_.slowc = 1;
simul(periods=100);
simul(periods=80);
rplot c;

103
tests/ferhat/fs2000.mod Normal file
View File

@ -0,0 +1,103 @@
// This file replicates the estimation of the CIA model from
// Frank Schorfheide (2000) "Loss function-based evaluation of DSGE models"
// Journal of Applied Econometrics, 15, 645-670.
// the data are the ones provided on Schorfheide's web site with the programs.
// http://www.econ.upenn.edu/~schorf/programs/dsgesel.ZIP
// You need to have fsdat.m in the same directory as this file.
// This file replicates:
// -the posterior mode as computed by Frank's Gauss programs
// -the parameter mean posterior estimates reported in the paper
// -the model probability (harmonic mean) reported in the paper
// This file was tested with dyn_mat_test_0218.zip
// the smooth shocks are probably stil buggy
//
// The equations are taken from J. Nason and T. Cogley (1994)
// "Testing the implications of long-run neutrality for monetary business
// cycle models" Journal of Applied Econometrics, 9, S37-S70.
// 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;
varexo e_a e_m;
parameters alp bet gam mst rho psi del;
alp = 0.33;
bet = 0.99;
gam = 0.003;
mst = 1.011;
rho = 0.7;
psi = 0.787;
del = 0.02;
model(sparse_dll,gcc_compiler,cutoff=1e-17);
dA = exp(gam+e_a);
log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m;
-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;
W = l/n;
-(psi/(1-psi))*(c*P/(1-n))+l/n = 0;
R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W;
1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0;
c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1);
P*c = m;
m-1+d = l;
e = exp(e_a);
y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a));
gy_obs = dA*y/y(-1);
gp_obs = (P/P(-1))*m(-1)/dA;
end;
initval;
k = 6;
m = mst;
P = 2.25;
c = 0.45;
e = 1;
W = 4;
R = 1.02;
d = 0.85;
n = 0.19;
l = 0.86;
y = 0.6;
gy_obs = exp(gam);
gp_obs = exp(-gam);
dA = exp(gam);
e_a=0;
e_m=0;
end;
/*shocks;
var e_a; stderr 0.014;
var e_m; stderr 0.005;
end;
*/
steady;
shocks;
var e_a;
periods 1;
values 0.16;
end;
simul(periods=80);
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);
*/

124
tests/ferhat/ireland.mod Normal file
View File

@ -0,0 +1,124 @@
var y a k c i h eoy eoc eoh oy oc oh;
varexo e eeoy eeoc eeoh;
parameters theta rho eta gam bet delta aa r11 r12 r13 r21 r22 r23 r31 r32 r33 scy shc shy;
bet = 0.99;
delta = 0.025;
theta = 0.2;
rho = 0.9959;
eta = 1.0051;
gam = 0.0045;
aa = 1.8;
r11 = 0.99;
r12 = 0;
r13 = 0;
r21 = 0;
r22 = 0.99;
r23 = 0;
r31 = 0;
r32 = 0;
r33 = 0.99;
scy = 0.0040;
shy = 0.0015;
shc = 0.0010;
model(sparse_dll,gcc_compiler);
exp(y) = exp(a)*exp(k(-1))^theta*exp(h)^(1-theta);
a = (1-rho)*aa+rho*a(-1)+e;
exp(y) = exp(c) + exp(i);
eta*exp(k) = (1-delta)*exp(k(-1))+exp(i);
gam*exp(c)*exp(h) = (1-theta)*exp(y);
eta/exp(c) = bet*(1/exp(c(+1)))*(theta*(exp(y(+1))/exp(k))+1-delta);
eoy = r11*eoy(-1) + r12*eoc(-1) + r13*eoh(-1) + eeoy;
eoc = r21*eoy(-1) + r22*eoc(-1) + r23*eoh(-1) + scy*eeoy+eeoc;
eoh = r31*eoy(-1) + r32*eoc(-1) + r33*eoh(-1) + shy*eeoy+shc*eeoc+eeoh;
oy = y + eoy;
oc = c + eoc;
oh = h + eoh;
end;
initval;
/*a = 1.7;
y = 8;
c = 8;
k = 10;
i = 5;
h = 4;
eoy = 0;
eoc = 0;
eoh = 0;
oy = y;
oc = c;
oh = h;
*/
e=0;
eeoy=0;
eeoc=0;
eeoh=0;
y= 7.99331700544506;
a= 1.8;
k= 9.59646163090336;
c= 7.83132048725623;
i= 6.09323152367607;
h= 5.34253084908048;
eoy= 0.0000;
eoc= 0.0000;
eoh= 0;
oy= 7.99331700544506;
oc= 7.83132048725623;
oh= 5.34253084908048;
end;
options_.dynatol=1e-12;
options_.maxit_=500;
options_.slowc=1;
steady(solve_algo=3);
options_.dynatol=4e-8;
//check;
shocks;
var e;
periods 1;
values -0.0002;
end;
simul(periods=60);
rplot y;
rplot k;
/*estimated_params;
theta , 0.22, 0.1, 0.5;
rho , 0.99, 0.7, 0.9999;
eta , 1.0051, 1, 1.03;
gam , 0.0045, 0.001, 0.01;
aa , 1.8, 0.1, 4;
r11 , 1.4187, -2, 2;
r12 , 0.2251, -2, 2;
r13 , -0.4441, -2, 2;
r21 , 0.0935, -2, 2;
r22 , 1.0236, -2, 2;
r23 , -0.0908, -2, 2;
r31 , 0.7775, -2, 2;
r32 , 0.3706, -2, 2;
r33 , 0.2398, -2, 2;
scy , 0.0040, -2, 2;
shy , 0.0015, -2, 2;
shc , 0.0010, -2, 2;
stderr e , 0.0056, 0, 0.2;
stderr eeoy , 0.0070, 0, 0.1;
stderr eeoc , 0.0069, 0, 0.1;
stderr eeoh , 0.0018, 0, 0.1;
end;
varobs oy oc oh;
observation_trends;
oy (log(eta));
oc (log(eta));
end;
estimation(datafile=idata,mode_compute=1,nograph);
*/

79
tests/ferhat/ls2003.mod Normal file
View File

@ -0,0 +1,79 @@
var y y_s R pie dq pie_s de A y_obs pie_obs R_obs;
varexo e_R e_q e_ys e_pies e_A;
parameters psi1 psi2 psi3 rho_R tau alpha rr k rho_q rho_A rho_ys rho_pies;
psi1 = 1.54;
psi2 = 0.25;
psi3 = 0.25;
rho_R = 0.5;
alpha = 0.3;
rr = 2.51;
k = 0.5;
tau = 0.5;
rho_q = 0.4;
rho_A = 0.2;
rho_ys = 0.9;
rho_pies = 0.7;
model(sparse_dll,gcc_compiler,cutoff=1e-17);
y = y(+1) - (tau +alpha*(2-alpha)*(1-tau))*(R-pie(+1))-alpha*(tau +alpha*(2-alpha)*(1-tau))*dq(+1) + alpha*(2-alpha)*((1-tau)/tau)*(y_s-y_s(+1))-A(+1);
pie = exp(-rr/400)*pie(+1)+alpha*exp(-rr/400)*dq(+1)-alpha*dq+(k/(tau+alpha*(2-alpha)*(1-tau)))*y+alpha*(2-alpha)*(1-tau)/(tau*(tau+alpha*(2-alpha)*(1-tau)))*y_s;
pie = de+(1-alpha)*dq+pie_s;
R = rho_R*R(-1)+(1-rho_R)*(psi1*pie+psi2*(y+alpha*(2-alpha)*((1-tau)/tau)*y_s)+psi3*de)+e_R;
dq = rho_q*dq(-1)+e_q;
y_s = rho_ys*y_s(-1)+e_ys;
pie_s = rho_pies*pie_s(-1)+e_pies;
A = rho_A*A(-1)+e_A;
y_obs = y-y(-1)+A;
pie_obs = 4*pie;
R_obs = 4*R;
end;
/*shocks;
var e_R = 1.25^2;
var e_q = 2.5^2;
var e_A = 1.89;
var e_ys = 1.89;
var e_pies = 1.89;
end;
varobs y_obs R_obs pie_obs dq de;
estimated_params;
psi1 , gamma_pdf,1.5,0.5;
psi2 , gamma_pdf,0.25,0.125;
psi3 , gamma_pdf,0.25,0.125;
rho_R ,beta_pdf,0.5,0.2;
alpha ,beta_pdf,0.3,0.1;
rr ,gamma_pdf,2.5,1;
k , gamma_pdf,0.5,0.25;
tau ,gamma_pdf,0.5,0.2;
rho_q ,beta_pdf,0.4,0.2;
rho_A ,beta_pdf,0.5,0.2;
rho_ys ,beta_pdf,0.8,0.1;
rho_pies,beta_pdf,0.7,0.15;
stderr e_R,inv_gamma_pdf,1.2533,0.6551;
stderr e_q,inv_gamma_pdf,2.5066,1.3103;
stderr e_A,inv_gamma_pdf,1.2533,0.6551;
stderr e_ys,inv_gamma_pdf,1.2533,0.6551;
stderr e_pies,inv_gamma_pdf,1.88,0.9827;
end;
estimation(datafile=data_ca1,first_obs=8,nobs=79,mh_nblocks=10,prefilter=1,mh_jscale=0.5,mh_replic=0);
*/
steady;
shocks;
var e_q;
periods 1;
values 0.5;
end;
simul(periods=80);
rplot A;
rplot pie;

4558
tests/ferhat/multimod.mod Normal file

File diff suppressed because it is too large Load Diff

View File

@ -9,7 +9,7 @@ bet=0.05;
aa=0.5;
model/*(SPARSE_DLL,GCC_COMPILER)*/;
model(SPARSE_DLL,GCC_COMPILER);
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;
@ -30,7 +30,10 @@ periods 1;
values 1.2;
end;
simul(periods=200);
oo_.dynatol=1e-19;
options_.maxit_ = 10;
simul(periods=6);
rplot c;
rplot k;

File diff suppressed because it is too large Load Diff