Skip to content

GitLab

dynare
dynare
Commit 97664607 authored by Ferhat Mihoubi's avatar Ferhat Mihoubi
Browse files
Options

- external functions are compatible with block decomposition and/or bytecode

parent befa1b96
Expand all Hide whitespace changes
Inline Side-by-side
Showing with 1167 additions and 87 deletions
mex/sources/bytecode/ErrorHandling.hh
View file @ 97664607
......@@ -139,6 +139,9 @@ double *u, *y, *ya;
double *steady_y, *steady_x;
double *g2, *g1, *r;
vector<mxArray*> jacobian_block, jacobian_other_endo_block, jacobian_exo_block, jacobian_det_exo_block;
map<unsigned int,double> TEF;
map<pair<unsigned int, unsigned int>, double > TEFD;
map<pair<unsigned int, pair<unsigned int, unsigned int> >, double > TEFDD;
ExpressionType EQN_type;
it_code_type it_code_expr;
......
mex/sources/bytecode/Interpreter.cc
View file @ 97664607
......@@ -160,6 +160,9 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate, int block_num, int si
#endif
EQN_type = TemporaryTerm;
EQN_equation = ((FNUMEXPR_ *) it_code->second)->get_equation();
#ifdef DEBUG
mexPrintf("EQN_equation=%d\n",EQN_equation);mexEvalString("drawnow;");
#endif
break;
case ModelEquation:
#ifdef DEBUG
......@@ -1058,6 +1061,339 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate, int block_num, int si
}
break;
case FPUSH:
break;
case FCALL:
{
#ifdef DEBUG
mexPrintf("------------------------------\n");
mexPrintf("CALL ");mexEvalString("drawnow;");
#endif
FCALL_ *fc = (FCALL_ *) it_code->second;
string function_name = fc->get_function_name();
#ifdef DEBUG
mexPrintf("function_name=%s ", function_name.c_str());mexEvalString("drawnow;");
#endif
unsigned int nb_input_arguments = fc->get_nb_input_arguments();
#ifdef DEBUG
mexPrintf("nb_input_arguments=%d ", nb_input_arguments);mexEvalString("drawnow;");
#endif
unsigned int nb_output_arguments = fc->get_nb_output_arguments();
#ifdef DEBUG
mexPrintf("nb_output_arguments=%d\n", nb_output_arguments);mexEvalString("drawnow;");
#endif
mxArray *output_arguments[3];
string arg_func_name = fc->get_arg_func_name();
#ifdef DEBUG
mexPrintf("arg_func_name.length() = %d\n",arg_func_name.length());
mexPrintf("arg_func_name.c_str() = %s\n",arg_func_name.c_str());
#endif
unsigned int nb_add_input_arguments = fc->get_nb_add_input_arguments();
external_function_type function_type = fc->get_function_type();
#ifdef DEBUG
mexPrintf("function_type=%d ExternalFunctionWithoutDerivative=%d\n",function_type, ExternalFunctionWithoutDerivative);
mexEvalString("drawnow;");
#endif
mxArray **input_arguments;
switch (function_type)
{
case ExternalFunctionWithoutDerivative:
case ExternalFunctionWithFirstDerivative:
case ExternalFunctionWithFirstandSecondDerivative:
{
input_arguments = (mxArray**)mxMalloc(nb_input_arguments * sizeof(mxArray*));
#ifdef DEBUG
mexPrintf("Stack.size()=%d\n",Stack.size());
mexEvalString("drawnow;");
#endif
for (unsigned int i = 0; i < nb_input_arguments ; i++)
{
mxArray *vv = mxCreateDoubleScalar(Stack.top());
input_arguments[nb_input_arguments - i - 1] = vv;
Stack.pop();
}
mexCallMATLAB(nb_output_arguments, output_arguments, nb_input_arguments, input_arguments, function_name.c_str());
double *rr = mxGetPr(output_arguments[0]);
Stack.push(*rr);
if (function_type == ExternalFunctionWithFirstDerivative || function_type == ExternalFunctionWithFirstandSecondDerivative)
{
unsigned int indx = fc->get_indx();
double *FD1 = mxGetPr(output_arguments[1]);
unsigned int rows = mxGetN(output_arguments[1]);
for(unsigned int i = 0; i < rows; i++)
TEFD[make_pair(indx, i)] = FD1[i];
}
if (function_type == ExternalFunctionWithFirstandSecondDerivative)
{
unsigned int indx = fc->get_indx();
double *FD2 = mxGetPr(output_arguments[2]);
unsigned int rows = mxGetM(output_arguments[2]);
unsigned int cols = mxGetN(output_arguments[2]);
unsigned int k = 0;
for (unsigned int j = 0; j < cols; j++)
for (unsigned int i = 0; i < rows; i++)
TEFDD[make_pair(indx, make_pair(i, j))] = FD2[k++];
}
}
break;
case ExternalFunctionNumericalFirstDerivative:
{
input_arguments = (mxArray**)mxMalloc((nb_input_arguments+1+nb_add_input_arguments) * sizeof(mxArray*));
mxArray *vv = mxCreateString(arg_func_name.c_str());
input_arguments[0] = vv;
vv = mxCreateDoubleScalar(fc->get_row());
input_arguments[1] = vv;
vv = mxCreateCellMatrix(1, nb_add_input_arguments);
for (unsigned int i = 0; i < nb_add_input_arguments; i++)
{
double rr = Stack.top();
#ifdef DEBUG
mexPrintf("i=%d rr = %f Stack.size()=%d\n",i, rr, Stack.size());
#endif
mxSetCell(vv, nb_add_input_arguments - (i+1), mxCreateDoubleScalar(rr));
Stack.pop();
}
input_arguments[nb_input_arguments+nb_add_input_arguments] = vv;
#ifdef DEBUG
mexCallMATLAB(0, NULL, 1, & input_arguments[0], "disp");
mexCallMATLAB(0, NULL, 1, & input_arguments[1], "disp");
mexCallMATLAB(0, NULL, 1, & input_arguments[2], "celldisp");
mexPrintf("OK\n");
mexEvalString("drawnow;");
#endif
nb_input_arguments = 3;
mexCallMATLAB(nb_output_arguments, output_arguments, nb_input_arguments, input_arguments, function_name.c_str());
double *rr = mxGetPr(output_arguments[0]);
#ifdef DEBUG
mexPrintf("*rr=%f\n",*rr);
#endif
Stack.push(*rr);
}
break;
case ExternalFunctionFirstDerivative:
{
input_arguments = (mxArray**)mxMalloc(nb_input_arguments * sizeof(mxArray*));
for (unsigned int i = 0; i < nb_input_arguments; i++)
{
mxArray *vv = mxCreateDoubleScalar(Stack.top());
input_arguments[(nb_input_arguments - 1) - i] = vv;
Stack.pop();
}
mexCallMATLAB(nb_output_arguments, output_arguments, nb_input_arguments, input_arguments, function_name.c_str());
unsigned int indx = fc->get_indx();
double *FD1 = mxGetPr(output_arguments[0]);
//mexPrint
unsigned int rows = mxGetN(output_arguments[0]);
for (unsigned int i = 0; i < rows; i++)
TEFD[make_pair(indx, i)] = FD1[i];
}
break;
case ExternalFunctionNumericalSecondDerivative:
{
input_arguments = (mxArray**)mxMalloc((nb_input_arguments+1+nb_add_input_arguments) * sizeof(mxArray*));
mxArray *vv = mxCreateString(arg_func_name.c_str());
input_arguments[0] = vv;
vv = mxCreateDoubleScalar(fc->get_row());
input_arguments[1] = vv;
vv = mxCreateDoubleScalar(fc->get_col());
input_arguments[2] = vv;
vv = mxCreateCellMatrix(1, nb_add_input_arguments);
for (unsigned int i = 0; i < nb_add_input_arguments; i++)
{
double rr = Stack.top();
mexPrintf("i=%d rr = %f\n",i, rr);
mxSetCell(vv, (nb_add_input_arguments - 1) - i, mxCreateDoubleScalar(rr));
Stack.pop();
}
input_arguments[nb_input_arguments+nb_add_input_arguments] = vv;
mexCallMATLAB(0, NULL, 1, & input_arguments[0], "disp");
mexCallMATLAB(0, NULL, 1, & input_arguments[1], "disp");
mexCallMATLAB(0, NULL, 1, & input_arguments[2], "celldisp");
mexPrintf("OK\n");
mexEvalString("drawnow;");
nb_input_arguments = 3;
mexCallMATLAB(nb_output_arguments, output_arguments, nb_input_arguments, input_arguments, function_name.c_str());
double *rr = mxGetPr(output_arguments[0]);
Stack.push(*rr);
}
break;
case ExternalFunctionSecondDerivative:
{
input_arguments = (mxArray**)mxMalloc(nb_input_arguments * sizeof(mxArray*));
for (unsigned int i = 0; i < nb_input_arguments ; i++)
{
mxArray *vv = mxCreateDoubleScalar(Stack.top());
input_arguments[i] = vv;
Stack.pop();
}
mexCallMATLAB(nb_output_arguments, output_arguments, nb_input_arguments, input_arguments, function_name.c_str());
unsigned int indx = fc->get_indx();
double *FD2 = mxGetPr(output_arguments[2]);
unsigned int rows = mxGetM(output_arguments[0]);
unsigned int cols = mxGetN(output_arguments[0]);
unsigned int k = 0;
for (unsigned int j = 0; j < cols; j++)
for (unsigned int i = 0; i < rows; i++)
TEFDD[make_pair(indx, make_pair(i, j))] = FD2[k++];
}
break;
}
/*if (f)
{
input_arguments = (mxArray**)mxMalloc((nb_input_arguments+1+nb_add_input_arguments) * sizeof(mxArray*));
//the called function is jacob_element or hess_element
mxArray *vv = mxCreateString(arg_func_name.c_str());
input_arguments[0] = vv;
vv = mxCreateDoubleScalar(fc->get_indx());
input_arguments[1] = vv;
start_input_arg += 2;
}
else
input_arguments = (mxArray**)mxMalloc(nb_input_arguments * sizeof(mxArray*));
for (unsigned int i = start_input_arg; i < nb_input_arguments + start_input_arg; i++)
{
mxArray *vv = mxCreateDoubleScalar(Stack.top());
input_arguments[i] = vv;
Stack.pop();
}
mexPrintf("nb_add_input_arguments=%d Stack.size()=%d\n",nb_add_input_arguments, Stack.size());mexEvalString("drawnow;");
if (arg_func_name.length() > 0)
{
mxArray *vv = mxCreateCellArray(1, &nb_add_input_arguments);
for (unsigned int i = 0; i < nb_add_input_arguments; i++)
{
double rr = Stack.top();
mexPrintf("i=%d rr = %f\n",i, rr);
mxSetCell(vv, i, mxCreateDoubleScalar(rr));
Stack.pop();
}
input_arguments[nb_input_arguments+nb_add_input_arguments] = vv;
mexCallMATLAB(0, NULL, 1, & input_arguments[0], "disp");
mexCallMATLAB(0, NULL, 1, & input_arguments[1], "disp");
mexCallMATLAB(0, NULL, 1, & input_arguments[2], "celldisp");
mexPrintf("OK\n");
mexEvalString("drawnow;");
nb_input_arguments = 3;
}
mexCallMATLAB(nb_output_arguments, output_arguments, nb_input_arguments, input_arguments, function_name.c_str());
double *rr = mxGetPr(output_arguments[0]);
Stack.push(*rr);
if (nb_output_arguments >= 2) //its the return of a TEF
{
unsigned int indx = fc->get_indx();
double *FD1 = mxGetPr(output_arguments[1]);
unsigned int rows = mxGetM(output_arguments[1]);
for(unsigned int i = 0; i < rows; i++)
TEFD[make_pair(indx, i)] = FD1[i];
if (nb_output_arguments == 3)
{
double *FD2 = mxGetPr(output_arguments[2]);
unsigned int rows = mxGetM(output_arguments[2]);
unsigned int cols = mxGetN(output_arguments[2]);
unsigned int k = 0;
for (unsigned int j = 0; j < cols; j++)
for (unsigned int i = 0; i < rows; i++)
TEFDD[make_pair(indx, make_pair(i, j))] = FD2[k++];
}
}
else*/
/*
#ifdef DEBUG
mexPrintf("Stack.size()=%d, *rr=%f\n",Stack.size(), *rr);
mexPrintf("done\n");
mexEvalString("drawnow;");
#endif*/
}
break;
case FSTPTEF:
var = ((FSTPTEF_ *) it_code->second)->get_number();
#ifdef DEBUG
mexPrintf("FSTPTEF\n");
mexPrintf("var=%d Stack.size()=%d\n",var, Stack.size());
#endif
TEF[var-1] = Stack.top();
#ifdef DEBUG
mexPrintf("FSTP TEF[var-1]=%f done\n",TEF[var-1]);
mexEvalString("drawnow;");
#endif
Stack.pop();
break;
case FLDTEF:
var = ((FLDTEF_ *) it_code->second)->get_number();
#ifdef DEBUG
mexPrintf("FLDTEF\n");
mexPrintf("var=%d Stack.size()=%d\n",var, Stack.size());
mexPrintf("FLD TEF[var-1]=%f done\n",TEF[var-1]);
mexEvalString("drawnow;");
#endif
Stack.push(TEF[var-1]);
break;
case FSTPTEFD:
{
unsigned int indx = ((FSTPTEFD_ *) it_code->second)->get_indx();
unsigned int row = ((FSTPTEFD_ *) it_code->second)->get_row();
#ifdef DEBUG
mexPrintf("FSTPTEFD\n");
mexPrintf("indx=%d Stack.size()=%d\n",indx, Stack.size());
#endif
TEFD[make_pair(indx, row-1)] = Stack.top();
#ifdef DEBUG
mexPrintf("FSTP TEFD[make_pair(indx, row)]=%f done\n",TEFD[make_pair(indx, row-1)]);
mexEvalString("drawnow;");
#endif
Stack.pop();
}
break;
case FLDTEFD:
{
unsigned int indx = ((FLDTEFD_ *) it_code->second)->get_indx();
unsigned int row = ((FLDTEFD_ *) it_code->second)->get_row();
#ifdef DEBUG
mexPrintf("FLDTEFD\n");
mexPrintf("indx=%d row=%d Stack.size()=%d\n",indx, row, Stack.size());
mexPrintf("FLD TEFD[make_pair(indx, row)]=%f done\n",TEFD[make_pair(indx, row-1)]);
mexEvalString("drawnow;");
#endif
Stack.push(TEFD[make_pair(indx, row-1)]);
}
break;
case FSTPTEFDD:
{
unsigned int indx = ((FSTPTEFDD_ *) it_code->second)->get_indx();
unsigned int row = ((FSTPTEFDD_ *) it_code->second)->get_row();
unsigned int col = ((FSTPTEFDD_ *) it_code->second)->get_col();
#ifdef DEBUG
mexPrintf("FSTPTEFD\n");
mexPrintf("indx=%d Stack.size()=%d\n",indx, Stack.size());
#endif
TEFDD[make_pair(indx, make_pair(row-1, col-1))] = Stack.top();
#ifdef DEBUG
mexPrintf("FSTP TEFDD[make_pair(indx, make_pair(row, col))]=%f done\n",TEFDD[make_pair(indx, make_pair(row, col))]);
mexEvalString("drawnow;");
#endif
Stack.pop();
}
break;
case FLDTEFDD:
{
unsigned int indx = ((FLDTEFDD_ *) it_code->second)->get_indx();
unsigned int row = ((FLDTEFDD_ *) it_code->second)->get_row();
unsigned int col = ((FSTPTEFDD_ *) it_code->second)->get_col();
#ifdef DEBUG
mexPrintf("FLDTEFD\n");
mexPrintf("indx=%d Stack.size()=%d\n",indx, Stack.size());
mexPrintf("FLD TEFD[make_pair(indx, make_pair(row, col))]=%f done\n",TEFDD[make_pair(indx, make_pair(row, col))]);
mexEvalString("drawnow;");
#endif
Stack.push(TEFDD[make_pair(indx, make_pair(row-1, col-1))]);
}
break;
case FCUML:
v1 = Stack.top();
Stack.pop();
......
preprocessor/CodeInterpreter.hh
View file @ 97664607
......@@ -93,7 +93,17 @@ enum Tags
FOK, //!< Used for debugging purpose - 21 (33)
FNUMEXPR //!< Store the expression type and references - 22 (34)
FNUMEXPR, //!< Store the expression type and references - 22 (34)
FCALL, //!< Call an external function - 23 (35)
FPUSH, //!< Push a double in the stack - 24 (36)
FPOP, //!< Pop a double from the stack - 25 (37)
FLDTEF, //!< Stores the result of an external function in the stack - 26 (38)
FSTPTEF, //!< Loads the result of an external function from the stack- 27 (39)
FLDTEFD, //!< Stores the result of an external function in the stack - 28 (40)
FSTPTEFD, //!< Loads the result of an external function from the stack- 29 (41)
FLDTEFDD, //!< Stores the result of an external function in the stack - 28 (42)
FSTPTEFDD //!< Loads the result of an external function from the stack- 29 (43)
};
......@@ -207,6 +217,17 @@ enum TrinaryOpcode
oNormpdf
};
enum external_function_type
{
ExternalFunctionWithoutDerivative,
ExternalFunctionWithFirstDerivative,
ExternalFunctionWithFirstandSecondDerivative,
ExternalFunctionNumericalFirstDerivative,
ExternalFunctionFirstDerivative,
ExternalFunctionNumericalSecondDerivative,
ExternalFunctionSecondDerivative
};
struct Block_contain_type
{
int Equation, Variable, Own_Derivative;
......@@ -362,6 +383,24 @@ public:
};
};
class FPUSH_ : public TagWithoutArgument
{
public:
inline FPUSH_() : TagWithoutArgument(FPUSH)
{
};
};
class FPOP_ : public TagWithoutArgument
{
public:
inline FPOP_() : TagWithoutArgument(FPOP)
{
};
};
class FDIMT_ : public TagWithOneArgument<unsigned int>
{
public:
......@@ -734,6 +773,136 @@ public:
}
};
class FLDTEF_ : public TagWithOneArgument<unsigned int>
{
public:
inline FLDTEF_() : TagWithOneArgument<unsigned int>::TagWithOneArgument(FLDTEF)
{
};
inline FLDTEF_(unsigned int number) : TagWithOneArgument<unsigned int>::TagWithOneArgument(FLDTEF, number)
{
};
inline unsigned int
get_number()
{
return arg1;
}
};
class FSTPTEF_ : public TagWithOneArgument<unsigned int>
{
public:
inline FSTPTEF_() : TagWithOneArgument<unsigned int>::TagWithOneArgument(FSTPTEF)
{
};
inline FSTPTEF_(unsigned int number) : TagWithOneArgument<unsigned int>::TagWithOneArgument(FSTPTEF, number)
{
};
inline unsigned int
get_number()
{
return arg1;
}
};
class FLDTEFD_ : public TagWithTwoArguments<unsigned int, unsigned int>
{
public:
inline FLDTEFD_() : TagWithTwoArguments<unsigned int, unsigned int>::TagWithTwoArguments(FLDTEFD)
{
};
inline FLDTEFD_(unsigned int indx, unsigned int row) : TagWithTwoArguments<unsigned int, unsigned int>::TagWithTwoArguments(FLDTEFD, indx, row)
{
};
inline unsigned int
get_indx()
{
return arg1;
};
inline unsigned int
get_row()
{
return arg2;
};
};
class FSTPTEFD_ : public TagWithTwoArguments<unsigned int, unsigned int>
{
public:
inline FSTPTEFD_() : TagWithTwoArguments<unsigned int, unsigned int>::TagWithTwoArguments(FSTPTEFD)
{
};
inline FSTPTEFD_(unsigned int indx, unsigned int row) : TagWithTwoArguments<unsigned int, unsigned int>::TagWithTwoArguments(FSTPTEFD, indx, row)
{
};
inline unsigned int
get_indx()
{
return arg1;
};
inline unsigned int
get_row()
{
return arg2;
};
};
class FLDTEFDD_ : public TagWithThreeArguments<unsigned int, unsigned int, unsigned int>
{
public:
inline FLDTEFDD_() : TagWithThreeArguments<unsigned int, unsigned int, unsigned int>::TagWithThreeArguments(FLDTEFDD)
{
};
inline FLDTEFDD_(unsigned int indx, unsigned int row, unsigned int col) : TagWithThreeArguments<unsigned int, unsigned int, unsigned int>::TagWithThreeArguments(FLDTEFDD, indx, row, col)
{
};
inline unsigned int
get_indx()
{
return arg1;
};
inline unsigned int
get_row()
{
return arg2;
};
inline unsigned int
get_col()
{
return arg3;
};
};
class FSTPTEFDD_ : public TagWithThreeArguments<unsigned int, unsigned int, unsigned int>
{
public:
inline FSTPTEFDD_() : TagWithThreeArguments<unsigned int, unsigned int, unsigned int>::TagWithThreeArguments(FSTPTEFDD)
{
};
inline FSTPTEFDD_(unsigned int indx, unsigned int row, unsigned int col) : TagWithThreeArguments<unsigned int, unsigned int, unsigned int>::TagWithThreeArguments(FSTPTEF, indx, row, col)
{
};
inline unsigned int
get_indx()
{
return arg1;
};
inline unsigned int
get_row()
{
return arg2;
};
inline unsigned int
get_col()
{
return arg3;
};
};
class FLDVS_ : public TagWithTwoArguments<uint8_t, unsigned int>
{
public:
......@@ -861,6 +1030,156 @@ public:
};
};
class FCALL_ : public TagWithFourArguments<unsigned int, unsigned int, string, unsigned int>
{
string func_name;
string arg_func_name;
unsigned int add_input_arguments, row, col;
external_function_type function_type;
public:
inline FCALL_() : TagWithFourArguments<unsigned int, unsigned int, string, unsigned int>::TagWithFourArguments(FCALL)
{
arg_func_name = "";
add_input_arguments = 0;
row = 0;
col = 0;
function_type = ExternalFunctionWithoutDerivative;
};
inline FCALL_(unsigned int nb_output_arguments, unsigned int nb_input_arguments, string f_name, unsigned int indx) :
TagWithFourArguments<unsigned int, unsigned int, string, unsigned int>::TagWithFourArguments(FCALL, nb_output_arguments, nb_input_arguments, f_name, indx)
{
arg_func_name = "";
add_input_arguments = 0;
row = 0;
col = 0;
function_type = ExternalFunctionWithoutDerivative;
func_name = f_name;
};
inline string
get_function_name()
{
//printf("get_function_name => func_name=%s\n",func_name.c_str());fflush(stdout);
return func_name;
};
inline unsigned int
get_nb_output_arguments()
{
return arg1;