diff --git a/src/DynamicModel.cc b/src/DynamicModel.cc
index 84aefaafcfd0bf90aeb00de059978e3935c44bd8..3eb18aa3d17dba318775dc12ac2ae66e4875f817 100644
--- a/src/DynamicModel.cc
+++ b/src/DynamicModel.cc
@@ -3972,7 +3972,58 @@ DynamicModel::writeDynamicFile(const string &basename, bool block, bool bytecode
else if (julia)
writeDynamicJuliaFile(basename);
else
- writeDynamicMFile(t_basename);
+ {
+ writeDynamicMFile(t_basename);
+ writeSetAuxiliaryVariables(t_basename, julia);
+ }
+}
+
+void
+DynamicModel::writeSetAuxiliaryVariables(const string &basename, const bool julia) const
+{
+ ostringstream output_func_body;
+ writeAuxVarRecursiveDefinitions(output_func_body, oMatlabDseries);
+
+ if (output_func_body.str().empty())
+ return;
+
+ string func_name = basename + "_set_auxiliary_series";
+ string filename = julia ? func_name + ".jl" : func_name + ".m";
+ string comment = julia ? "#" : "%";
+
+ ofstream output;
+ output.open(filename.c_str(), ios::out | ios::binary);
+ if (!output.is_open())
+ {
+ cerr << "ERROR: Can't open file " << filename << " for writing" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ output << "function ds = " << func_name + "(ds, params)" << endl
+ << comment << endl
+ << comment << " Status : Computes Auxiliary variables of the dynamic model and returns a dseries" << endl
+ << comment << endl
+ << comment << " Warning : this file is generated automatically by Dynare" << endl
+ << comment << " from model file (.mod)" << endl << endl
+ << output_func_body.str();
+
+ output.close();
+}
+
+void
+DynamicModel::writeAuxVarRecursiveDefinitions(ostream &output, ExprNodeOutputType output_type) const
+{
+ deriv_node_temp_terms_t tef_terms;
+ temporary_terms_t temporary_terms;
+ for (int i = 0; i < (int) aux_equations.size(); i++)
+ if (dynamic_cast<ExprNode *>(aux_equations[i])->containsExternalFunction())
+ dynamic_cast<ExprNode *>(aux_equations[i])->writeExternalFunctionOutput(output, output_type,
+ temporary_terms, tef_terms);
+ for (int i = 0; i < (int) aux_equations.size(); i++)
+ {
+ dynamic_cast<ExprNode *>(aux_equations[i])->writeOutput(output, output_type, temporary_terms, tef_terms);
+ output << ";" << endl;
+ }
}
void
@@ -4848,6 +4899,39 @@ DynamicModel::substituteLeadLagInternal(aux_var_t type, bool deterministic_model
aux_equations[i] = substeq;
}
+ // Substitute in diff_aux_equations
+ // Without this loop, the auxiliary equations in equations
+ // will diverge from those in diff_aux_equations
+ for (int i = 0; i < (int) diff_aux_equations.size(); i++)
+ {
+ expr_t subst;
+ switch (type)
+ {
+ case avEndoLead:
+ subst = diff_aux_equations[i]->substituteEndoLeadGreaterThanTwo(subst_table,
+ neweqs, deterministic_model);
+ break;
+ case avEndoLag:
+ subst = diff_aux_equations[i]->substituteEndoLagGreaterThanTwo(subst_table, neweqs);
+ break;
+ case avExoLead:
+ subst = diff_aux_equations[i]->substituteExoLead(subst_table, neweqs, deterministic_model);
+ break;
+ case avExoLag:
+ subst = diff_aux_equations[i]->substituteExoLag(subst_table, neweqs);
+ break;
+ case avDiffForward:
+ subst = diff_aux_equations[i]->differentiateForwardVars(subset, subst_table, neweqs);
+ break;
+ default:
+ cerr << "DynamicModel::substituteLeadLagInternal: impossible case" << endl;
+ exit(EXIT_FAILURE);
+ }
+ BinaryOpNode *substeq = dynamic_cast<BinaryOpNode *>(subst);
+ assert(substeq != NULL);
+ diff_aux_equations[i] = substeq;
+ }
+
// Add new equations
for (int i = 0; i < (int) neweqs.size(); i++)
addEquation(neweqs[i], -1);
@@ -4897,20 +4981,29 @@ DynamicModel::substituteAdl()
}
void
-DynamicModel::substituteDiff()
+DynamicModel::substituteDiff(StaticModel &static_model)
{
- ExprNode::subst_table_t subst_table;
- vector<BinaryOpNode *> neweqs;
+ // Find diff Nodes
+ diff_table_t diff_table;
+ for (map<int, expr_t>::iterator it = local_variables_table.begin();
+ it != local_variables_table.end(); it++)
+ it->second->findDiffNodes(static_model, diff_table);
+
+ for (int i = 0; i < (int) equations.size(); i++)
+ equations[i]->findDiffNodes(static_model, diff_table);
// Substitute in model local variables
+ vector<BinaryOpNode *> neweqs;
+ ExprNode::subst_table_t diff_subst_table;
for (map<int, expr_t>::iterator it = local_variables_table.begin();
it != local_variables_table.end(); it++)
- it->second = it->second->substituteDiff(subst_table, neweqs);
+ it->second = it->second->substituteDiff(static_model, diff_table, diff_subst_table, neweqs);
// Substitute in equations
for (int i = 0; i < (int) equations.size(); i++)
{
- BinaryOpNode *substeq = dynamic_cast<BinaryOpNode *>(equations[i]->substituteDiff(subst_table, neweqs));
+ BinaryOpNode *substeq = dynamic_cast<BinaryOpNode *>(equations[i]->
+ substituteDiff(static_model, diff_table, diff_subst_table, neweqs));
assert(substeq != NULL);
equations[i] = substeq;
}
@@ -4919,10 +5012,18 @@ DynamicModel::substituteDiff()
for (int i = 0; i < (int) neweqs.size(); i++)
addEquation(neweqs[i], -1);
- if (subst_table.size() > 0)
+ copy(neweqs.begin(), neweqs.end(), back_inserter(diff_aux_equations));
+
+ if (diff_subst_table.size() > 0)
cout << "Substitution of Diff operator: added " << neweqs.size() << " auxiliary variables and equations." << endl;
}
+void
+DynamicModel::combineDiffAuxEquations()
+{
+ copy(diff_aux_equations.begin(), diff_aux_equations.end(), back_inserter(aux_equations));
+}
+
void
DynamicModel::substituteExpectation(bool partial_information_model)
{
diff --git a/src/DynamicModel.hh b/src/DynamicModel.hh
index afcd72b9bb34ecc2aadb03360f2445ac85dba70e..ad41e03f5f9c7b2dead60ab4a4bb761296cf3a91 100644
--- a/src/DynamicModel.hh
+++ b/src/DynamicModel.hh
@@ -116,6 +116,9 @@ private:
//! Writes the code of the model in virtual machine bytecode
void writeModelEquationsCode(string &file_name, const string &bin_basename, const map_idx_t &map_idx) const;
+ void writeSetAuxiliaryVariables(const string &basename, const bool julia) const;
+ void writeAuxVarRecursiveDefinitions(ostream &output, ExprNodeOutputType output_type) const;
+
//! Computes jacobian and prepares for equation normalization
/*! Using values from initval/endval blocks and parameter initializations:
- computes the jacobian for the model w.r. to contemporaneous variables
@@ -392,7 +395,10 @@ public:
void substituteAdl();
//! Substitutes diff operator
- void substituteDiff();
+ void substituteDiff(StaticModel &static_model);
+
+ //! Adds contents of diff_aux_equations to the back of aux_equations
+ void combineDiffAuxEquations();
//! Fill var_expectation_functions_to_write
void fillVarExpectationFunctionsToWrite();
diff --git a/src/ExprNode.cc b/src/ExprNode.cc
index 43095bcfcffc0a35acca8ef61a716a4d808d62f3..f6393c7d49a6ab442f4f40abde4cb67ecd378af6 100644
--- a/src/ExprNode.cc
+++ b/src/ExprNode.cc
@@ -451,6 +451,12 @@ NumConstNode::maxLead() const
return 0;
}
+int
+NumConstNode::maxLag() const
+{
+ return 0;
+}
+
expr_t
NumConstNode::decreaseLeadsLags(int n) const
{
@@ -499,8 +505,13 @@ NumConstNode::substituteAdl() const
return const_cast<NumConstNode *>(this);
}
+void
+NumConstNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+}
+
expr_t
-NumConstNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+NumConstNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
return const_cast<NumConstNode *>(this);
}
@@ -839,6 +850,11 @@ VariableNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
case oCSteadyStateFile:
output << "ys_[" << tsid << "]";
break;
+ case oMatlabDseries:
+ output << "ds." << datatree.symbol_table.getName(symb_id);
+ if (lag != 0)
+ output << LEFT_ARRAY_SUBSCRIPT(output_type) << lag << RIGHT_ARRAY_SUBSCRIPT(output_type);
+ break;
default:
cerr << "VariableNode::writeOutput: should not reach this point" << endl;
exit(EXIT_FAILURE);
@@ -891,6 +907,11 @@ VariableNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
case oCSteadyStateFile:
output << "exo_[" << i - 1 << "]";
break;
+ case oMatlabDseries:
+ output << "ds." << datatree.symbol_table.getName(symb_id);
+ if (lag != 0)
+ output << LEFT_ARRAY_SUBSCRIPT(output_type) << lag << RIGHT_ARRAY_SUBSCRIPT(output_type);
+ break;
default:
cerr << "VariableNode::writeOutput: should not reach this point" << endl;
exit(EXIT_FAILURE);
@@ -943,6 +964,11 @@ VariableNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
case oCSteadyStateFile:
output << "exo_[" << i - 1 << "]";
break;
+ case oMatlabDseries:
+ output << "ds." << datatree.symbol_table.getName(symb_id);
+ if (lag != 0)
+ output << LEFT_ARRAY_SUBSCRIPT(output_type) << lag << RIGHT_ARRAY_SUBSCRIPT(output_type);
+ break;
default:
cerr << "VariableNode::writeOutput: should not reach this point" << endl;
exit(EXIT_FAILURE);
@@ -1282,14 +1308,35 @@ VariableNode::maxLead() const
}
}
+int
+VariableNode::maxLag() const
+{
+ switch (type)
+ {
+ case eEndogenous:
+ return -lag;
+ case eExogenous:
+ return -lag;
+ case eModelLocalVariable:
+ return datatree.local_variables_table[symb_id]->maxLag();
+ default:
+ return 0;
+ }
+}
+
expr_t
VariableNode::substituteAdl() const
{
return const_cast<VariableNode *>(this);
}
+void
+VariableNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+}
+
expr_t
-VariableNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+VariableNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
return const_cast<VariableNode *>(this);
}
@@ -2788,6 +2835,12 @@ UnaryOpNode::maxLead() const
return arg->maxLead();
}
+int
+UnaryOpNode::maxLag() const
+{
+ return arg->maxLag();
+}
+
expr_t
UnaryOpNode::substituteAdl() const
{
@@ -2829,38 +2882,101 @@ UnaryOpNode::isDiffPresent() const
return arg->isDiffPresent();
}
+void
+UnaryOpNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+ if (op_code != oDiff)
+ return;
+
+ arg->findDiffNodes(static_datatree, diff_table);
+
+ expr_t sthis = this->toStatic(static_datatree);
+ int arg_max_lag = -arg->maxLag();
+ diff_table_t::iterator it = diff_table.find(sthis);
+ if (it != diff_table.end())
+ {
+ for (map<int, expr_t>::const_iterator it1 = it->second.begin();
+ it1 != it->second.end(); it1++)
+ if (arg == it1->second)
+ return;
+ it->second[arg_max_lag] = const_cast<UnaryOpNode *>(this);
+ }
+ else
+ diff_table[sthis][arg_max_lag] = const_cast<UnaryOpNode *>(this);
+}
+
expr_t
-UnaryOpNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+UnaryOpNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
if (op_code != oDiff)
{
- expr_t argsubst = arg->substituteDiff(subst_table, neweqs);
+ expr_t argsubst = arg->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
return buildSimilarUnaryOpNode(argsubst, datatree);
}
- subst_table_t::iterator it = subst_table.find(const_cast<UnaryOpNode *>(this));
- if (it != subst_table.end())
- return const_cast<VariableNode *>(it->second);
+ subst_table_t::const_iterator sit = subst_table.find(this);
+ if (sit != subst_table.end())
+ return const_cast<VariableNode *>(sit->second);
- expr_t argsubst = arg->substituteDiff(subst_table, neweqs);
- assert(argsubst != NULL);
+ expr_t sthis = dynamic_cast<UnaryOpNode *>(this->toStatic(static_datatree));
+ diff_table_t::iterator it = diff_table.find(sthis);
+ assert(it != diff_table.end() && it->second[-arg->maxLag()] == this);
- int symb_id = -1;
- VariableNode *vn = dynamic_cast<VariableNode *>(argsubst);
- if (vn != NULL)
- symb_id = datatree.symbol_table.addDiffAuxiliaryVar(argsubst->idx, argsubst, vn->get_symb_id(), vn->get_lag());
- else
- symb_id = datatree.symbol_table.addDiffAuxiliaryVar(argsubst->idx, argsubst);
-
- expr_t newAuxVar = datatree.AddVariable(symb_id, 0);
-
- //AUX_DIFF_IDX = argsubst - argsubst(-1)
- neweqs.push_back(dynamic_cast<BinaryOpNode *>(datatree.AddEqual(newAuxVar,
- datatree.AddMinus(argsubst, argsubst->decreaseLeadsLags(1)))));
-
- assert(dynamic_cast<VariableNode *>(newAuxVar) != NULL);
- subst_table[this] = dynamic_cast<VariableNode *>(newAuxVar);
- return newAuxVar;
+ int last_arg_max_lag;
+ VariableNode *last_aux_var;
+ for (map<int, expr_t>::reverse_iterator rit = it->second.rbegin();
+ rit != it->second.rend(); rit++)
+ {
+ expr_t argsubst = dynamic_cast<UnaryOpNode *>(rit->second)->
+ get_arg()->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
+ int symb_id;
+ VariableNode *vn = dynamic_cast<VariableNode *>(argsubst);
+ if (rit == it->second.rbegin())
+ {
+ if (vn != NULL)
+ symb_id = datatree.symbol_table.addDiffAuxiliaryVar(argsubst->idx, argsubst, vn->get_symb_id(), vn->get_lag());
+ else
+ symb_id = datatree.symbol_table.addDiffAuxiliaryVar(argsubst->idx, argsubst);
+
+ // make originating aux var & equation
+ last_arg_max_lag = rit->first;
+ last_aux_var = datatree.AddVariable(symb_id, 0);
+ //ORIG_AUX_DIFF = argsubst - argsubst(-1)
+ neweqs.push_back(dynamic_cast<BinaryOpNode *>(datatree.AddEqual(last_aux_var,
+ datatree.AddMinus(argsubst,
+ argsubst->decreaseLeadsLags(1)))));
+ subst_table[rit->second] = dynamic_cast<VariableNode *>(last_aux_var);
+ }
+ else
+ {
+ // just add equation of form: AUX_DIFF = ORIG_AUX_DIFF(last_arg_max_lag - rit->first)
+ VariableNode *new_aux_var;
+ for (int i = last_arg_max_lag; i > rit->first; i--)
+ {
+ if (vn == NULL)
+ if (i == last_arg_max_lag)
+ symb_id = datatree.symbol_table.addDiffAuxiliaryVar(argsubst->idx, argsubst);
+ else
+ symb_id = datatree.symbol_table.addDiffAuxiliaryVar(new_aux_var->idx, new_aux_var);
+ else
+ if (i == last_arg_max_lag)
+ symb_id = datatree.symbol_table.addDiffAuxiliaryVar(argsubst->idx, argsubst,
+ vn->get_symb_id(), i - 1);
+ else
+ symb_id = datatree.symbol_table.addDiffAuxiliaryVar(new_aux_var->idx, new_aux_var,
+ vn->get_symb_id(), i - 1);
+
+
+ new_aux_var = datatree.AddVariable(symb_id, 0);
+ neweqs.push_back(dynamic_cast<BinaryOpNode *>(datatree.AddEqual(new_aux_var,
+ last_aux_var->decreaseLeadsLags(1))));
+ last_aux_var = new_aux_var;
+ }
+ subst_table[rit->second] = dynamic_cast<VariableNode *>(new_aux_var);
+ last_arg_max_lag = rit->first;
+ }
+ }
+ return const_cast<VariableNode *>(subst_table[this]);
}
expr_t
@@ -4355,6 +4471,12 @@ BinaryOpNode::maxLead() const
return max(arg1->maxLead(), arg2->maxLead());
}
+int
+BinaryOpNode::maxLag() const
+{
+ return max(arg1->maxLag(), arg2->maxLag());
+}
+
expr_t
BinaryOpNode::decreaseLeadsLags(int n) const
{
@@ -4491,11 +4613,18 @@ BinaryOpNode::substituteAdl() const
return buildSimilarBinaryOpNode(arg1subst, arg2subst, datatree);
}
+void
+BinaryOpNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+ arg1->findDiffNodes(static_datatree, diff_table);
+ arg2->findDiffNodes(static_datatree, diff_table);
+}
+
expr_t
-BinaryOpNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+BinaryOpNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
- expr_t arg1subst = arg1->substituteDiff(subst_table, neweqs);
- expr_t arg2subst = arg2->substituteDiff(subst_table, neweqs);
+ expr_t arg1subst = arg1->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
+ expr_t arg2subst = arg2->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
return buildSimilarBinaryOpNode(arg1subst, arg2subst, datatree);
}
@@ -5260,6 +5389,12 @@ TrinaryOpNode::maxLead() const
return max(arg1->maxLead(), max(arg2->maxLead(), arg3->maxLead()));
}
+int
+TrinaryOpNode::maxLag() const
+{
+ return max(arg1->maxLag(), max(arg2->maxLag(), arg3->maxLag()));
+}
+
expr_t
TrinaryOpNode::decreaseLeadsLags(int n) const
{
@@ -5346,12 +5481,20 @@ TrinaryOpNode::substituteAdl() const
return buildSimilarTrinaryOpNode(arg1subst, arg2subst, arg3subst, datatree);
}
+void
+TrinaryOpNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+ arg1->findDiffNodes(static_datatree, diff_table);
+ arg2->findDiffNodes(static_datatree, diff_table);
+ arg3->findDiffNodes(static_datatree, diff_table);
+}
+
expr_t
-TrinaryOpNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+TrinaryOpNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
- expr_t arg1subst = arg1->substituteDiff(subst_table, neweqs);
- expr_t arg2subst = arg2->substituteDiff(subst_table, neweqs);
- expr_t arg3subst = arg3->substituteDiff(subst_table, neweqs);
+ expr_t arg1subst = arg1->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
+ expr_t arg2subst = arg2->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
+ expr_t arg3subst = arg3->substituteDiff(static_datatree, diff_table, subst_table, neweqs);
return buildSimilarTrinaryOpNode(arg1subst, arg2subst, arg3subst, datatree);
}
@@ -5634,6 +5777,16 @@ AbstractExternalFunctionNode::maxLead() const
return val;
}
+int
+AbstractExternalFunctionNode::maxLag() const
+{
+ int val = 0;
+ for (vector<expr_t>::const_iterator it = arguments.begin();
+ it != arguments.end(); it++)
+ val = max(val, (*it)->maxLag());
+ return val;
+}
+
expr_t
AbstractExternalFunctionNode::decreaseLeadsLags(int n) const
{
@@ -5706,12 +5859,19 @@ AbstractExternalFunctionNode::substituteAdl() const
return buildSimilarExternalFunctionNode(arguments_subst, datatree);
}
+void
+AbstractExternalFunctionNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+ for (vector<expr_t>::const_iterator it = arguments.begin(); it != arguments.end(); it++)
+ (*it)->findDiffNodes(static_datatree, diff_table);
+}
+
expr_t
-AbstractExternalFunctionNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+AbstractExternalFunctionNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
vector<expr_t> arguments_subst;
for (vector<expr_t>::const_iterator it = arguments.begin(); it != arguments.end(); it++)
- arguments_subst.push_back((*it)->substituteDiff(subst_table, neweqs));
+ arguments_subst.push_back((*it)->substituteDiff(static_datatree, diff_table, subst_table, neweqs));
return buildSimilarExternalFunctionNode(arguments_subst, datatree);
}
@@ -7112,6 +7272,12 @@ VarExpectationNode::maxLead() const
return 0;
}
+int
+VarExpectationNode::maxLag() const
+{
+ return 0;
+}
+
expr_t
VarExpectationNode::decreaseLeadsLags(int n) const
{
@@ -7229,8 +7395,13 @@ VarExpectationNode::substituteAdl() const
return const_cast<VarExpectationNode *>(this);
}
+void
+VarExpectationNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+}
+
expr_t
-VarExpectationNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+VarExpectationNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
return const_cast<VarExpectationNode *>(this);
}
@@ -7530,6 +7701,12 @@ PacExpectationNode::maxLead() const
return 0;
}
+int
+PacExpectationNode::maxLag() const
+{
+ return 0;
+}
+
expr_t
PacExpectationNode::decreaseLeadsLags(int n) const
{
@@ -7646,8 +7823,13 @@ PacExpectationNode::substituteAdl() const
return const_cast<PacExpectationNode *>(this);
}
+void
+PacExpectationNode::findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const
+{
+}
+
expr_t
-PacExpectationNode::substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
+PacExpectationNode::substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
return const_cast<PacExpectationNode *>(this);
}
diff --git a/src/ExprNode.hh b/src/ExprNode.hh
index 7fe371cc4368230694808cdbe7d39929f8a88731..1ea6dbb5a48020dd1078e70b9776a724862ec312 100644
--- a/src/ExprNode.hh
+++ b/src/ExprNode.hh
@@ -34,6 +34,7 @@ using namespace std;
class DataTree;
class VariableNode;
+class UnaryOpNode;
class BinaryOpNode;
class PacExpectationNode;
@@ -57,6 +58,9 @@ typedef map<int, double> eval_context_t;
//! Type for tracking first/second derivative functions that have already been written as temporary terms
typedef map<pair<int, vector<expr_t> >, int> deriv_node_temp_terms_t;
+//! Type for the substitution map used in the process of substitutitng diff expressions
+typedef map<expr_t, map<int, expr_t> > diff_table_t;
+
//! Possible types of output when writing ExprNode(s)
enum ExprNodeOutputType
{
@@ -79,7 +83,8 @@ enum ExprNodeOutputType
oJuliaDynamicSteadyStateOperator, //!< Julia code, dynamic model, inside a steady state operator
oSteadyStateFile, //!< Matlab code, in the generated steady state file
oCSteadyStateFile, //!< C code, in the generated steady state file
- oJuliaSteadyStateFile //!< Julia code, in the generated steady state file
+ oJuliaSteadyStateFile, //!< Julia code, in the generated steady state file
+ oMatlabDseries //!< Matlab code for dseries
};
#define IS_MATLAB(output_type) ((output_type) == oMatlabStaticModel \
@@ -89,7 +94,8 @@ enum ExprNodeOutputType
|| (output_type) == oMatlabDynamicModelSparse \
|| (output_type) == oMatlabDynamicSteadyStateOperator \
|| (output_type) == oMatlabDynamicSparseSteadyStateOperator \
- || (output_type) == oSteadyStateFile)
+ || (output_type) == oSteadyStateFile \
+ || (output_type) == oMatlabDseries)
#define IS_JULIA(output_type) ((output_type) == oJuliaStaticModel \
|| (output_type) == oJuliaDynamicModel \
@@ -339,6 +345,10 @@ class ExprNode
/*! A negative value means that the expression contains only lagged variables */
virtual int maxLead() const = 0;
+ //! Returns the relative period of the most backward term in this expression
+ /*! A negative value means that the expression contains only leaded variables */
+ virtual int maxLag() const = 0;
+
//! Returns a new expression where all the leads/lags have been shifted backwards by the same amount
/*!
Only acts on endogenous, exogenous, exogenous det
@@ -350,8 +360,7 @@ class ExprNode
//! Type for the substitution map used in the process of creating auxiliary vars for leads >= 2
typedef map<const ExprNode *, const VariableNode *> subst_table_t;
- //! Type for the substitution map used in the process of substituting adl/diff expressions
- // typedef map<const ExprNode *, const UnaryOpNode *> subst_table_diff_t;
+ //! Type for the substitution map used in the process of substituting adl expressions
typedef map<const ExprNode *, const expr_t> subst_table_adl_t;
//! Creates auxiliary endo lead variables corresponding to this expression
@@ -466,7 +475,8 @@ class ExprNode
virtual expr_t substituteAdl() const = 0;
//! Substitute diff operator
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const = 0;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const = 0;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const = 0;
//! Substitute pac_expectation operator
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table) = 0;
@@ -544,6 +554,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual expr_t substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const;
virtual expr_t substituteEndoLagGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -551,7 +562,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual expr_t decreaseLeadsLagsPredeterminedVariables() const;
virtual expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -625,6 +637,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual expr_t substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const;
virtual expr_t substituteEndoLagGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -632,7 +645,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual expr_t decreaseLeadsLagsPredeterminedVariables() const;
virtual expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -726,6 +740,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual expr_t substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const;
//! Creates another UnaryOpNode with the same opcode, but with a possibly different datatree and argument
@@ -735,7 +750,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual expr_t decreaseLeadsLagsPredeterminedVariables() const;
virtual expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -845,6 +861,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual expr_t substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const;
//! Creates another BinaryOpNode with the same opcode, but with a possibly different datatree and arguments
@@ -854,7 +871,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual expr_t decreaseLeadsLagsPredeterminedVariables() const;
virtual expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -941,6 +959,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual expr_t substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const;
//! Creates another TrinaryOpNode with the same opcode, but with a possibly different datatree and arguments
@@ -950,7 +969,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual expr_t decreaseLeadsLagsPredeterminedVariables() const;
virtual expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -1039,6 +1059,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual expr_t substituteEndoLeadGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool deterministic_model) const;
virtual expr_t substituteEndoLagGreaterThanTwo(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
@@ -1046,7 +1067,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual expr_t buildSimilarExternalFunctionNode(vector<expr_t> &alt_args, DataTree &alt_datatree) const = 0;
virtual expr_t decreaseLeadsLagsPredeterminedVariables() const;
@@ -1220,6 +1242,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual void prepareForDerivation();
virtual expr_t computeDerivative(int deriv_id);
@@ -1233,7 +1256,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual pair<int, expr_t> normalizeEquation(int symb_id_endo, vector<pair<int, pair<expr_t, expr_t> > > &List_of_Op_RHS) const;
virtual void compile(ostream &CompileCode, unsigned int &instruction_number,
@@ -1295,6 +1319,7 @@ public:
virtual int maxEndoLag() const;
virtual int maxExoLag() const;
virtual int maxLead() const;
+ virtual int maxLag() const;
virtual expr_t decreaseLeadsLags(int n) const;
virtual void prepareForDerivation();
virtual expr_t computeDerivative(int deriv_id);
@@ -1308,7 +1333,8 @@ public:
virtual expr_t substituteExoLag(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substituteExpectation(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs, bool partial_information_model) const;
virtual expr_t substituteAdl() const;
- virtual expr_t substituteDiff(subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
+ virtual void findDiffNodes(DataTree &static_datatree, diff_table_t &diff_table) const;
+ virtual expr_t substituteDiff(DataTree &static_datatree, diff_table_t &diff_table, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const;
virtual expr_t substitutePacExpectation(map<const PacExpectationNode *, const BinaryOpNode *> &subst_table);
virtual pair<int, expr_t> normalizeEquation(int symb_id_endo, vector<pair<int, pair<expr_t, expr_t> > > &List_of_Op_RHS) const;
virtual void compile(ostream &CompileCode, unsigned int &instruction_number,
diff --git a/src/ModFile.cc b/src/ModFile.cc
index eb4a374b2b73e53277f46c15769ee4030dbff03e..b5fd81d87ba1cd0f32bf056b72a0c45d89c5abef 100644
--- a/src/ModFile.cc
+++ b/src/ModFile.cc
@@ -39,6 +39,7 @@ ModFile::ModFile(WarningConsolidation &warnings_arg)
orig_ramsey_dynamic_model(symbol_table, num_constants, external_functions_table),
static_model(symbol_table, num_constants, external_functions_table),
steady_state_model(symbol_table, num_constants, external_functions_table, static_model),
+ diff_static_model(symbol_table, num_constants, external_functions_table),
linear(false), block(false), byte_code(false), use_dll(false), no_static(false),
differentiate_forward_vars(false), nonstationary_variables(false),
param_used_with_lead_lag(false), warnings(warnings_arg)
@@ -362,7 +363,7 @@ ModFile::transformPass(bool nostrict, bool stochastic, bool compute_xrefs, const
}
// Create auxiliary variable and equations for Diff operator
- dynamic_model.substituteDiff();
+ dynamic_model.substituteDiff(diff_static_model);
// Var Model
map<string, pair<SymbolList, int> > var_model_info_var_expectation;
@@ -478,6 +479,8 @@ ModFile::transformPass(bool nostrict, bool stochastic, bool compute_xrefs, const
dynamic_model.substituteEndoLagGreaterThanTwo(true);
}
+ dynamic_model.combineDiffAuxEquations();
+
if (differentiate_forward_vars)
dynamic_model.differentiateForwardVars(differentiate_forward_vars_subset);
diff --git a/src/ModFile.hh b/src/ModFile.hh
index ac678b5b52ba79970f086d25b72a0c50878e3410..08ee6c4213b839a854b705ab1b0acf70e43ae3e1 100644
--- a/src/ModFile.hh
+++ b/src/ModFile.hh
@@ -72,6 +72,8 @@ public:
StaticModel static_model;
//! Static model, as declared in the "steady_state_model" block if present
SteadyStateModel steady_state_model;
+ //! Static model used for mapping arguments of diff operator
+ StaticModel diff_static_model;
//! Option linear
bool linear;
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index be6b66cbbf7646be6c35e1fa1e2414eadf836bd8..312169a1f90bfc7d3af1d8d3f2f14f76cfe6f1ed 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -60,6 +60,12 @@ protected:
//! Only stores generated auxiliary equations, in an order meaningful for evaluation
deque<BinaryOpNode *> aux_equations;
+ //! Temporarily stores aux equations for diff operator
+ //! This is a hack to make diff aux vars work: they must be substituted before we consider VARs
+ //! But leads lags can't be substituted until after we consider VARs
+ //! The diff_aux_vars may depend on aux_vars created for leads lags, as in
+ //! diff(x(-1)) => x(-1) - x(-2) => x(-1) - aux_var(-1); aux_var = x(-1);
+ deque<BinaryOpNode *> diff_aux_equations;
//! Stores equation tags
vector<pair<int, pair<string, string> > > equation_tags;
diff --git a/src/StaticModel.cc b/src/StaticModel.cc
index 11a033a867280b3ce747d9ef8319d11c4bd58d43..1ef7bad5a0e1e0a7ace36d9577f1f00fd97216ba 100644
--- a/src/StaticModel.cc
+++ b/src/StaticModel.cc
@@ -2162,6 +2162,8 @@ StaticModel::writeSetAuxiliaryVariables(const string &basename, const bool julia
<< comment << " Warning : this file is generated automatically by Dynare" << endl
<< comment << " from model file (.mod)" << endl << endl
<< output_func_body.str();
+
+ output.close();
}
void
diff --git a/src/SymbolTable.cc b/src/SymbolTable.cc
index ad61e3749c2b089a88d80af2b86d7bfd2cd7365d..3220685771f1013dafdb1bf21f7b23f37647f6e2 100644
--- a/src/SymbolTable.cc
+++ b/src/SymbolTable.cc
@@ -732,7 +732,7 @@ SymbolTable::addDiffAuxiliaryVar(int index, expr_t expr_arg, int orig_symb_id, i
int
SymbolTable::addDiffAuxiliaryVar(int index, expr_t expr_arg) throw (FrozenException)
{
- return addDiffAuxiliaryVar(index, expr_arg, 0, 0);
+ return addDiffAuxiliaryVar(index, expr_arg, -1, 0);
}
int