diff --git a/mex/sources/bytecode/Interpreter.cc b/mex/sources/bytecode/Interpreter.cc
index 2c53681eca43762498b67c0e3ce29b3a0cdd9f5e..122f161fe5149697eb43d24faf6239fa868dbff0 100644
--- a/mex/sources/bytecode/Interpreter.cc
+++ b/mex/sources/bytecode/Interpreter.cc
@@ -242,8 +242,8 @@ Interpreter::print_expression(it_code_type it_code, bool evaluate)
stack<string> Stack;
stack<double> Stackf;
ostringstream tmp_out, tmp_out2;
- string v1, v2;
- double v1f, v2f;
+ string v1, v2, v3;
+ double v1f, v2f, v3f;
bool go_on = true;
double ll;
ExpressionType equation_type;
@@ -1024,10 +1024,46 @@ Interpreter::print_expression(it_code_type it_code, bool evaluate)
tmp_out << "sqrt(" << v1 << ")";
Stack.push(tmp_out.str());
break;
+ case oErf:
+ Stackf.push(sqrt(v1f));
+ tmp_out.str("");
+ tmp_out << "erf(" << v1 << ")";
+ Stack.push(tmp_out.str());
+ break;
default:
;
}
break;
+ case FTRINARY:
+ op = ((FTRINARY_ *) it_code->second)->get_op_type();
+ v3 = Stack.top();
+ Stack.pop();
+ v2 = Stack.top();
+ Stack.pop();
+ v1 = Stack.top();
+ Stack.pop();
+ v3f = Stackf.top();
+ Stackf.pop();
+ v2f = Stackf.top();
+ Stackf.pop();
+ v1f = Stackf.top();
+ Stackf.pop();
+ switch (op)
+ {
+ case oNormcdf:
+ Stackf.push(0.5*(1+erf((v1f-v2f)/v3f/M_SQRT2)));
+ tmp_out.str("");
+ tmp_out << "normcdf(" << v1 << ", " << v2 << ", " << v3 << ")";
+ Stack.push(tmp_out.str());
+ break;
+ case oNormpdf:
+ Stackf.push(1/(v3f*sqrt(2*M_PI)*exp(pow((v1f-v2f)/v3f,2)/2)));
+ tmp_out.str("");
+ tmp_out << "normpdf(" << v1 << ", " << v2 << ", " << v3 << ")";
+ Stack.push(tmp_out.str());
+ break;
+ }
+ break;
case FCUML:
case FENDBLOCK:
case FOK:
@@ -1049,7 +1085,7 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate, int block_num)
{
int var, lag = 0, op;
ostringstream tmp_out;
- double v1, v2;
+ double v1, v2, v3;
bool go_on = true;
double ll;
@@ -1727,12 +1763,44 @@ Interpreter::compute_block_time(int Per_u_, bool evaluate, int block_num)
Stack.push(sqrt(v1));
#ifdef DEBUG
tmp_out << " |sqrt(" << v1 << ")|";
+#endif
+ break;
+ case oErf:
+ Stack.push(erf(v1));
+#ifdef DEBUG
+ tmp_out << " |erf(" << v1 << ")|";
#endif
break;
default:
;
}
break;
+ case FTRINARY:
+ op = ((FTRINARY_ *) it_code->second)->get_op_type();
+ v3 = Stack.top();
+ Stack.pop();
+ v2 = Stack.top();
+ Stack.pop();
+ v1 = Stack.top();
+ Stack.pop();
+ switch (op)
+ {
+ case oNormcdf:
+ //mexPrintf("normcdf(v1=%f, v2=%f, v3=%f)=%f\n", v1, v2, v3, 0.5*(1+erf((v1-v2)/v3/M_SQRT2)));
+ Stack.push(0.5*(1+erf((v1-v2)/v3/M_SQRT2)));
+#ifdef DEBUG
+ tmp_out << " |normcdf(" << v1 << ", " << v2 << ", " << v3 << ")|";
+#endif
+ break;
+ case oNormpdf:
+ Stack.push(1/(v3*sqrt(2*M_PI)*exp(pow((v1-v2)/v3,2)/2)));
+#ifdef DEBUG
+ tmp_out << " |normpdf(" << v1 << ", " << v2 << ", " << v3 << ")|";
+#endif
+ break;
+ }
+ break;
+
case FCUML:
v1 = Stack.top();
Stack.pop();
diff --git a/preprocessor/CodeInterpreter.hh b/preprocessor/CodeInterpreter.hh
index b542b4b52168379150ae27905c16642e1f536a0a..3d8952cd940d3676ae89ac02bf426e8f0d49eeb3 100644
--- a/preprocessor/CodeInterpreter.hh
+++ b/preprocessor/CodeInterpreter.hh
@@ -86,6 +86,7 @@ enum Tags
FUNARY, //!< A Unary operator - 14
FBINARY, //!< A binary operator - 15
+ FTRINARY, //!< A trinary operator - 15'
FCUML, //!< Cumulates the result - 16
@@ -610,6 +611,22 @@ public:
};
};
+class FTRINARY_ : public TagWithOneArgument<uint8_t>
+{
+public:
+ inline FTRINARY_() : TagWithOneArgument<uint8_t>::TagWithOneArgument(FTRINARY)
+ {
+ };
+ inline FTRINARY_(const int op_type_arg) : TagWithOneArgument<uint8_t>::TagWithOneArgument(FTRINARY, op_type_arg)
+ {
+ };
+ inline uint8_t
+ get_op_type()
+ {
+ return arg1;
+ };
+};
+
class FOK_ : public TagWithOneArgument<int>
{
public:
@@ -1178,6 +1195,13 @@ public:
tags_liste.push_back(make_pair(FBINARY, code));
code += sizeof(FBINARY_);
break;
+ case FTRINARY:
+# ifdef DEBUGL
+ mexPrintf("FTRINARY\n");
+# endif
+ tags_liste.push_back(make_pair(FTRINARY, code));
+ code += sizeof(FTRINARY_);
+ break;
case FOK:
# ifdef DEBUGL
mexPrintf("FOK\n");
diff --git a/preprocessor/ExprNode.cc b/preprocessor/ExprNode.cc
index 687f1b8cc726ea06c8e87ca6bde89338068e1a61..ea87a50ebf2e96b7b5be370b3562484d1d5bae03 100644
--- a/preprocessor/ExprNode.cc
+++ b/preprocessor/ExprNode.cc
@@ -237,6 +237,18 @@ ExprNode::createExoLeadAuxiliaryVarForMyself(subst_table_t &subst_table, vector<
return dynamic_cast<VariableNode *>(substexpr);
}
+bool
+ExprNode::isNumConstNodeEqualTo(double value) const
+{
+ return false;
+}
+
+bool
+ExprNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ return false;
+}
+
NumConstNode::NumConstNode(DataTree &datatree_arg, int id_arg) :
ExprNode(datatree_arg),
id(id_arg)
@@ -385,6 +397,21 @@ VariableNode::VariableNode(DataTree &datatree_arg, int symb_id_arg, int lag_arg)
&& (lag == 0 || (type != eModelLocalVariable && type != eModFileLocalVariable)));
}
+bool
+NumConstNode::isNumConstNodeEqualTo(double value) const
+{
+ if (datatree.num_constants.getDouble(id) == value)
+ return true;
+ else
+ return false;
+}
+
+bool
+NumConstNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ return false;
+}
+
void
VariableNode::prepareForDerivation()
{
@@ -999,6 +1026,21 @@ VariableNode::substituteExpectation(subst_table_t &subst_table, vector<BinaryOpN
return const_cast<VariableNode *>(this);
}
+bool
+VariableNode::isNumConstNodeEqualTo(double value) const
+{
+ return false;
+}
+
+bool
+VariableNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ if (type == type_arg && datatree.symbol_table.getTypeSpecificID(symb_id) == variable_id && lag == lag_arg)
+ return true;
+ else
+ return false;
+}
+
UnaryOpNode::UnaryOpNode(DataTree &datatree_arg, UnaryOpcode op_code_arg, const NodeID arg_arg, const int expectation_information_set_arg, const string &expectation_information_set_name_arg) :
ExprNode(datatree_arg),
arg(arg_arg),
@@ -1820,6 +1862,18 @@ UnaryOpNode::substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNo
}
}
+bool
+UnaryOpNode::isNumConstNodeEqualTo(double value) const
+{
+ return false;
+}
+
+bool
+UnaryOpNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ return false;
+}
+
BinaryOpNode::BinaryOpNode(DataTree &datatree_arg, const NodeID arg1_arg,
BinaryOpcode op_code_arg, const NodeID arg2_arg) :
ExprNode(datatree_arg),
@@ -2810,6 +2864,18 @@ BinaryOpNode::substituteExpectation(subst_table_t &subst_table, vector<BinaryOpN
return buildSimilarBinaryOpNode(arg1subst, arg2subst, datatree);
}
+bool
+BinaryOpNode::isNumConstNodeEqualTo(double value) const
+{
+ return false;
+}
+
+bool
+BinaryOpNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ return false;
+}
+
TrinaryOpNode::TrinaryOpNode(DataTree &datatree_arg, const NodeID arg1_arg,
TrinaryOpcode op_code_arg, const NodeID arg2_arg, const NodeID arg3_arg) :
ExprNode(datatree_arg),
@@ -3077,8 +3143,8 @@ TrinaryOpNode::compile(ostream &CompileCode, bool lhs_rhs, const temporary_terms
arg1->compile(CompileCode, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic);
arg2->compile(CompileCode, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic);
arg3->compile(CompileCode, lhs_rhs, temporary_terms, map_idx, dynamic, steady_dynamic);
- FBINARY_ fbinary(op_code);
- fbinary.write(CompileCode);
+ FTRINARY_ ftrinary(op_code);
+ ftrinary.write(CompileCode);
}
void
@@ -3306,6 +3372,18 @@ TrinaryOpNode::substituteExpectation(subst_table_t &subst_table, vector<BinaryOp
return buildSimilarTrinaryOpNode(arg1subst, arg2subst, arg3subst, datatree);
}
+bool
+TrinaryOpNode::isNumConstNodeEqualTo(double value) const
+{
+ return false;
+}
+
+bool
+TrinaryOpNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ return false;
+}
+
ExternalFunctionNode::ExternalFunctionNode(DataTree &datatree_arg,
int symb_id_arg,
const vector<NodeID> &arguments_arg) :
@@ -3632,6 +3710,18 @@ ExternalFunctionNode::getIndxInTefTerms(int the_symb_id, deriv_node_temp_terms_t
throw UnknownFunctionNameAndArgs();
}
+bool
+ExternalFunctionNode::isNumConstNodeEqualTo(double value) const
+{
+ return false;
+}
+
+bool
+ExternalFunctionNode::isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const
+{
+ return false;
+}
+
FirstDerivExternalFunctionNode::FirstDerivExternalFunctionNode(DataTree &datatree_arg,
int top_level_symb_id_arg,
const vector<NodeID> &arguments_arg,
diff --git a/preprocessor/ExprNode.hh b/preprocessor/ExprNode.hh
index bfe3f6c31a72308e06a27a13983b8ff6b14ad2ef..52a6b04bc5ad9c4063488c0143c33c13489c40c8 100644
--- a/preprocessor/ExprNode.hh
+++ b/preprocessor/ExprNode.hh
@@ -333,6 +333,20 @@ public:
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const = 0;
+ //! Return true if the nodeID is a numerical constant equal to value and false otherwise
+ /*!
+ \param[in] value of the numerical constante
+ \param[out] the boolean equal to true if NodeId is a constant equal to value
+ */
+ virtual bool isNumConstNodeEqualTo(double value) const = 0;
+
+ //! Return true if the nodeID is a variable withe a type equal to type_arg, a specific variable id aqual to varfiable_id and a lag equal to lag_arg and false otherwise
+ /*!
+ \param[in] the type (type_arg), specifique variable id (variable_id and the lag (lag_arg)
+ \param[out] the boolean equal to true if NodeId is the variable
+ */
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const = 0;
+
};
//! Object used to compare two nodes (using their indexes)
@@ -378,6 +392,8 @@ public:
virtual NodeID substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual NodeID substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const;
+ virtual bool isNumConstNodeEqualTo(double value) const;
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const;
};
//! Symbol or variable node
@@ -421,6 +437,8 @@ public:
virtual NodeID substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual NodeID substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const;
+ virtual bool isNumConstNodeEqualTo(double value) const;
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const;
};
//! Unary operator node
@@ -482,6 +500,8 @@ public:
virtual NodeID substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual NodeID substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const;
+ virtual bool isNumConstNodeEqualTo(double value) const;
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const;
};
//! Binary operator node
@@ -548,6 +568,8 @@ public:
virtual NodeID substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual NodeID substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const;
+ virtual bool isNumConstNodeEqualTo(double value) const;
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const;
};
//! Trinary operator node
@@ -596,6 +618,8 @@ public:
virtual NodeID substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual NodeID substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const;
+ virtual bool isNumConstNodeEqualTo(double value) const;
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const;
};
//! External function node
@@ -649,6 +673,8 @@ public:
virtual NodeID substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual NodeID buildSimilarExternalFunctionNode(vector<NodeID> &alt_args, DataTree &alt_datatree) const;
virtual NodeID decreaseLeadsLagsPredeterminedVariables() const;
+ virtual bool isNumConstNodeEqualTo(double value) const;
+ virtual bool isVariableNodeEqualTo(SymbolType type_arg, int variable_id, int lag_arg) const;
};
class FirstDerivExternalFunctionNode : public ExternalFunctionNode
diff --git a/preprocessor/ModelTree.cc b/preprocessor/ModelTree.cc
index 2bcce6150f6d7b632f24300c6f016a455a71d320..afd4f8c1c3b4cd21095f5f62462b485e1c0a6380 100644
--- a/preprocessor/ModelTree.cc
+++ b/preprocessor/ModelTree.cc
@@ -414,10 +414,7 @@ t_equation_type_and_normalized_equation
ModelTree::equationTypeDetermination(vector<BinaryOpNode *> &equations, map<pair<int, pair<int, int> >, NodeID> &first_order_endo_derivatives, vector<int> &Index_Var_IM, vector<int> &Index_Equ_IM, int mfs)
{
NodeID lhs, rhs;
- ostringstream tmp_output;
BinaryOpNode *eq_node;
- ostringstream tmp_s;
- temporary_terms_type temporary_terms;
EquationType Equation_Simulation_Type;
t_equation_type_and_normalized_equation V_Equation_Simulation_Type(equations.size());
for (unsigned int i = 0; i < equations.size(); i++)
@@ -429,10 +426,6 @@ ModelTree::equationTypeDetermination(vector<BinaryOpNode *> &equations, map<pair
lhs = eq_node->get_arg1();
rhs = eq_node->get_arg2();
Equation_Simulation_Type = E_SOLVE;
- tmp_s.str("");
- tmp_output.str("");
- lhs->writeOutput(tmp_output, oMatlabDynamicModelSparse, temporary_terms);
- 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)));
pair<bool, NodeID> res;
if (derivative != first_order_endo_derivatives.end())
@@ -441,9 +434,7 @@ ModelTree::equationTypeDetermination(vector<BinaryOpNode *> &equations, map<pair
derivative->second->collectEndogenous(result);
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);
- if (tmp_output.str() == tmp_s.str() and tt.str() == "1")
+ if (lhs->isVariableNodeEqualTo(eEndogenous, Index_Var_IM[i], 0) && derivative->second->isNumConstNodeEqualTo(1))
{
Equation_Simulation_Type = E_EVALUATE;
}