The local variables are implemented with byte code
git-svn-id: https://www.dynare.org/svn/dynare/trunk@3015 ac1d8469-bf42-47a9-8791-bf33cf982152time-shift
parent
541129c1be
commit
ece194b602
|
@ -131,7 +131,7 @@ Interpreter::log1(double a)
|
||||||
|
|
||||||
|
|
||||||
void
|
void
|
||||||
Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException)*/
|
Interpreter::compute_block_time(int Per_u_, bool evaluate, int block_num) /*throw(EvalException)*/
|
||||||
{
|
{
|
||||||
int var, lag, op;
|
int var, lag, op;
|
||||||
ostringstream tmp_out;
|
ostringstream tmp_out;
|
||||||
|
@ -179,7 +179,7 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException
|
||||||
Stack.push(x[it_+lag+var*nb_row_xd]);
|
Stack.push(x[it_+lag+var*nb_row_xd]);
|
||||||
break;
|
break;
|
||||||
default:
|
default:
|
||||||
mexPrintf("Unknown variable type\n");
|
mexPrintf("FLDV: Unknown variable type\n");
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case FLDSV :
|
case FLDSV :
|
||||||
|
@ -214,8 +214,12 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException
|
||||||
var=get_code_int;
|
var=get_code_int;
|
||||||
Stack.push(x[var]);
|
Stack.push(x[var]);
|
||||||
break;
|
break;
|
||||||
|
case eModelLocalVariable:
|
||||||
|
/*mexPrintf("FLDSV a local variable in Block %d Stack.size()=%d",block_num, Stack.size());
|
||||||
|
mexPrintf(" value=%f\n",Stack.top());*/
|
||||||
|
break;
|
||||||
default:
|
default:
|
||||||
mexPrintf("Unknown variable type\n");
|
mexPrintf("FLDSV: Unknown variable type\n");
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case FLDVS :
|
case FLDVS :
|
||||||
|
@ -251,7 +255,7 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException
|
||||||
Stack.push(x[var]);
|
Stack.push(x[var]);
|
||||||
break;
|
break;
|
||||||
default:
|
default:
|
||||||
mexPrintf("Unknown variable type\n");
|
mexPrintf("FLDVS: Unknown variable type\n");
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case FLDT :
|
case FLDT :
|
||||||
|
@ -343,7 +347,7 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException
|
||||||
Stack.pop();
|
Stack.pop();
|
||||||
break;
|
break;
|
||||||
default:
|
default:
|
||||||
mexPrintf("Unknown vraibale type\n");
|
mexPrintf("FSTPV: Unknown vraibale type\n");
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case FSTPSV :
|
case FSTPSV :
|
||||||
|
@ -376,7 +380,7 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException
|
||||||
Stack.pop();
|
Stack.pop();
|
||||||
break;
|
break;
|
||||||
default:
|
default:
|
||||||
mexPrintf("Unknown vraibale type\n");
|
mexPrintf("FSTPSV: Unknown vraibale type\n");
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case FSTPT :
|
case FSTPT :
|
||||||
|
@ -711,7 +715,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
{
|
{
|
||||||
case EVALUATE_FORWARD :
|
case EVALUATE_FORWARD :
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
begining=get_code_pointer;
|
begining=get_code_pointer;
|
||||||
|
@ -720,7 +724,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
|
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
@ -729,7 +733,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
r=(double*)mxMalloc(size*sizeof(double));
|
r=(double*)mxMalloc(size*sizeof(double));
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
{
|
{
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[Block_Contain[j].Variable] += r[j];
|
y[Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -741,7 +745,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
|
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -755,7 +759,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
r=(double*)mxMalloc(size*sizeof(double));
|
r=(double*)mxMalloc(size*sizeof(double));
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
{
|
{
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[Block_Contain[j].Variable] += r[j];
|
y[Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -767,7 +771,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
|
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -776,7 +780,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
break;
|
break;
|
||||||
case EVALUATE_BACKWARD :
|
case EVALUATE_BACKWARD :
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
begining=get_code_pointer;
|
begining=get_code_pointer;
|
||||||
|
@ -784,7 +788,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
@ -793,7 +797,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
r=(double*)mxMalloc(size*sizeof(double));
|
r=(double*)mxMalloc(size*sizeof(double));
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
{
|
{
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[Block_Contain[j].Variable] += r[j];
|
y[Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -804,7 +808,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0,true);
|
compute_block_time(0,true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -818,7 +822,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
r=(double*)mxMalloc(size*sizeof(double));
|
r=(double*)mxMalloc(size*sizeof(double));
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
{
|
{
|
||||||
compute_block_time(0, true);
|
compute_block_time(0, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[Block_Contain[j].Variable] += r[j];
|
y[Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -829,7 +833,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0,true);
|
compute_block_time(0,true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -848,7 +852,7 @@ Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool stea
|
||||||
Per_u_=(it_-y_kmin)*u_count_int;
|
Per_u_=(it_-y_kmin)*u_count_int;
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
compute_block_time(Per_u_, true);
|
compute_block_time(Per_u_, true, block_num);
|
||||||
for(int j=0; j<size; j++)
|
for(int j=0; j<size; j++)
|
||||||
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
y[it_*y_size+Block_Contain[j].Variable] += r[j];
|
||||||
}
|
}
|
||||||
|
@ -875,7 +879,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
case EVALUATE_FORWARD :
|
case EVALUATE_FORWARD :
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
begining=get_code_pointer;
|
begining=get_code_pointer;
|
||||||
|
@ -883,13 +887,13 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
case EVALUATE_BACKWARD :
|
case EVALUATE_BACKWARD :
|
||||||
if(steady_state)
|
if(steady_state)
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
begining=get_code_pointer;
|
begining=get_code_pointer;
|
||||||
|
@ -897,7 +901,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
@ -913,7 +917,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
y[Block_Contain[0].Variable] += -r[0]/g1[0];
|
y[Block_Contain[0].Variable] += -r[0]/g1[0];
|
||||||
double rr;
|
double rr;
|
||||||
rr=r[0];
|
rr=r[0];
|
||||||
|
@ -923,6 +927,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
if (!cvg)
|
if (!cvg)
|
||||||
{
|
{
|
||||||
mexPrintf("Convergence not achieved in block %d, after %d iterations\n",Block_Count,iter);
|
mexPrintf("Convergence not achieved in block %d, after %d iterations\n",Block_Count,iter);
|
||||||
|
mexPrintf("r[0]= %f\n",r[0]);
|
||||||
/*mexEvalString("st=fclose('all');clear all;");
|
/*mexEvalString("st=fclose('all');clear all;");
|
||||||
mexErrMsgTxt("End of simulate");*/
|
mexErrMsgTxt("End of simulate");*/
|
||||||
return false;
|
return false;
|
||||||
|
@ -939,7 +944,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
|
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
|
||||||
double rr;
|
double rr;
|
||||||
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
|
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
|
||||||
|
@ -972,7 +977,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
y[Block_Contain[0].Variable] += -r[0]/g1[0];
|
y[Block_Contain[0].Variable] += -r[0]/g1[0];
|
||||||
double rr;
|
double rr;
|
||||||
rr=r[0];
|
rr=r[0];
|
||||||
|
@ -998,7 +1003,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
{
|
{
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
|
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
|
||||||
double rr;
|
double rr;
|
||||||
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
|
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
|
||||||
|
@ -1049,7 +1054,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
res2=0;
|
res2=0;
|
||||||
res1=0;
|
res1=0;
|
||||||
max_res=0;
|
max_res=0;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
/*if (isnan(res1)||isinf(res1))
|
/*if (isnan(res1)||isinf(res1))
|
||||||
{
|
{
|
||||||
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
||||||
|
@ -1091,7 +1096,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
iter = 0;
|
iter = 0;
|
||||||
res1=res2=max_res=0;max_res_idx=0;
|
res1=res2=max_res=0;max_res_idx=0;
|
||||||
error_not_printed = true;
|
error_not_printed = true;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
cvg=false;
|
cvg=false;
|
||||||
result = simulate_NG(Block_Count, symbol_table_endo_nbr, 0, 0, 0, size, false, cvg, iter, true);
|
result = simulate_NG(Block_Count, symbol_table_endo_nbr, 0, 0, 0, size, false, cvg, iter, true);
|
||||||
}
|
}
|
||||||
|
@ -1115,7 +1120,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
res2=0;
|
res2=0;
|
||||||
res1=0;
|
res1=0;
|
||||||
max_res=0;
|
max_res=0;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
/*if (isnan(res1)||isinf(res1))
|
/*if (isnan(res1)||isinf(res1))
|
||||||
{
|
{
|
||||||
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
||||||
|
@ -1161,7 +1166,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
iter = 0;
|
iter = 0;
|
||||||
res1=res2=max_res=0;max_res_idx=0;
|
res1=res2=max_res=0;max_res_idx=0;
|
||||||
error_not_printed = true;
|
error_not_printed = true;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
cvg=false;
|
cvg=false;
|
||||||
result = simulate_NG(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, false, cvg, iter, false);
|
result = simulate_NG(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, false, cvg, iter, false);
|
||||||
}
|
}
|
||||||
|
@ -1203,7 +1208,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
res2=0;
|
res2=0;
|
||||||
res1=0;
|
res1=0;
|
||||||
max_res=0;
|
max_res=0;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
/*if (isnan(res1)||isinf(res1))
|
/*if (isnan(res1)||isinf(res1))
|
||||||
{
|
{
|
||||||
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
||||||
|
@ -1244,7 +1249,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
iter = 0;
|
iter = 0;
|
||||||
res1=res2=max_res=0;max_res_idx=0;
|
res1=res2=max_res=0;max_res_idx=0;
|
||||||
error_not_printed = true;
|
error_not_printed = true;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
cvg=false;
|
cvg=false;
|
||||||
result = simulate_NG(Block_Count, symbol_table_endo_nbr, 0, 0, 0, size, false, cvg, iter, true);
|
result = simulate_NG(Block_Count, symbol_table_endo_nbr, 0, 0, 0, size, false, cvg, iter, true);
|
||||||
}
|
}
|
||||||
|
@ -1266,7 +1271,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
res2=0;
|
res2=0;
|
||||||
res1=0;
|
res1=0;
|
||||||
max_res=0;
|
max_res=0;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
/*if (isnan(res1)||isinf(res1))
|
/*if (isnan(res1)||isinf(res1))
|
||||||
{
|
{
|
||||||
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
||||||
|
@ -1310,7 +1315,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
error_not_printed = true;
|
error_not_printed = true;
|
||||||
compute_block_time(0, false);
|
compute_block_time(0, false, block_num);
|
||||||
cvg=false;
|
cvg=false;
|
||||||
result = simulate_NG(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, false, cvg, iter, false);
|
result = simulate_NG(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, false, cvg, iter, false);
|
||||||
}
|
}
|
||||||
|
@ -1363,7 +1368,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
Per_u_=(it_-y_kmin)*u_count_int;
|
Per_u_=(it_-y_kmin)*u_count_int;
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
compute_block_time(Per_u_, false);
|
compute_block_time(Per_u_, false, block_num);
|
||||||
if (isnan(res1)||isinf(res1))
|
if (isnan(res1)||isinf(res1))
|
||||||
{
|
{
|
||||||
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
|
||||||
|
@ -1408,7 +1413,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
|
||||||
Per_u_=(it_-y_kmin)*u_count_int;
|
Per_u_=(it_-y_kmin)*u_count_int;
|
||||||
Per_y_=it_*y_size;
|
Per_y_=it_*y_size;
|
||||||
set_code_pointer(begining);
|
set_code_pointer(begining);
|
||||||
compute_block_time(Per_u_, false);
|
compute_block_time(Per_u_, false, block_num);
|
||||||
for (i=0; i< size; i++)
|
for (i=0; i< size; i++)
|
||||||
{
|
{
|
||||||
double rr;
|
double rr;
|
||||||
|
|
|
@ -68,7 +68,7 @@ class Interpreter : SparseMatrix
|
||||||
protected :
|
protected :
|
||||||
double pow1(double a, double b);
|
double pow1(double a, double b);
|
||||||
double log1(double a);
|
double log1(double a);
|
||||||
void compute_block_time(int Per_u_, bool evaluate);
|
void compute_block_time(int Per_u_, bool evaluate, int block_num);
|
||||||
void evaluate_a_block(int size, int type, string bin_basename, bool steady_state, int block_num);
|
void evaluate_a_block(int size, int type, string bin_basename, bool steady_state, int block_num);
|
||||||
bool simulate_a_block(int size, int type, string file_name, string bin_basename, bool Gaussian_Elimination, bool steady_state, int block_num);
|
bool simulate_a_block(int size, int type, string file_name, string bin_basename, bool Gaussian_Elimination, bool steady_state, int block_num);
|
||||||
double *T;
|
double *T;
|
||||||
|
|
|
@ -120,6 +120,9 @@ BlockTriangular::Compute_Normalization(bool *IM, int equation_number, int prolog
|
||||||
|
|
||||||
typedef adjacency_list<vecS, vecS, undirectedS> BipartiteGraph;
|
typedef adjacency_list<vecS, vecS, undirectedS> BipartiteGraph;
|
||||||
|
|
||||||
|
|
||||||
|
if(verbose == 2)
|
||||||
|
cout << "trying to normlized even in singular case\n";
|
||||||
/*
|
/*
|
||||||
Vertices 0 to n-1 are for endogenous (using type specific ID)
|
Vertices 0 to n-1 are for endogenous (using type specific ID)
|
||||||
Vertices n to 2*n-1 are for equations (using equation no.)
|
Vertices n to 2*n-1 are for equations (using equation no.)
|
||||||
|
@ -148,7 +151,7 @@ BlockTriangular::Compute_Normalization(bool *IM, int equation_number, int prolog
|
||||||
{
|
{
|
||||||
if (verbose == 1)
|
if (verbose == 1)
|
||||||
{
|
{
|
||||||
cerr << "WARNING: Could not normalize dynamic model. Variable "
|
cout << "WARNING: Could not normalize dynamic model. Variable "
|
||||||
<< symbol_table.getName(symbol_table.getID(eEndogenous, it - mate_map.begin()))
|
<< symbol_table.getName(symbol_table.getID(eEndogenous, it - mate_map.begin()))
|
||||||
<< " is not in the maximum cardinality matching. Trying to find a singular normalization." << endl;
|
<< " is not in the maximum cardinality matching. Trying to find a singular normalization." << endl;
|
||||||
//exit(EXIT_FAILURE);
|
//exit(EXIT_FAILURE);
|
||||||
|
@ -169,6 +172,7 @@ BlockTriangular::Compute_Normalization(bool *IM, int equation_number, int prolog
|
||||||
memcpy(SIM, IM, equation_number*equation_number*sizeof(bool));
|
memcpy(SIM, IM, equation_number*equation_number*sizeof(bool));
|
||||||
for (int i = 0; i < n; i++)
|
for (int i = 0; i < n; i++)
|
||||||
{
|
{
|
||||||
|
//printf("match equation %4d with variable %4d \n", mate_map[i] - n, i);
|
||||||
Index_Equ_IM[i + prologue] = Index_Equ_IM_tmp[mate_map[i] - n + prologue];
|
Index_Equ_IM[i + prologue] = Index_Equ_IM_tmp[mate_map[i] - n + prologue];
|
||||||
for (int k = 0; k < n; k++)
|
for (int k = 0; k < n; k++)
|
||||||
IM[(i+prologue)*equation_number+k +prologue] = SIM[(mate_map[i]-n+prologue)*equation_number+k +prologue];
|
IM[(i+prologue)*equation_number+k +prologue] = SIM[(mate_map[i]-n+prologue)*equation_number+k +prologue];
|
||||||
|
@ -337,6 +341,7 @@ BlockTriangular::Allocate_Block(int size, int *count_Equ, int count_Block, Block
|
||||||
bool *Cur_IM;
|
bool *Cur_IM;
|
||||||
bool *IM, OK;
|
bool *IM, OK;
|
||||||
int Lag_Endo, Lead_Endo, Lag_Exo, Lead_Exo, Lag_Other_Endo, Lead_Other_Endo;
|
int Lag_Endo, Lead_Endo, Lag_Exo, Lead_Exo, Lag_Other_Endo, Lead_Other_Endo;
|
||||||
|
//cout << "block " << count_Block << " size " << size << " SimType=" << BlockSim(SimType) << "\n";
|
||||||
|
|
||||||
ModelBlock->Periods = periods;
|
ModelBlock->Periods = periods;
|
||||||
ModelBlock->Block_List[count_Block].is_linear = true;
|
ModelBlock->Block_List[count_Block].is_linear = true;
|
||||||
|
@ -725,35 +730,35 @@ BlockTriangular::Equation_Type_determination(vector<BinaryOpNode *> &equations,
|
||||||
lhs->writeOutput(tmp_output, oMatlabDynamicModelSparse, temporary_terms);
|
lhs->writeOutput(tmp_output, oMatlabDynamicModelSparse, temporary_terms);
|
||||||
tmp_s << "y(it_, " << Index_Var_IM[i]+1 << ")";
|
tmp_s << "y(it_, " << Index_Var_IM[i]+1 << ")";
|
||||||
map<pair<int, pair<int, int> >, NodeID>::iterator derivative = first_order_endo_derivatives.find(make_pair(eq, make_pair(var, 0)));
|
map<pair<int, pair<int, int> >, NodeID>::iterator derivative = first_order_endo_derivatives.find(make_pair(eq, make_pair(var, 0)));
|
||||||
set<pair<int, int> > result;
|
|
||||||
pair<bool, NodeID> res;
|
pair<bool, NodeID> res;
|
||||||
derivative->second->collectEndogenous(result);
|
if(derivative != first_order_endo_derivatives.end())
|
||||||
set<pair<int, int> >::const_iterator d_endo_variable = result.find(make_pair(var, 0));
|
|
||||||
//Determine whether the equation could be evaluated rather than to be solved
|
|
||||||
ostringstream tt("");
|
|
||||||
derivative->second->writeOutput(tt, oMatlabDynamicModelSparse, temporary_terms);
|
|
||||||
//cout << tt.str().c_str() << " tmp_output.str()=" << tmp_output.str() << " tmp_s.str()=" << tmp_s.str();
|
|
||||||
if (tmp_output.str() == tmp_s.str() and tt.str()=="1")
|
|
||||||
{
|
{
|
||||||
Equation_Simulation_Type = E_EVALUATE;
|
set<pair<int, int> > result;
|
||||||
//cout << " E_EVALUATE ";
|
derivative->second->collectEndogenous(result);
|
||||||
}
|
set<pair<int, int> >::const_iterator d_endo_variable = result.find(make_pair(var, 0));
|
||||||
else
|
//Determine whether the equation could be evaluated rather than to be solved
|
||||||
{
|
ostringstream tt("");
|
||||||
vector<pair<int, pair<NodeID, NodeID> > > List_of_Op_RHS;
|
derivative->second->writeOutput(tt, oMatlabDynamicModelSparse, temporary_terms);
|
||||||
res = equations[eq]->normalizeEquation(var, List_of_Op_RHS);
|
if (tmp_output.str() == tmp_s.str() and tt.str()=="1")
|
||||||
if(mfs==2)
|
|
||||||
{
|
{
|
||||||
if(d_endo_variable == result.end() && res.second)
|
Equation_Simulation_Type = E_EVALUATE;
|
||||||
Equation_Simulation_Type = E_EVALUATE_S;
|
|
||||||
}
|
}
|
||||||
else if(mfs==3)
|
else
|
||||||
{
|
{
|
||||||
if(res.second) // The equation could be solved analytically
|
vector<pair<int, pair<NodeID, NodeID> > > List_of_Op_RHS;
|
||||||
Equation_Simulation_Type = E_EVALUATE_S;
|
res = equations[eq]->normalizeEquation(var, List_of_Op_RHS);
|
||||||
|
if(mfs==2)
|
||||||
|
{
|
||||||
|
if(d_endo_variable == result.end() && res.second)
|
||||||
|
Equation_Simulation_Type = E_EVALUATE_S;
|
||||||
|
}
|
||||||
|
else if(mfs==3)
|
||||||
|
{
|
||||||
|
if(res.second) // The equation could be solved analytically
|
||||||
|
Equation_Simulation_Type = E_EVALUATE_S;
|
||||||
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
//cout << " " << c_Equation_Type(Equation_Simulation_Type) << endl;
|
|
||||||
V_Equation_Simulation_Type[eq] = make_pair(Equation_Simulation_Type, dynamic_cast<BinaryOpNode *>(res.second));
|
V_Equation_Simulation_Type[eq] = make_pair(Equation_Simulation_Type, dynamic_cast<BinaryOpNode *>(res.second));
|
||||||
}
|
}
|
||||||
return (V_Equation_Simulation_Type);
|
return (V_Equation_Simulation_Type);
|
||||||
|
@ -1000,10 +1005,10 @@ BlockTriangular::Normalize_and_BlockDecompose(bool *IM, Model_Block *ModelBlock,
|
||||||
Compute_Normalization(IM, n, prologue, epilogue, 2, IM_0, Index_Equ_IM);
|
Compute_Normalization(IM, n, prologue, epilogue, 2, IM_0, Index_Equ_IM);
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
SIM0 = (bool *) malloc(n * n * sizeof(bool));
|
SIM0 = (bool *) malloc(n * n * sizeof(bool));
|
||||||
memcpy(SIM0, IM_0, n*n*sizeof(bool));
|
memcpy(SIM0, IM_0, n*n*sizeof(bool));
|
||||||
Prologue_Epilogue(IM, prologue, epilogue, n, Index_Var_IM, Index_Equ_IM, SIM0);
|
Prologue_Epilogue(IM, prologue, epilogue, n, Index_Var_IM, Index_Equ_IM, SIM0);
|
||||||
|
|
||||||
free(SIM0);
|
free(SIM0);
|
||||||
|
|
||||||
V_Equation_Type = Equation_Type_determination(equations, first_order_endo_derivatives, Index_Var_IM, Index_Equ_IM, mfs);
|
V_Equation_Type = Equation_Type_determination(equations, first_order_endo_derivatives, Index_Var_IM, Index_Equ_IM, mfs);
|
||||||
|
|
|
@ -2424,12 +2424,16 @@ DynamicModel::toStaticDll(StaticDllModel &static_model) const
|
||||||
// Convert model local variables (need to be done first)
|
// Convert model local variables (need to be done first)
|
||||||
for (map<int, NodeID>::const_iterator it = local_variables_table.begin();
|
for (map<int, NodeID>::const_iterator it = local_variables_table.begin();
|
||||||
it != local_variables_table.end(); it++)
|
it != local_variables_table.end(); it++)
|
||||||
|
{
|
||||||
static_model.AddLocalVariable(symbol_table.getName(it->first), it->second->toStatic(static_model));
|
static_model.AddLocalVariable(symbol_table.getName(it->first), it->second->toStatic(static_model));
|
||||||
|
}
|
||||||
|
|
||||||
// Convert equations
|
// Convert equations
|
||||||
for (vector<BinaryOpNode *>::const_iterator it = equations.begin();
|
for (vector<BinaryOpNode *>::const_iterator it = equations.begin();
|
||||||
it != equations.end(); it++)
|
it != equations.end(); it++)
|
||||||
|
{
|
||||||
static_model.addEquation((*it)->toStatic(static_model));
|
static_model.addEquation((*it)->toStatic(static_model));
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
void
|
void
|
||||||
|
|
|
@ -52,12 +52,10 @@ ExprNode::getDerivative(int deriv_id)
|
||||||
{
|
{
|
||||||
if (!preparedForDerivation)
|
if (!preparedForDerivation)
|
||||||
prepareForDerivation();
|
prepareForDerivation();
|
||||||
|
|
||||||
// Return zero if derivative is necessarily null (using symbolic a priori)
|
// Return zero if derivative is necessarily null (using symbolic a priori)
|
||||||
set<int>::const_iterator it = non_null_derivatives.find(deriv_id);
|
set<int>::const_iterator it = non_null_derivatives.find(deriv_id);
|
||||||
if (it == non_null_derivatives.end())
|
if (it == non_null_derivatives.end())
|
||||||
return datatree.Zero;
|
return datatree.Zero;
|
||||||
|
|
||||||
// If derivative is stored in cache, use the cached value, otherwise compute it (and cache it)
|
// If derivative is stored in cache, use the cached value, otherwise compute it (and cache it)
|
||||||
map<int, NodeID>::const_iterator it2 = derivatives.find(deriv_id);
|
map<int, NodeID>::const_iterator it2 = derivatives.find(deriv_id);
|
||||||
if (it2 != derivatives.end())
|
if (it2 != derivatives.end())
|
||||||
|
@ -308,7 +306,7 @@ VariableNode::prepareForDerivation()
|
||||||
{
|
{
|
||||||
if (preparedForDerivation)
|
if (preparedForDerivation)
|
||||||
return;
|
return;
|
||||||
|
int der_id;
|
||||||
preparedForDerivation = true;
|
preparedForDerivation = true;
|
||||||
|
|
||||||
// Fill in non_null_derivatives
|
// Fill in non_null_derivatives
|
||||||
|
@ -319,12 +317,25 @@ VariableNode::prepareForDerivation()
|
||||||
case eExogenousDet:
|
case eExogenousDet:
|
||||||
case eParameter:
|
case eParameter:
|
||||||
// For a variable or a parameter, the only non-null derivative is with respect to itself
|
// For a variable or a parameter, the only non-null derivative is with respect to itself
|
||||||
non_null_derivatives.insert(datatree.getDerivID(symb_id, lag));
|
//non_null_derivatives.insert(datatree.getDerivID(symb_id, lag));
|
||||||
|
der_id = datatree.getDerivID(symb_id, lag);
|
||||||
|
non_null_derivatives = set<int>(&der_id, &der_id+1);
|
||||||
break;
|
break;
|
||||||
case eModelLocalVariable:
|
case eModelLocalVariable:
|
||||||
datatree.local_variables_table[symb_id]->prepareForDerivation();
|
datatree.local_variables_table[symb_id]->prepareForDerivation();
|
||||||
// Non null derivatives are those of the value of the local parameter
|
// Non null derivatives are those of the value of the local parameter
|
||||||
non_null_derivatives = datatree.local_variables_table[symb_id]->non_null_derivatives;
|
non_null_derivatives = datatree.local_variables_table[symb_id]->non_null_derivatives;
|
||||||
|
/*for(set<int>::const_iterator it=datatree.local_variables_table[symb_id]->non_null_derivatives.begin(); it != datatree.local_variables_table[symb_id]->non_null_derivatives.end(); it++)
|
||||||
|
non_null_derivatives.insert(*it);*/
|
||||||
|
//cout << "symb_id=" << symb_id << " " << datatree.symbol_table.getName(symb_id) << " " << datatree.symbol_table.maxID();
|
||||||
|
//non_null_derivatives = datatree.local_variables_table[symb_id]->non_null_derivatives;
|
||||||
|
/*cout << " non_null_derivatives.size()= " << non_null_derivatives.size() << endl;
|
||||||
|
if(non_null_derivatives.size()>0)
|
||||||
|
for(set<int>::const_iterator it=non_null_derivatives.begin(); it != non_null_derivatives.end(); it++)
|
||||||
|
{
|
||||||
|
cout << "symb_id=" << symb_id << " *it = " << *it << "\n" ;
|
||||||
|
//cout << datatree.symbol_table.getName(staticdllmodel.getSymbIDByDerivID(*it)) << "(" << staticdllmodel.getLagByDerivID(*it) << ")\n";
|
||||||
|
}*/
|
||||||
break;
|
break;
|
||||||
case eModFileLocalVariable:
|
case eModFileLocalVariable:
|
||||||
// Such a variable is never derived
|
// Such a variable is never derived
|
||||||
|
@ -349,6 +360,7 @@ VariableNode::computeDerivative(int deriv_id)
|
||||||
else
|
else
|
||||||
return datatree.Zero;
|
return datatree.Zero;
|
||||||
case eModelLocalVariable:
|
case eModelLocalVariable:
|
||||||
|
//datatree.local_variables_table[symb_id]->writeOutput(cout);
|
||||||
return datatree.local_variables_table[symb_id]->getDerivative(deriv_id);
|
return datatree.local_variables_table[symb_id]->getDerivative(deriv_id);
|
||||||
case eModFileLocalVariable:
|
case eModFileLocalVariable:
|
||||||
cerr << "ModFileLocalVariable is not derivable" << endl;
|
cerr << "ModFileLocalVariable is not derivable" << endl;
|
||||||
|
@ -628,7 +640,7 @@ VariableNode::compile(ostream &CompileCode, bool lhs_rhs, const temporary_terms_
|
||||||
break;
|
break;
|
||||||
case eModelLocalVariable:
|
case eModelLocalVariable:
|
||||||
case eModFileLocalVariable:
|
case eModFileLocalVariable:
|
||||||
//cout << "eModelLocalVariable=" << symb_id << "\n";
|
//cout << "eModelLocalVariable=" << tsid << "lhs_rhs=" << lhs_rhs << "\n";
|
||||||
datatree.local_variables_table[symb_id]->compile(CompileCode, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic);
|
datatree.local_variables_table[symb_id]->compile(CompileCode, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic);
|
||||||
break;
|
break;
|
||||||
case eUnknownFunction:
|
case eUnknownFunction:
|
||||||
|
|
|
@ -241,7 +241,7 @@ public:
|
||||||
typedef map<const ExprNode *, const VariableNode *> subst_table_t;
|
typedef map<const ExprNode *, const VariableNode *> subst_table_t;
|
||||||
|
|
||||||
//! Creates auxiliary lead variables corresponding to this expression
|
//! Creates auxiliary lead variables corresponding to this expression
|
||||||
/*!
|
/*!
|
||||||
If maximum endogenous lead >= 3, this method will also create intermediary auxiliary var, and will add the equations of the form aux1 = aux2(+1) to the substitution table.
|
If maximum endogenous lead >= 3, this method will also create intermediary auxiliary var, and will add the equations of the form aux1 = aux2(+1) to the substitution table.
|
||||||
\pre This expression is assumed to have maximum endogenous lead >= 2
|
\pre This expression is assumed to have maximum endogenous lead >= 2
|
||||||
\param[in,out] subst_table The table to which new auxiliary variables and their correspondance will be added
|
\param[in,out] subst_table The table to which new auxiliary variables and their correspondance will be added
|
||||||
|
|
|
@ -53,13 +53,15 @@ ModelTree::writeDerivative(ostream &output, int eq, int symb_id, int lag,
|
||||||
void
|
void
|
||||||
ModelTree::computeJacobian(const set<int> &vars)
|
ModelTree::computeJacobian(const set<int> &vars)
|
||||||
{
|
{
|
||||||
|
int i = 0;
|
||||||
for(set<int>::const_iterator it = vars.begin();
|
for(set<int>::const_iterator it = vars.begin();
|
||||||
it != vars.end(); it++)
|
it != vars.end(); it++, i++)
|
||||||
for (int eq = 0; eq < (int) equations.size(); eq++)
|
for (int eq = 0; eq < (int) equations.size(); eq++)
|
||||||
{
|
{
|
||||||
NodeID d1 = equations[eq]->getDerivative(*it);
|
NodeID d1 = equations[eq]->getDerivative(*it);
|
||||||
if (d1 == Zero)
|
if (d1 == Zero)
|
||||||
continue;
|
continue;
|
||||||
|
//printf("eq=%4d var=%4d i=%4d\n",eq,*it,i);
|
||||||
first_derivatives[make_pair(eq, *it)] = d1;
|
first_derivatives[make_pair(eq, *it)] = d1;
|
||||||
++NNZDerivatives[0];
|
++NNZDerivatives[0];
|
||||||
}
|
}
|
||||||
|
|
|
@ -47,7 +47,8 @@ void
|
||||||
StaticDllModel::compileDerivative(ofstream &code_file, int eq, int symb_id, int lag, map_idx_type &map_idx) const
|
StaticDllModel::compileDerivative(ofstream &code_file, int eq, int symb_id, int lag, map_idx_type &map_idx) const
|
||||||
{
|
{
|
||||||
//first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, getDerivID(symb_id, lag)));
|
//first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, getDerivID(symb_id, lag)));
|
||||||
first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, getDerivID(symbol_table.getID(eEndogenous, symb_id), lag)));
|
//first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, getDerivID(symbol_table.getID(eEndogenous, symb_id), lag)));
|
||||||
|
first_derivatives_type::const_iterator it = first_derivatives.find(make_pair(eq, symbol_table.getID(eEndogenous, symb_id)));
|
||||||
if (it != first_derivatives.end())
|
if (it != first_derivatives.end())
|
||||||
(it->second)->compile(code_file, false, temporary_terms, map_idx, false, false);
|
(it->second)->compile(code_file, false, temporary_terms, map_idx, false, false);
|
||||||
else
|
else
|
||||||
|
@ -645,8 +646,8 @@ StaticDllModel::evaluateJacobian(const eval_context_type &eval_context, jacob_ma
|
||||||
it != first_derivatives.end(); it++)
|
it != first_derivatives.end(); it++)
|
||||||
{
|
{
|
||||||
//cout << "it->first.second=" << it->first.second << " variable_table.getSymbolID(it->first.second)=" << variable_table.getSymbolID(it->first.second) << " Type=" << variable_table.getType(it->first.second) << " eEndogenous=" << eEndogenous << " eExogenous=" << eExogenous << " variable_table.getLag(it->first.second)=" << variable_table.getLag(it->first.second) << "\n";
|
//cout << "it->first.second=" << it->first.second << " variable_table.getSymbolID(it->first.second)=" << variable_table.getSymbolID(it->first.second) << " Type=" << variable_table.getType(it->first.second) << " eEndogenous=" << eEndogenous << " eExogenous=" << eExogenous << " variable_table.getLag(it->first.second)=" << variable_table.getLag(it->first.second) << "\n";
|
||||||
if (getTypeByDerivID(it->first.second) == eEndogenous)
|
/*if (getTypeByDerivID(it->first.second) == eEndogenous)
|
||||||
{
|
{*/
|
||||||
NodeID Id = it->second;
|
NodeID Id = it->second;
|
||||||
double val = 0;
|
double val = 0;
|
||||||
try
|
try
|
||||||
|
@ -655,14 +656,17 @@ StaticDllModel::evaluateJacobian(const eval_context_type &eval_context, jacob_ma
|
||||||
}
|
}
|
||||||
catch (ExprNode::EvalException &e)
|
catch (ExprNode::EvalException &e)
|
||||||
{
|
{
|
||||||
cout << "evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(getSymbIDByDerivID(it->first.second)) << "(" << getLagByDerivID(it->first.second) << ") [" << getSymbIDByDerivID(it->first.second) << "] !" << endl;
|
//cout << "evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(getSymbIDByDerivID(it->first.second)) << "(" << getLagByDerivID(it->first.second) << ") [" << getSymbIDByDerivID(it->first.second) << "] !" << endl;
|
||||||
|
cout << "evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(it->first.second) << "(" << 0 << ") [" << it->first.second << "] !" << endl;
|
||||||
Id->writeOutput(cout, oMatlabStaticModelSparse, temporary_terms);
|
Id->writeOutput(cout, oMatlabStaticModelSparse, temporary_terms);
|
||||||
cout << "\n";
|
cout << "\n";
|
||||||
cerr << "StaticDllModel::evaluateJacobian: evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(getSymbIDByDerivID(it->first.second)) << "(" << getLagByDerivID(it->first.second) << ")!" << endl;
|
//cerr << "StaticDllModel::evaluateJacobian: evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(getSymbIDByDerivID(it->first.second)) << "(" << getLagByDerivID(it->first.second) << ")!" << endl;
|
||||||
|
cerr << "StaticDllModel::evaluateJacobian: evaluation of Jacobian failed for equation " << it->first.first+1 << " and variable " << symbol_table.getName(it->first.second) << "(" << 0 << ")!" << endl;
|
||||||
}
|
}
|
||||||
int eq=it->first.first;
|
int eq=it->first.first;
|
||||||
int var = symbol_table.getTypeSpecificID(getSymbIDByDerivID(it->first.second));///symbol_table.getID(eEndogenous,it->first.second);//variable_table.getSymbolID(it->first.second);
|
//int var = symbol_table.getTypeSpecificID(getSymbIDByDerivID(it->first.second));///symbol_table.getID(eEndogenous,it->first.second);//variable_table.getSymbolID(it->first.second);
|
||||||
int k1 = getLagByDerivID(it->first.second);
|
int var = symbol_table.getTypeSpecificID(it->first.second);///symbol_table.getID(eEndogenous,it->first.second);//variable_table.getSymbolID(it->first.second);
|
||||||
|
int k1 = 0;//getLagByDerivID(it->first.second);
|
||||||
if (a_variable_lag!=k1)
|
if (a_variable_lag!=k1)
|
||||||
{
|
{
|
||||||
IM=block_triangular.incidencematrix.Get_IM(k1, eEndogenous);
|
IM=block_triangular.incidencematrix.Get_IM(k1, eEndogenous);
|
||||||
|
@ -673,22 +677,26 @@ StaticDllModel::evaluateJacobian(const eval_context_type &eval_context, jacob_ma
|
||||||
j++;
|
j++;
|
||||||
(*j_m)[make_pair(eq,var)]+=val;
|
(*j_m)[make_pair(eq,var)]+=val;
|
||||||
}
|
}
|
||||||
if (IM[eq*symbol_table.endo_nbr()+var] && (fabs(val) < cutoff))
|
/*if (IM[eq*symbol_table.endo_nbr()+var] && (fabs(val) < cutoff))
|
||||||
{
|
{
|
||||||
if (block_triangular.bt_verbose)
|
//if (block_triangular.bt_verbose)
|
||||||
cout << "the coefficient related to variable " << var << " with lag " << k1 << " in equation " << eq << " is equal to " << val << " and is set to 0 in the incidence matrix (size=" << symbol_table.endo_nbr() << ")\n";
|
cout << "the coefficient related to variable " << var << " with lag " << k1 << " in equation " << eq << " is equal to " << val << " and is set to 0 in the incidence matrix (size=" << symbol_table.endo_nbr() << ")\n";
|
||||||
block_triangular.incidencematrix.unfill_IM(eq, var, k1, eEndogenous);
|
block_triangular.incidencematrix.unfill_IM(eq, var, k1, eEndogenous);
|
||||||
i++;
|
i++;
|
||||||
}
|
}*/
|
||||||
}
|
/*}*/
|
||||||
}
|
}
|
||||||
//Get ride of the elements of the incidence matrix equal to Zero
|
//Get ride of the elements of the incidence matrix equal to Zero
|
||||||
IM=block_triangular.incidencematrix.Get_IM(0, eEndogenous);
|
IM=block_triangular.incidencematrix.Get_IM(0, eEndogenous);
|
||||||
for (int i=0;i<symbol_table.endo_nbr();i++)
|
/*for (int i=0;i<symbol_table.endo_nbr();i++)
|
||||||
for (int j=0;j<symbol_table.endo_nbr();j++)
|
for (int j=0;j<symbol_table.endo_nbr();j++)
|
||||||
if (IM[i*symbol_table.endo_nbr()+j])
|
if (IM[i*symbol_table.endo_nbr()+j])
|
||||||
if (first_derivatives.find(make_pair(i,getDerivID(symbol_table.getID(eEndogenous, j), 0)))==first_derivatives.end())
|
//if (first_derivatives.find(make_pair(i,getDerivID(symbol_table.getID(eEndogenous, j), 0)))==first_derivatives.end())
|
||||||
block_triangular.incidencematrix.unfill_IM(i, j, 0, eEndogenous);
|
if (first_derivatives.find(make_pair(i,symbol_table.getID(eEndogenous, j)))==first_derivatives.end())
|
||||||
|
{
|
||||||
|
block_triangular.incidencematrix.unfill_IM(i, j, 0, eEndogenous);
|
||||||
|
cout << "eliminating equation " << i << " and variable " << j << " in the incidence matrix\n";
|
||||||
|
}*/
|
||||||
if (i>0)
|
if (i>0)
|
||||||
{
|
{
|
||||||
cout << i << " elements among " << first_derivatives.size() << " 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";
|
||||||
|
@ -711,7 +719,8 @@ StaticDllModel::BlockLinear(Model_Block *ModelBlock)
|
||||||
int eq=ModelBlock->Block_List[j].IM_lead_lag[ll].Equ_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];
|
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.getID(var,0)));
|
//first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,variable_table.getID(var,0)));
|
||||||
first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,getDerivID(symbol_table.getID(eEndogenous, var),0)));
|
//first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,getDerivID(symbol_table.getID(eEndogenous, var),0)));
|
||||||
|
first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,symbol_table.getID(eEndogenous, var)));
|
||||||
if (it!= first_derivatives.end())
|
if (it!= first_derivatives.end())
|
||||||
{
|
{
|
||||||
NodeID Id = it->second;
|
NodeID Id = it->second;
|
||||||
|
@ -741,7 +750,8 @@ StaticDllModel::BlockLinear(Model_Block *ModelBlock)
|
||||||
int eq=ModelBlock->Block_List[j].IM_lead_lag[m].Equ_Index[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 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.getID(var,k1)));
|
//first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,variable_table.getID(var,k1)));
|
||||||
first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,getDerivID(symbol_table.getID(eEndogenous, var),k1)));
|
//first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,getDerivID(symbol_table.getID(eEndogenous, var),k1)));
|
||||||
|
first_derivatives_type::const_iterator it=first_derivatives.find(make_pair(eq,symbol_table.getID(eEndogenous, var)));
|
||||||
NodeID Id = it->second;
|
NodeID Id = it->second;
|
||||||
if (it!= first_derivatives.end())
|
if (it!= first_derivatives.end())
|
||||||
{
|
{
|
||||||
|
@ -775,14 +785,17 @@ StaticDllModel::collect_first_order_derivatives_endogenous()
|
||||||
for (first_derivatives_type::iterator it2 = first_derivatives.begin();
|
for (first_derivatives_type::iterator it2 = first_derivatives.begin();
|
||||||
it2 != first_derivatives.end(); it2++)
|
it2 != first_derivatives.end(); it2++)
|
||||||
{
|
{
|
||||||
if (getTypeByDerivID(it2->first.second)==eEndogenous)
|
//if (getTypeByDerivID(it2->first.second)==eEndogenous)
|
||||||
{
|
/*if (it2->first.second)==eEndogenous)
|
||||||
|
{*/
|
||||||
int eq = it2->first.first;
|
int eq = it2->first.first;
|
||||||
int var=symbol_table.getTypeSpecificID(getSymbIDByDerivID(it2->first.second));
|
//int var=symbol_table.getTypeSpecificID(getSymbIDByDerivID(it2->first.second));
|
||||||
int lag=getLagByDerivID(it2->first.second);
|
int var=symbol_table.getTypeSpecificID(it2->first.second);
|
||||||
|
//int lag=getLagByDerivID(it2->first.second);
|
||||||
|
int lag = 0;
|
||||||
//if (lag==0)
|
//if (lag==0)
|
||||||
endo_derivatives[make_pair(eq, make_pair(var, lag))] = it2->second;
|
endo_derivatives[make_pair(eq, make_pair(var, lag))] = it2->second;
|
||||||
}
|
//}
|
||||||
}
|
}
|
||||||
return endo_derivatives;
|
return endo_derivatives;
|
||||||
}
|
}
|
||||||
|
@ -796,26 +809,27 @@ StaticDllModel::computingPass(const eval_context_type &eval_context, bool no_tmp
|
||||||
|
|
||||||
// Compute derivatives w.r. to all endogenous, and possibly exogenous and exogenous deterministic
|
// Compute derivatives w.r. to all endogenous, and possibly exogenous and exogenous deterministic
|
||||||
set<int> vars;
|
set<int> vars;
|
||||||
for (deriv_id_table_t::const_iterator it = deriv_id_table.begin();
|
/*for (deriv_id_table_t::const_iterator it = deriv_id_table.begin();
|
||||||
it != deriv_id_table.end(); it++)
|
it != deriv_id_table.end(); it++)
|
||||||
{
|
{
|
||||||
SymbolType type = symbol_table.getType(it->first.first);
|
SymbolType type = symbol_table.getType(it->first.first);
|
||||||
if (type == eEndogenous)
|
if (type == eEndogenous)
|
||||||
vars.insert(it->second);
|
vars.insert(it->second);
|
||||||
}
|
}*/
|
||||||
|
for(int i = 0; i < symbol_table.endo_nbr(); i++)
|
||||||
|
vars.insert(symbol_table.getID(eEndogenous, i));
|
||||||
|
|
||||||
// Launch computations
|
// Launch computations
|
||||||
cout << "Computing static model derivatives:" << endl
|
cout << "Computing static model derivatives:" << endl
|
||||||
<< " - order 1" << endl;
|
<< " - order 1" << endl;
|
||||||
computeJacobian(vars);
|
computeJacobian(vars);
|
||||||
//cout << "mode=" << mode << " eSparseDLLMode=" << eSparseDLLMode << " eSparseMode=" << eSparseMode << "\n";
|
//cout << "mode=" << mode << " eSparseDLLMode=" << eSparseDLLMode << " eSparseMode=" << eSparseMode << "\n";
|
||||||
|
|
||||||
BuildIncidenceMatrix();
|
BuildIncidenceMatrix();
|
||||||
|
|
||||||
|
|
||||||
jacob_map j_m;
|
jacob_map j_m;
|
||||||
evaluateJacobian(eval_context, &j_m, true);
|
evaluateJacobian(eval_context, &j_m, true);
|
||||||
|
|
||||||
|
|
||||||
if (block_triangular.bt_verbose)
|
if (block_triangular.bt_verbose)
|
||||||
{
|
{
|
||||||
cout << "The gross incidence matrix \n";
|
cout << "The gross incidence matrix \n";
|
||||||
|
@ -914,9 +928,11 @@ StaticDllModel::computeChainRuleJacobian(Model_Block *ModelBlock)
|
||||||
for(int i = 0; i < ModelBlock->Block_List[blck].Nb_Recursives; i++)
|
for(int i = 0; i < ModelBlock->Block_List[blck].Nb_Recursives; i++)
|
||||||
{
|
{
|
||||||
if (ModelBlock->Block_List[blck].Equation_Type[i] == E_EVALUATE_S)
|
if (ModelBlock->Block_List[blck].Equation_Type[i] == E_EVALUATE_S)
|
||||||
recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = ModelBlock->Block_List[blck].Equation_Normalized[i];
|
//recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = ModelBlock->Block_List[blck].Equation_Normalized[i];
|
||||||
|
recursive_variables[symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i])] = ModelBlock->Block_List[blck].Equation_Normalized[i];
|
||||||
else
|
else
|
||||||
recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = equations[ModelBlock->Block_List[blck].Equation[i]];
|
//recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = equations[ModelBlock->Block_List[blck].Equation[i]];
|
||||||
|
recursive_variables[symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i])] = equations[ModelBlock->Block_List[blck].Equation[i]];
|
||||||
}
|
}
|
||||||
map<pair<pair<int, pair<int, int> >, pair<int, int> >, int> Derivatives = block_triangular.get_Derivatives(ModelBlock, blck);
|
map<pair<pair<int, pair<int, int> >, pair<int, int> >, int> Derivatives = block_triangular.get_Derivatives(ModelBlock, blck);
|
||||||
|
|
||||||
|
@ -934,18 +950,20 @@ StaticDllModel::computeChainRuleJacobian(Model_Block *ModelBlock)
|
||||||
int varr = it_l.second.second;
|
int varr = it_l.second.second;
|
||||||
if(Deriv_type == 0)
|
if(Deriv_type == 0)
|
||||||
{
|
{
|
||||||
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = first_derivatives[make_pair(eqr, getDerivID(symbol_table.getID(eEndogenous, varr), lag))];
|
//first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = first_derivatives[make_pair(eqr, getDerivID(symbol_table.getID(eEndogenous, varr), lag))];
|
||||||
|
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = first_derivatives[make_pair(eqr, symbol_table.getID(eEndogenous, varr))];
|
||||||
}
|
}
|
||||||
else if (Deriv_type == 1)
|
else if (Deriv_type == 1)
|
||||||
{
|
{
|
||||||
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = ModelBlock->Block_List[blck].Equation_Normalized[eq]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
|
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = ModelBlock->Block_List[blck].Equation_Normalized[eq]->getChainRuleDerivative(symbol_table.getID(eEndogenous, varr), recursive_variables);
|
||||||
}
|
}
|
||||||
else if (Deriv_type == 2)
|
else if (Deriv_type == 2)
|
||||||
{
|
{
|
||||||
if(ModelBlock->Block_List[blck].Equation_Type[eq] == E_EVALUATE_S && eq<ModelBlock->Block_List[blck].Nb_Recursives)
|
if(ModelBlock->Block_List[blck].Equation_Type[eq] == E_EVALUATE_S && eq<ModelBlock->Block_List[blck].Nb_Recursives)
|
||||||
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = ModelBlock->Block_List[blck].Equation_Normalized[eq]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
|
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = ModelBlock->Block_List[blck].Equation_Normalized[eq]->getChainRuleDerivative(symbol_table.getID(eEndogenous, varr), recursive_variables);
|
||||||
else
|
else
|
||||||
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = equations[eqr]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
|
//first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = equations[eqr]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), lag), recursive_variables);
|
||||||
|
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, lag))] = equations[eqr]->getChainRuleDerivative(symbol_table.getID(eEndogenous, varr), recursive_variables);
|
||||||
}
|
}
|
||||||
ModelBlock->Block_List[blck].Chain_Rule_Derivatives->push_back(make_pair( make_pair(lag, make_pair(eq, var)), make_pair(eqr, varr)));
|
ModelBlock->Block_List[blck].Chain_Rule_Derivatives->push_back(make_pair( make_pair(lag, make_pair(eq, var)), make_pair(eqr, varr)));
|
||||||
}
|
}
|
||||||
|
@ -957,9 +975,11 @@ StaticDllModel::computeChainRuleJacobian(Model_Block *ModelBlock)
|
||||||
for(int i = 0; i < ModelBlock->Block_List[blck].Nb_Recursives; i++)
|
for(int i = 0; i < ModelBlock->Block_List[blck].Nb_Recursives; i++)
|
||||||
{
|
{
|
||||||
if (ModelBlock->Block_List[blck].Equation_Type[i] == E_EVALUATE_S)
|
if (ModelBlock->Block_List[blck].Equation_Type[i] == E_EVALUATE_S)
|
||||||
recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = ModelBlock->Block_List[blck].Equation_Normalized[i];
|
//recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = ModelBlock->Block_List[blck].Equation_Normalized[i];
|
||||||
|
recursive_variables[symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i])] = ModelBlock->Block_List[blck].Equation_Normalized[i];
|
||||||
else
|
else
|
||||||
recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = equations[ModelBlock->Block_List[blck].Equation[i]];
|
//recursive_variables[getDerivID(symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i]), 0)] = equations[ModelBlock->Block_List[blck].Equation[i]];
|
||||||
|
recursive_variables[symbol_table.getID(eEndogenous, ModelBlock->Block_List[blck].Variable[i])] = equations[ModelBlock->Block_List[blck].Equation[i]];
|
||||||
}
|
}
|
||||||
for(int eq = ModelBlock->Block_List[blck].Nb_Recursives; eq < ModelBlock->Block_List[blck].Size; eq++)
|
for(int eq = ModelBlock->Block_List[blck].Nb_Recursives; eq < ModelBlock->Block_List[blck].Size; eq++)
|
||||||
{
|
{
|
||||||
|
@ -967,7 +987,8 @@ StaticDllModel::computeChainRuleJacobian(Model_Block *ModelBlock)
|
||||||
for(int var = ModelBlock->Block_List[blck].Nb_Recursives; var < ModelBlock->Block_List[blck].Size; var++)
|
for(int var = ModelBlock->Block_List[blck].Nb_Recursives; var < ModelBlock->Block_List[blck].Size; var++)
|
||||||
{
|
{
|
||||||
int varr = ModelBlock->Block_List[blck].Variable[var];
|
int varr = ModelBlock->Block_List[blck].Variable[var];
|
||||||
NodeID d1 = equations[eqr]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), 0), recursive_variables);
|
//NodeID d1 = equations[eqr]->getChainRuleDerivative(getDerivID(symbol_table.getID(eEndogenous, varr), 0), recursive_variables);
|
||||||
|
NodeID d1 = equations[eqr]->getChainRuleDerivative(symbol_table.getID(eEndogenous, varr), recursive_variables);
|
||||||
if (d1 == Zero)
|
if (d1 == Zero)
|
||||||
continue;
|
continue;
|
||||||
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, 0))] = d1;
|
first_chain_rule_derivatives[make_pair(eqr, make_pair(varr, 0))] = d1;
|
||||||
|
|
|
@ -250,6 +250,8 @@ StaticModel::computeNormalization()
|
||||||
|
|
||||||
#if 1
|
#if 1
|
||||||
bool check = checked_edmonds_maximum_cardinality_matching(g, &mate_map[0]);
|
bool check = checked_edmonds_maximum_cardinality_matching(g, &mate_map[0]);
|
||||||
|
for (int i = 0; i < n; i++)
|
||||||
|
printf("match equation %4d with variable %4d \n", mate_map[i] - n, i);
|
||||||
#else // Alternative way to compute normalization, by giving an initial matching using natural normalizations
|
#else // Alternative way to compute normalization, by giving an initial matching using natural normalizations
|
||||||
fill(mate_map.begin(), mate_map.end(), graph_traits<BipartiteGraph>::null_vertex());
|
fill(mate_map.begin(), mate_map.end(), graph_traits<BipartiteGraph>::null_vertex());
|
||||||
|
|
||||||
|
@ -404,11 +406,11 @@ StaticModel::computeSortedBlockDecomposition()
|
||||||
for(int i = 0; i < n; i++)
|
for(int i = 0; i < n; i++)
|
||||||
blocks[unordered2ordered[endo2block[i]]].insert(i);
|
blocks[unordered2ordered[endo2block[i]]].insert(i);
|
||||||
|
|
||||||
#ifdef DEBUG
|
//#ifdef DEBUG
|
||||||
cout << "Found " << m << " blocks" << endl;
|
cout << "Found " << m << " blocks" << endl;
|
||||||
for(int i = 0; i < m; i++)
|
for(int i = 0; i < m; i++)
|
||||||
cout << " Block " << i << " of size " << blocks[i].size() << endl;
|
cout << " Block " << i << " of size " << blocks[i].size() << endl;
|
||||||
#endif
|
//#endif
|
||||||
}
|
}
|
||||||
|
|
||||||
void
|
void
|
||||||
|
|
Loading…
Reference in New Issue