stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

New option: "evaluate' to evaluate the static or the dynamic model 

git-svn-id: https://www.dynare.org/svn/dynare/trunk@2919 ac1d8469-bf42-47a9-8791-bf33cf982152
time-shift
ferhat 2009-09-11 17:06:54 +00:00
parent e627f62441
commit 153c7b2542
3 changed files with 324 additions and 83 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

297
mex/sources/bytecode/Interpreter.cc
View File

@ -128,7 +128,7 @@ Interpreter::log1(double a)
void
Interpreter::compute_block_time(int Per_u_) /*throw(EvalException)*/
Interpreter::compute_block_time(int Per_u_, bool evaluate) /*throw(EvalException)*/
{
int var, lag, op;
ostringstream tmp_out;
@ -154,7 +154,10 @@ Interpreter::compute_block_time(int Per_u_) /*throw(EvalException)*/
case eEndogenous :
var=get_code_int;
lag=get_code_int;
Stack.push(y[(it_+lag)*y_size+var]);
if(evaluate)
Stack.push(ya[(it_+lag)*y_size+var]);
else
Stack.push(y[(it_+lag)*y_size+var]);
#ifdef DEBUG
tmp_out << " y[" << it_+lag << ", " << var << "](" << y[(it_+lag)*y_size+var] << ")";
#endif
@ -189,7 +192,10 @@ Interpreter::compute_block_time(int Per_u_) /*throw(EvalException)*/
break;
case eEndogenous :
var=get_code_int;
Stack.push(y[var]);
if(evaluate)
Stack.push(ya[var]);
else
Stack.push(y[var]);
#ifdef DEBUG
tmp_out << " y[" << var << "](" << y[var] << ")";
#endif
@ -640,6 +646,178 @@ Interpreter::compute_block_time(int Per_u_) /*throw(EvalException)*/
}
}
void
Interpreter::evaluate_a_block(int size, int type, string bin_basename, bool steady_state, int block_num)
{
char *begining;
int max_lag_plus_max_lead_plus_1;
int symbol_table_endo_nbr;
int Block_List_Max_Lag;
int Block_List_Max_Lead;
int u_count_int;
bool is_linear;
//mexPrintf("Block = %d type = %d, size = %d, steady_state = %d\n",block_num, type, size, steady_state);
//mexPrintf("y_kmin=%d periods=%d\n",y_kmin, periods);
//mexEvalString("drawnow;");
if(type == SOLVE_FORWARD_COMPLETE || type == SOLVE_BACKWARD_COMPLETE || type == SOLVE_TWO_BOUNDARIES_SIMPLE || type == SOLVE_TWO_BOUNDARIES_COMPLETE)
{
is_linear=get_code_bool;
max_lag_plus_max_lead_plus_1=get_code_int;
symbol_table_endo_nbr=get_code_int;
Block_List_Max_Lag=get_code_int;
Block_List_Max_Lead=get_code_int;
u_count_int=get_code_int;
}
switch (type)
{
case EVALUATE_FORWARD :
if(steady_state)
compute_block_time(0, true);
else
{
begining=get_code_pointer;
for (it_=y_kmin;it_<periods+y_kmin;it_++)
{
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0, true);
}
}
break;
case SOLVE_FORWARD_SIMPLE :
g1=(double*)mxMalloc(size*size*sizeof(double));
r=(double*)mxMalloc(size*sizeof(double));
if(steady_state)
{
compute_block_time(0, true);
for(int j=0; j<size; j++)
y[Block_Contain[j].Variable] += r[j];
}
else
{
begining=get_code_pointer;
for (it_=y_kmin;it_<periods+y_kmin;it_++)
{
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0, true);
for(int j=0; j<size; j++)
y[it_*y_size+Block_Contain[j].Variable] += r[j];
}
}
mxFree(g1);
mxFree(r);
break;
case SOLVE_FORWARD_COMPLETE :
fixe_u(&u, u_count_int, u_count_int);
Read_SparseMatrix(bin_basename, size, 1, 0, 0, steady_state, false);
r=(double*)mxMalloc(size*sizeof(double));
if(steady_state)
{
compute_block_time(0, true);
for(int j=0; j<size; j++)
y[Block_Contain[j].Variable] += r[j];
}
else
{
begining=get_code_pointer;
for (it_=y_kmin;it_<periods+y_kmin;it_++)
{
//mexPrintf("it_=%d, y_size=%d\n",it_, y_size);
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0, true);
for(int j=0; j<size; j++)
y[it_*y_size+Block_Contain[j].Variable] += r[j];
}
}
mxFree(r);
break;
case EVALUATE_BACKWARD :
if(steady_state)
compute_block_time(0, true);
else
{
begining=get_code_pointer;
for (it_=periods+y_kmin-1;it_>=y_kmin;it_--)
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0, true);
}
}
break;
case SOLVE_BACKWARD_SIMPLE :
g1=(double*)mxMalloc(size*size*sizeof(double));
r=(double*)mxMalloc(size*sizeof(double));
if(steady_state)
{
compute_block_time(0, true);
for(int j=0; j<size; j++)
y[Block_Contain[j].Variable] += r[j];
}
else
{
begining=get_code_pointer;
for (it_=periods+y_kmin-1;it_>=y_kmin;it_--)
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0,true);
for(int j=0; j<size; j++)
y[it_*y_size+Block_Contain[j].Variable] += r[j];
}
}
mxFree(g1);
mxFree(r);
break;
case SOLVE_BACKWARD_COMPLETE :
fixe_u(&u, u_count_int, u_count_int);
Read_SparseMatrix(bin_basename, size, 1, 0, 0, steady_state, false);
r=(double*)mxMalloc(size*sizeof(double));
if(steady_state)
{
compute_block_time(0, true);
for(int j=0; j<size; j++)
y[Block_Contain[j].Variable] += r[j];
}
else
{
begining=get_code_pointer;
for (it_=periods+y_kmin-1;it_>=y_kmin;it_--)
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0,true);
for(int j=0; j<size; j++)
y[it_*y_size+Block_Contain[j].Variable] += r[j];
}
}
mxFree(r);
break;
case SOLVE_TWO_BOUNDARIES_SIMPLE :
case SOLVE_TWO_BOUNDARIES_COMPLETE:
fixe_u(&u, u_count_int, u_count_int);
Read_SparseMatrix(bin_basename, size, periods, y_kmin, y_kmax, steady_state, true);
u_count=u_count_int*(periods+y_kmax+y_kmin);
r=(double*)mxMalloc(size*sizeof(double));
begining=get_code_pointer;
for (it_=y_kmin;it_<periods+y_kmin;it_++)
{
Per_u_=(it_-y_kmin)*u_count_int;
Per_y_=it_*y_size;
set_code_pointer(begining);
compute_block_time(Per_u_, true);
for(int j=0; j<size; j++)
y[it_*y_size+Block_Contain[j].Variable] += r[j];
}
mxFree(r);
break;
}
}
bool
Interpreter::simulate_a_block(int size,int type, string file_name, string bin_basename, bool Gaussian_Elimination, bool steady_state, int block_num)
{
@ -659,7 +837,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
{
case EVALUATE_FORWARD :
if(steady_state)
compute_block_time(0);
compute_block_time(0, false);
else
{
begining=get_code_pointer;
@ -667,13 +845,13 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0);
compute_block_time(0, false);
}
}
break;
case EVALUATE_BACKWARD :
if(steady_state)
compute_block_time(0);
compute_block_time(0, false);
else
{
begining=get_code_pointer;
@ -681,7 +859,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0);
compute_block_time(0, false);
}
}
break;
@ -697,50 +875,50 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0);
y[Block_Contain[0].Variable] += -r[0]/g1[0];
compute_block_time(0, false);
y[Block_Contain[0].Variable] += -r[0]/g1[0];
double rr;
rr=r[0];
cvg=((rr*rr)<solve_tolf);
iter++;
}
if (!cvg)
{
mexPrintf("Convergence not achieved in block %d, after %d iterations\n",Block_Count,iter);
/*mexEvalString("st=fclose('all');clear all;");
mexErrMsgTxt("End of simulate");*/
return false;
}
}
else
{
for (it_=y_kmin;it_<periods+y_kmin;it_++)
{
cvg=false;
iter=0;
Per_y_=it_*y_size;
while (!(cvg||(iter>maxit_)))
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0, false);
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
double rr;
rr=r[0];
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
rr=r[0]/(1+y[Per_y_+Block_Contain[0].Variable]);
else
rr=r[0];
cvg=((rr*rr)<solve_tolf);
iter++;
}
if (!cvg)
{
mexPrintf("Convergence not achieved in block %d, after %d iterations\n",Block_Count,iter);
/*mexEvalString("st=fclose('all');clear all;");
mexErrMsgTxt("End of simulate");*/
return false;
mexPrintf("Convergence not achieved in block %d, at time %d after %d iterations\n",Block_Count,it_,iter);
mexEvalString("st=fclose('all');clear all;");
mexErrMsgTxt("End of simulate");
}
}
else
{
for (it_=y_kmin;it_<periods+y_kmin;it_++)
{
cvg=false;
iter=0;
Per_y_=it_*y_size;
while (!(cvg||(iter>maxit_)))
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0);
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
double rr;
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
rr=r[0]/(1+y[Per_y_+Block_Contain[0].Variable]);
else
rr=r[0];
cvg=((rr*rr)<solve_tolf);
iter++;
}
if (!cvg)
{
mexPrintf("Convergence not achieved in block %d, at time %d after %d iterations\n",Block_Count,it_,iter);
mexEvalString("st=fclose('all');clear all;");
mexErrMsgTxt("End of simulate");
}
}
}
}
}
mxFree(g1);
mxFree(r);
break;
@ -756,7 +934,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0);
compute_block_time(0, false);
y[Block_Contain[0].Variable] += -r[0]/g1[0];
double rr;
rr=r[0];
@ -782,7 +960,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
{
set_code_pointer(begining);
Per_y_=it_*y_size;
compute_block_time(0);
compute_block_time(0, false);
y[Per_y_+Block_Contain[0].Variable] += -r[0]/g1[0];
double rr;
if(fabs(1+y[Per_y_+Block_Contain[0].Variable])>eps)
@ -833,7 +1011,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
res2=0;
res1=0;
max_res=0;
compute_block_time(0);
compute_block_time(0, false);
/*if (isnan(res1)||isinf(res1))
{
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
@ -875,7 +1053,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
iter = 0;
res1=res2=max_res=0;max_res_idx=0;
error_not_printed = true;
compute_block_time(0);
compute_block_time(0, false);
cvg=false;
result = simulate_NG(Block_Count, symbol_table_endo_nbr, 0, 0, 0, size, false, cvg, iter, true);
}
@ -899,7 +1077,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
res2=0;
res1=0;
max_res=0;
compute_block_time(0);
compute_block_time(0, false);
/*if (isnan(res1)||isinf(res1))
{
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
@ -945,7 +1123,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
iter = 0;
res1=res2=max_res=0;max_res_idx=0;
error_not_printed = true;
compute_block_time(0);
compute_block_time(0, false);
cvg=false;
result = simulate_NG(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, false, cvg, iter, false);
}
@ -987,7 +1165,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
res2=0;
res1=0;
max_res=0;
compute_block_time(0);
compute_block_time(0, false);
/*if (isnan(res1)||isinf(res1))
{
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
@ -1028,7 +1206,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
iter = 0;
res1=res2=max_res=0;max_res_idx=0;
error_not_printed = true;
compute_block_time(0);
compute_block_time(0, false);
cvg=false;
result = simulate_NG(Block_Count, symbol_table_endo_nbr, 0, 0, 0, size, false, cvg, iter, true);
}
@ -1050,7 +1228,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
res2=0;
res1=0;
max_res=0;
compute_block_time(0);
compute_block_time(0, false);
/*if (isnan(res1)||isinf(res1))
{
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
@ -1094,7 +1272,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
set_code_pointer(begining);
Per_y_=it_*y_size;
error_not_printed = true;
compute_block_time(0);
compute_block_time(0, false);
cvg=false;
result = simulate_NG(Block_Count, symbol_table_endo_nbr, it_, y_kmin, y_kmax, size, false, cvg, iter, false);
}
@ -1147,7 +1325,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
Per_u_=(it_-y_kmin)*u_count_int;
Per_y_=it_*y_size;
set_code_pointer(begining);
compute_block_time(Per_u_);
compute_block_time(Per_u_, false);
if (isnan(res1)||isinf(res1))
{
memcpy(y, y_save, y_size*sizeof(double)*(periods+y_kmax+y_kmin));
@ -1192,7 +1370,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
Per_u_=(it_-y_kmin)*u_count_int;
Per_y_=it_*y_size;
set_code_pointer(begining);
compute_block_time(Per_u_);
compute_block_time(Per_u_, false);
for (i=0; i< size; i++)
{
double rr;
@ -1226,7 +1404,7 @@ Interpreter::simulate_a_block(int size,int type, string file_name, string bin_ba
}
bool
Interpreter::compute_blocks(string file_name, string bin_basename, bool steady_state)
Interpreter::compute_blocks(string file_name, string bin_basename, bool steady_state, bool evaluate)
{
ifstream CompiledCode;
bool result = true;
@ -1283,7 +1461,10 @@ Interpreter::compute_blocks(string file_name, string bin_basename, bool steady_s
lBlock_Contain.Own_Derivative=get_code_int;
Block_Contain.push_back(lBlock_Contain);
}
result = simulate_a_block(lBlock.size,lBlock.type, file_name, bin_basename,true, steady_state, Block_Count);
if(evaluate)
evaluate_a_block(lBlock.size,lBlock.type, bin_basename, steady_state, Block_Count);
else
result = simulate_a_block(lBlock.size,lBlock.type, file_name, bin_basename,true, steady_state, Block_Count);
if(!result)
go_on = false;
break;

7
mex/sources/bytecode/Interpreter.hh
View File

@ -68,8 +68,9 @@ class Interpreter : SparseMatrix
protected :
double pow1(double a, double b);
double log1(double a);
void compute_block_time(int Per_u_);
bool simulate_a_block(int size,int type, string file_name, string bin_basename, bool Gaussian_Elimination, bool steady_state, int block_num);
void compute_block_time(int Per_u_, bool evaluate);
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);
double *T;
vector<Block_contain_type> Block_Contain;
vector<Block_type> Block;
@ -92,7 +93,7 @@ class Interpreter : SparseMatrix
Interpreter(double *params_arg, double *y_arg, double *ya_arg, double *x_arg, double *direction_arg, int y_size_arg, int nb_row_x_arg,
int nb_row_xd_arg, int periods_arg, int y_kmin_arg, int y_kmax_arg, int maxit_arg_, double solve_tolf_arg, int size_o_direction_arg,
double slowc_arg, int y_decal_arg, double markowitz_c_arg, string &filename_arg);
bool compute_blocks(string file_name, string bin_basename, bool steady_state);
bool compute_blocks(string file_name, string bin_basename, bool steady_state, bool evaluate);
};

103
mex/sources/bytecode/bytecode.cc
View File

@ -46,12 +46,21 @@ max(int a, int b)
using namespace std;
#include <sstream>
#include "mex_interface.hh"
string
Get_Argument(const char *argv)
{
//mexPrintf("number=%d\n",number);
string f(argv);
return f;
}
int
main( int argc, const char* argv[] )
{
FILE *fid;
bool steady_state = false;
bool evaluate = false;
printf("argc=%d\n",argc);
if(argc<2)
{
@ -70,11 +79,23 @@ main( int argc, const char* argv[] )
string file_name(argv[1]);
if(argc>2)
for(i=2;i<argc; i++)
{
string f(argv[1]);
if(f == "steady_state")
if(Get_Argument(argv[i])=="steady_state")
steady_state = true;
else if(Get_Argument(argv[i])=="dynamic")
steady_state = false;
else if(Get_Argument(argv[i])=="evaluate")
evaluate = true;
else
{
mexPrintf("Unknown argument : ");
mexEvalString("st=fclose('all');clear all;");
string f;
f = Get_Argument(argv[i]);
f.append("\n");
mexErrMsgTxt(f.c_str());
}
}
fid = fopen(tmp_out.str().c_str(),"r");
int periods;
@ -147,9 +168,10 @@ main( int argc, const char* argv[] )
clock_t t0= clock();
Interpreter interprete(params, y, ya, x, direction, y_size, nb_row_x, nb_row_xd, periods, y_kmin, y_kmax, maxit_, solve_tolf, size_of_direction, slowc, y_decal, markowitz_c, file_name);
string f(file_name);
interprete.compute_blocks(f, f, steady_state);
interprete.compute_blocks(f, f, steady_state, evaluate);
clock_t t1= clock();
mexPrintf("Simulation Time=%f milliseconds\n",1000.0*(double(t1)-double(t0))/double(CLOCKS_PER_SEC));
if(!evaluate)
mexPrintf("Simulation Time=%f milliseconds\n",1000.0*(double(t1)-double(t0))/double(CLOCKS_PER_SEC));
if(x)
mxFree(x);
if(y)
@ -162,28 +184,51 @@ main( int argc, const char* argv[] )
}
#else
string
Get_Argument(const mxArray *prhs)
{
//mexPrintf("number=%d\n",number);
const mxArray *mxa = prhs;
int buflen=mxGetM(mxa) * mxGetN(mxa) + 1;
char *first_argument;
first_argument=(char*)mxCalloc(buflen, sizeof(char));
int status = mxGetString(mxa, first_argument, buflen);
if (status != 0)
mexWarnMsgTxt("Not enough space. The first argument is truncated.");
string f(first_argument);
mxFree(first_argument);
return f;
}
/* The gateway routine */
void
mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
mxArray *M_, *oo_, *options_;
int i, row_y, col_y, row_x, col_x, nb_row_xd;
int steady_row_y, steady_col_y, steady_row_x, steady_col_x, steady_nb_row_xd;
int y_kmin=0, y_kmax=0, y_decal=0, periods=1;
double * pind ;
double *direction;
bool steady_state = false;
if(nrhs>0)
bool evaluate = false;
for(i=0;i<nrhs; i++)
{
const mxArray *mxa = prhs[0];
int buflen=mxGetM(mxa) * mxGetN(mxa) + 1;
char *first_argument;
first_argument=(char*)mxCalloc(buflen, sizeof(char));
int status = mxGetString(mxa, first_argument, buflen);
if (status != 0)
mexWarnMsgTxt("Not enough space. The first argument is truncated.");
string f(first_argument);
if(f == "steady_state")
if(Get_Argument(prhs[i])=="steady_state")
steady_state = true;
else if(Get_Argument(prhs[i])=="dynamic")
steady_state = false;
else if(Get_Argument(prhs[i])=="evaluate")
evaluate = true;
else
{
mexPrintf("Unknown argument : ");
mexEvalString("st=fclose('all');clear all;");
string f;
f = Get_Argument(prhs[i]);
f.append("\n");
mexErrMsgTxt(f.c_str());
}
}
M_ = mexGetVariable("global","M_");
if (M_ == NULL )
@ -209,7 +254,7 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
}
//mexPrintf("ok0\n");
double * params = mxGetPr(mxGetFieldByNumber(M_, 0, mxGetFieldNumber(M_,"params")));
double *yd, *xd;
double *yd, *xd, *steady_yd, *steady_xd ;
if(!steady_state)
{
yd= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"endo_simul")));
@ -224,6 +269,14 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
y_kmax=int(floor(*(mxGetPr(mxGetFieldByNumber(M_, 0, mxGetFieldNumber(M_,"maximum_lead"))))));
y_decal=max(0,y_kmin-int(floor(*(mxGetPr(mxGetFieldByNumber(M_, 0, mxGetFieldNumber(M_,"maximum_endo_lag")))))));
periods=int(floor(*(mxGetPr(mxGetFieldByNumber(options_, 0, mxGetFieldNumber(options_,"periods"))))));
steady_yd= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"steady_state")));
steady_row_y=mxGetM(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"steady_state")));
steady_col_y=mxGetN(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"steady_state")));;
steady_xd= mxGetPr(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"exo_steady_state")));
steady_row_x=mxGetM(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"exo_steady_state")));
steady_col_x=mxGetN(mxGetFieldByNumber(oo_, 0, mxGetFieldNumber(oo_,"exo_steady_state")));
steady_nb_row_xd=int(floor(*(mxGetPr(mxGetFieldByNumber(M_, 0, mxGetFieldNumber(M_,"exo_det_nbr"))))));
}
else
{
@ -262,27 +315,33 @@ mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
for (i=0;i<row_x*col_x;i++)
x[i]=double(xd[i]);
for (i=0;i<row_y*col_y;i++)
y[i]=double(yd[i]);
{
y[i] = double(yd[i]);
ya[i] = double(yd[i]);
}
int y_size=row_y;
int nb_row_x=row_x;
/*int it_ = y_kmin;
for (int j = 0; j < y_size; j++)
mexPrintf(" variable %d at time %d and %d = %f\n", j+1, it_, it_+1, y[j+it_*y_size]);*/
clock_t t0= clock();
Interpreter interprete(params, y, ya, x, direction, y_size, nb_row_x, nb_row_xd, periods, y_kmin, y_kmax, maxit_, solve_tolf, size_of_direction, slowc, y_decal, markowitz_c, file_name);
string f(fname);
bool result = interprete.compute_blocks(f, f, steady_state);
bool result = interprete.compute_blocks(f, f, steady_state, evaluate);
clock_t t1= clock();
if(!steady_state)
if(!steady_state && !evaluate)
mexPrintf("Simulation Time=%f milliseconds\n",1000.0*(double(t1)-double(t0))/double(CLOCKS_PER_SEC));
if (nlhs>0)
{
plhs[0] = mxCreateDoubleMatrix(row_y, col_y, mxREAL);
pind = mxGetPr(plhs[0]);
for (i=0;i<row_y*col_y;i++)
pind[i]=y[i];
if(evaluate)
for (i=0;i<row_y*col_y;i++)
pind[i]=y[i]-ya[i];
else
for (i=0;i<row_y*col_y;i++)
pind[i]=y[i];
if(nlhs>1)
{
plhs[1] = mxCreateDoubleMatrix(1, 1, mxREAL);