diff --git a/src/DynamicModel.cc b/src/DynamicModel.cc index 5717f017..f4c7a317 100644 --- a/src/DynamicModel.cc +++ b/src/DynamicModel.cc @@ -3568,15 +3568,15 @@ DynamicModel::fillTrendComponentModelTable() const trend_component_model_table.setVals(eqnums, trend_eqnums, lhsr, lhs_expr_tr); } -pair, expr_t>>, map, expr_t>>> +pair, expr_t>>, map, expr_t>>> DynamicModel::computeErrorComponentMatrices(const ExprNode::subst_table_t &diff_subst_table) const { - map, expr_t>> A0r, A0starr; + map, expr_t>> A0r, A0starr; for (const auto &[model_name, eqns] : trend_component_model_table.getEqNums()) { int i = 0; - map, expr_t> A0, A0star; + map, expr_t> A0, A0star; vector target_lhs = trend_component_model_table.getTargetLhs(model_name); vector nontarget_eqnums = trend_component_model_table.getNonTargetEqNums(model_name); vector undiff_nontarget_lhs = getUndiffLHSForPac(model_name, diff_subst_table); diff --git a/src/DynamicModel.hh b/src/DynamicModel.hh index e0776c00..420dbe2a 100644 --- a/src/DynamicModel.hh +++ b/src/DynamicModel.hh @@ -294,7 +294,7 @@ private: //! Compute error component matrices of trend component_models /*! Returns a pair (A0r, A0starr) */ - pair, expr_t>>, map, expr_t>>> computeErrorComponentMatrices(const ExprNode::subst_table_t &diff_subst_table) const; + pair, expr_t>>, map, expr_t>>> computeErrorComponentMatrices(const ExprNode::subst_table_t &diff_subst_table) const; /* For a VAR model, given the symbol ID of a LHS variable, and a (negative) lag, returns all the corresponding deriv_ids (by properly dealing with two diff --git a/src/ExprNode.cc b/src/ExprNode.cc index b7975f66..48654d7e 100644 --- a/src/ExprNode.cc +++ b/src/ExprNode.cc @@ -290,8 +290,8 @@ void ExprNode::fillErrorCorrectionRow(int eqn, const vector &nontarget_lhs, const vector &target_lhs, - map, expr_t> &A0, - map, expr_t> &A0star) const + map, expr_t> &A0, + map, expr_t> &A0star) const { vector> terms; decomposeAdditiveTerms(terms, 1); @@ -330,6 +330,11 @@ ExprNode::fillErrorCorrectionRow(int eqn, auto [orig_vid, orig_lag] = datatree.symbol_table.unrollDiffLeadLagChain(var_id, lag); if (find(target_lhs.begin(), target_lhs.end(), orig_vid) == target_lhs.end()) { + if (orig_lag != -1) + { + cerr << "ERROR in trend component model: variables in the error correction term should appear with a lag of -1" << endl; + exit(EXIT_FAILURE); + } // This an LHS variable, so fill A0 if (constant != 1) { @@ -342,13 +347,13 @@ ExprNode::fillErrorCorrectionRow(int eqn, exit(EXIT_FAILURE); } int colidx = static_cast(distance(nontarget_lhs.begin(), find(nontarget_lhs.begin(), nontarget_lhs.end(), orig_vid))); - if (A0.find({eqn, -orig_lag, colidx}) != A0.end()) + if (A0.find({eqn, colidx}) != A0.end()) { cerr << "ExprNode::fillErrorCorrection: Error filling A0 matrix: " - << "lag/symb_id encountered more than once in equation" << endl; + << "symb_id encountered more than once in equation" << endl; exit(EXIT_FAILURE); } - A0[{eqn, -orig_lag, colidx}] = datatree.AddVariable(m.first); + A0[{eqn, colidx}] = datatree.AddVariable(m.first); } else { @@ -357,7 +362,7 @@ ExprNode::fillErrorCorrectionRow(int eqn, expr_t e = datatree.AddTimes(datatree.AddVariable(m.first), datatree.AddPossiblyNegativeConstant(-constant)); if (param_id != -1) e = datatree.AddTimes(e, datatree.AddVariable(param_id)); - if (auto coor = tuple(eqn, -orig_lag, colidx); A0star.find(coor) == A0star.end()) + if (auto coor = make_pair(eqn, colidx); A0star.find(coor) == A0star.end()) A0star[coor] = e; else A0star[coor] = datatree.AddPlus(e, A0star[coor]); diff --git a/src/ExprNode.hh b/src/ExprNode.hh index a2ffa392..6b9c0051 100644 --- a/src/ExprNode.hh +++ b/src/ExprNode.hh @@ -669,8 +669,8 @@ public: //! Fills the EC matrix structure void fillErrorCorrectionRow(int eqn, const vector &nontarget_lhs, const vector &target_lhs, - map, expr_t> &A0, - map, expr_t> &A0star) const; + map, expr_t> &A0, + map, expr_t> &A0star) const; //! Replaces variables found in BinaryOpNode::findConstantEquations() with their constant values virtual expr_t replaceVarsInEquation(map &table) const = 0; diff --git a/src/SubModel.cc b/src/SubModel.cc index 55e7d269..b6b764a1 100644 --- a/src/SubModel.cc +++ b/src/SubModel.cc @@ -113,8 +113,8 @@ TrendComponentModelTable::setAR(map, expr_t>> A } void -TrendComponentModelTable::setA0(map, expr_t>> A0_arg, - map, expr_t>> A0star_arg) +TrendComponentModelTable::setA0(map, expr_t>> A0_arg, + map, expr_t>> A0star_arg) { A0 = move(A0_arg); A0star = move(A0star_arg); @@ -361,30 +361,24 @@ TrendComponentModelTable::writeOutput(const string &basename, ostream &output) c ar_ec_output << ";" << endl; } - int a0_lag = 0; - for (const auto &[key, expr] : A0.at(name)) - a0_lag = max(a0_lag, get<1>(key)); ar_ec_output << endl << " % A0" << endl - << " A0 = zeros(" << nontarget_lhs_vec.size() << ", " << nontarget_lhs_vec.size() << ", " << a0_lag << ");" << endl; + << " A0 = zeros(" << nontarget_lhs_vec.size() << ", " << nontarget_lhs_vec.size() << ");" << endl; for (const auto &[key, expr] : A0.at(name)) { - auto [eqn, lag, colidx] = key; - ar_ec_output << " A0(" << eqn + 1 << ", " << colidx + 1 << ", " << lag << ") = "; + auto [eqn, colidx] = key; + ar_ec_output << " A0(" << eqn + 1 << ", " << colidx + 1 << ") = "; expr->writeOutput(ar_ec_output, ExprNodeOutputType::matlabDynamicModel); ar_ec_output << ";" << endl; } - int a0star_lag = 0; - for (const auto &[key, expr] : A0star.at(name)) - a0star_lag = max(a0star_lag, get<1>(key)); ar_ec_output << endl << " % A0star" << endl - << " A0star = zeros(" << nontarget_lhs_vec.size() << ", " << target_lhs_vec.size() << ", " << a0star_lag << ");" << endl; + << " A0star = zeros(" << nontarget_lhs_vec.size() << ", " << target_lhs_vec.size() << ");" << endl; for (const auto &[key, expr] : A0star.at(name)) { - auto [eqn, lag, colidx] = key; - ar_ec_output << " A0star(" << eqn + 1 << ", " << colidx + 1 << ", " << lag << ") = "; + auto [eqn, colidx] = key; + ar_ec_output << " A0star(" << eqn + 1 << ", " << colidx + 1 << ") = "; expr->writeOutput(ar_ec_output, ExprNodeOutputType::matlabDynamicModel); ar_ec_output << ";" << endl; } diff --git a/src/SubModel.hh b/src/SubModel.hh index c035fbae..3bb42633 100644 --- a/src/SubModel.hh +++ b/src/SubModel.hh @@ -53,7 +53,7 @@ private: map, expr_t>> AR; // name -> (eqn, lag, lhs_symb_id) -> expr_t /* Note that A0 in the trend-component model context is not the same thing as in the structural VAR context. */ - map, expr_t>> A0, A0star; // name -> (eqn, lag, col) -> expr_t + map, expr_t>> A0, A0star; // name -> (eqn, col) -> expr_t public: explicit TrendComponentModelTable(SymbolTable &symbol_table_arg); @@ -91,8 +91,8 @@ public: void setOrigDiffVar(map> orig_diff_var_arg); void setTargetVar(map> target_vars_arg); void setAR(map, expr_t>> AR_arg); - void setA0(map, expr_t>> A0_arg, - map, expr_t>> A0star_arg); + void setA0(map, expr_t>> A0_arg, + map, expr_t>> A0star_arg); //! Write output of this class void writeOutput(const string &basename, ostream &output) const;