From ea6fb40db7f0d6a8bd9638e89f392a0b9782399c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Tue, 26 Oct 2021 18:06:26 +0200
Subject: [PATCH] =?UTF-8?q?PAC:=20new=20=E2=80=9Cpac=5Ftarget=5Finfo?=
=?UTF-8?q?=E2=80=9D=20block=20and=20=E2=80=9Cpac=5Ftarget=5Fnonstationary?=
=?UTF-8?q?=E2=80=9D=20operator?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit
Ref. Madysson/estimation-codes#5
---
src/CodeInterpreter.hh | 8 +
src/DataTree.cc | 15 ++
src/DataTree.hh | 20 +-
src/DynamicModel.cc | 93 ++++++++
src/DynamicModel.hh | 22 ++
src/DynareBison.yy | 53 +++++
src/DynareFlex.ll | 25 ++-
src/ExprNode.cc | 261 +++++++++++++++++++++++
src/ExprNode.hh | 60 +++++-
src/ModFile.cc | 6 +-
src/ModelTree.cc | 9 +
src/ModelTree.hh | 4 +
src/ParsingDriver.cc | 60 ++++++
src/ParsingDriver.hh | 13 ++
src/SubModel.cc | 466 +++++++++++++++++++++++++++++++++++------
src/SubModel.hh | 38 +++-
src/SymbolTable.cc | 20 ++
src/SymbolTable.hh | 5 +-
18 files changed, 1099 insertions(+), 79 deletions(-)
diff --git a/src/CodeInterpreter.hh b/src/CodeInterpreter.hh
index 899fb014..efdb0367 100644
--- a/src/CodeInterpreter.hh
+++ b/src/CodeInterpreter.hh
@@ -241,6 +241,14 @@ enum class PriorDistributions
weibull = 8
};
+enum class PacTargetKind
+ {
+ unspecified, // Must be the first one, because it’s the default initializer
+ ll,
+ dl,
+ dd
+ };
+
struct Block_contain_type
{
int Equation, Variable, Own_Derivative;
diff --git a/src/DataTree.cc b/src/DataTree.cc
index d460ad08..654ff9d5 100644
--- a/src/DataTree.cc
+++ b/src/DataTree.cc
@@ -93,6 +93,7 @@ DataTree::operator=(const DataTree &d)
external_function_node_map.clear();
var_expectation_node_map.clear();
pac_expectation_node_map.clear();
+ pac_target_nonstationary_node_map.clear();
first_deriv_external_function_node_map.clear();
second_deriv_external_function_node_map.clear();
@@ -704,6 +705,20 @@ DataTree::AddPacExpectation(const string &model_name)
return p;
}
+expr_t
+DataTree::AddPacTargetNonstationary(const string &model_name)
+{
+ if (auto it = pac_target_nonstationary_node_map.find(model_name);
+ it != pac_target_nonstationary_node_map.end())
+ return it->second;
+
+ auto sp = make_unique<PacTargetNonstationaryNode>(*this, node_list.size(), model_name);
+ auto p = sp.get();
+ node_list.push_back(move(sp));
+ pac_target_nonstationary_node_map[model_name] = p;
+ return p;
+}
+
BinaryOpNode *
DataTree::AddEqual(expr_t iArg1, expr_t iArg2)
{
diff --git a/src/DataTree.hh b/src/DataTree.hh
index e7baa641..565a0189 100644
--- a/src/DataTree.hh
+++ b/src/DataTree.hh
@@ -81,6 +81,10 @@ private:
using pac_expectation_node_map_t = map<string, PacExpectationNode *>;
pac_expectation_node_map_t pac_expectation_node_map;
+ // model_name -> PacTargetNonstationaryNode
+ using pac_target_nonstationary_node_map_t = map<string, PacTargetNonstationaryNode *>;
+ pac_target_nonstationary_node_map_t pac_target_nonstationary_node_map;
+
// (arguments, deriv_idx, symb_id) -> FirstDerivExternalFunctionNode
using first_deriv_external_function_node_map_t = map<tuple<vector<expr_t>, int, int>, FirstDerivExternalFunctionNode *>;
first_deriv_external_function_node_map_t first_deriv_external_function_node_map;
@@ -101,13 +105,6 @@ protected:
//! Internal implementation of ParamUsedWithLeadLag()
bool ParamUsedWithLeadLagInternal() const;
- /*! Takes a MATLAB/Octave package name (possibly with several levels nested using dots),
- and returns the name of the corresponding filesystem directory (which
- is created by the function if it does not exist).
- In practice the package nesting is used for the planner_objective (stored
- inside +objective subdir). */
- static string packageDir(const string &package);
-
/* Writes the contents of “new_contents” to the file “filename”. However, if
the file already exists and would not be modified by this operation, then do
nothing. */
@@ -260,6 +257,8 @@ public:
expr_t AddVarExpectation(const string &model_name);
//! Adds pac_expectation command to model tree
expr_t AddPacExpectation(const string &model_name);
+ //! Adds a pac_target_nonstationary node to model tree
+ expr_t AddPacTargetNonstationary(const string &model_name);
//! Adds a model local variable with its value
void AddLocalVariable(int symb_id, expr_t value) noexcept(false);
//! Adds an external function node
@@ -344,6 +343,13 @@ public:
{
no_commutativity = true;
}
+
+ /*! Takes a MATLAB/Octave package name (possibly with several levels nested using dots),
+ and returns the name of the corresponding filesystem directory (which
+ is created by the function if it does not exist).
+ In practice the package nesting is used for the planner_objective (stored
+ inside +objective subdir). */
+ static string packageDir(const string &package);
};
inline expr_t
diff --git a/src/DynamicModel.cc b/src/DynamicModel.cc
index 39a6ce52..6a4885da 100644
--- a/src/DynamicModel.cc
+++ b/src/DynamicModel.cc
@@ -4080,6 +4080,80 @@ DynamicModel::computePacBackwardExpectationSubstitution(const string &name,
pac_expectation_substitution[name] = AddVariable(expect_var_id);
}
+void
+DynamicModel::computePacBackwardExpectationSubstitutionWithComponents(const string &name,
+ const vector<int> &lhs,
+ int max_lag,
+ const string &aux_model_type,
+ vector<PacModelTable::target_component_t> &pac_target_components,
+ map<string, expr_t> &pac_expectation_substitution)
+{
+ auto create_aux_param = [&](const string ¶m_name)
+ {
+ try
+ {
+ return symbol_table.addSymbol(param_name, SymbolType::parameter);
+ }
+ catch (SymbolTable::AlreadyDeclaredException)
+ {
+ cerr << "ERROR: the variable/parameter '" << param_name << "' conflicts with some auxiliary parameter that will be generated for the '" << name << "' PAC model. Please rename that parameter." << endl;
+ exit(EXIT_FAILURE);
+ }
+ };
+
+ expr_t substexpr = Zero;
+ int component_idx = 1;
+
+ for (auto &[component, growth, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth, growth_info] : pac_target_components)
+ {
+ string name_component = name + "_component" + to_string(component_idx);
+
+ // Create the linear combination of the variables from the auxiliary model
+ expr_t auxdef = Zero;
+ if (aux_model_type == "var")
+ {
+ /* If the auxiliary model is a VAR, add a parameter corresponding
+ to the constant. */
+ int new_param_symb_id = create_aux_param("h_" + name_component + "_constant");
+ h_indices.push_back(new_param_symb_id);
+ auxdef = AddPlus(auxdef, AddVariable(new_param_symb_id));
+ }
+ for (int i = 1; i < max_lag + 1; i++)
+ for (auto lhsit : lhs)
+ {
+ int new_param_symb_id = create_aux_param("h_" + name_component + "_var_" + symbol_table.getName(lhsit) + "_lag_" + to_string(i));
+ h_indices.push_back(new_param_symb_id);
+ auxdef = AddPlus(auxdef, AddTimes(AddVariable(new_param_symb_id),
+ AddVariable(lhsit, -i)));
+ }
+
+ // If needed, add the growth neutrality correction for this component
+ if (growth)
+ {
+ growth_neutrality_param = create_aux_param(name_component + "_pac_growth_neutrality_correction");
+ auxdef = AddPlus(auxdef, AddTimes(growth, AddVariable(growth_neutrality_param)));
+ }
+ else
+ growth_neutrality_param = -1;
+
+ // Create the auxiliary variable for this component
+ int aux_id = symbol_table.addPacExpectationAuxiliaryVar(auxname, auxdef);
+ expr_t auxvar = AddVariable(aux_id);
+
+ // Add the equation defining the auxiliary variable for this component
+ expr_t neweq = AddEqual(auxvar, auxdef);
+ addEquation(neweq, -1);
+ addAuxEquation(neweq);
+
+ // Update the expression to be substituted for the pac_expectation operator
+ substexpr = AddPlus(substexpr, AddTimes(coeff, auxvar));
+
+ component_idx++;
+ }
+
+ pac_expectation_substitution[name] = substexpr;
+}
+
void
DynamicModel::substitutePacExpectation(const map<string, expr_t> &pac_expectation_substitution,
const map<string, string> &pac_eq_name)
@@ -4093,6 +4167,13 @@ DynamicModel::substitutePacExpectation(const map<string, expr_t> &pac_expectatio
}
}
+void
+DynamicModel::substitutePacTargetNonstationary(const string &pac_model_name, expr_t substexpr)
+{
+ for (auto &eq : equations)
+ eq = dynamic_cast<BinaryOpNode *>(eq->substitutePacTargetNonstationary(pac_model_name, substexpr));
+}
+
void
DynamicModel::computingPass(bool jacobianExo, int derivsOrder, int paramsDerivsOrder,
const eval_context_t &eval_context, bool no_tmp_terms, bool block, bool use_dll,
@@ -6482,3 +6563,15 @@ DynamicModel::simplifyEquations()
findConstantEquationsWithoutMcpTag(subst_table);
}
}
+
+void
+DynamicModel::checkNoRemainingPacTargetNonstationary() const
+{
+ for (size_t eq = 0; eq < equations.size(); eq++)
+ if (equations[eq]->containsPacTargetNonstationary())
+ {
+ cerr << "ERROR: in equation " << equation_tags.getTagValueByEqnAndKey(eq, "name")
+ << ", the pac_target_nonstationary operator does not match a corresponding 'pac_target_info' block" << endl;
+ exit(EXIT_FAILURE);
+ }
+}
diff --git a/src/DynamicModel.hh b/src/DynamicModel.hh
index 05584e42..6907bcda 100644
--- a/src/DynamicModel.hh
+++ b/src/DynamicModel.hh
@@ -582,10 +582,28 @@ public:
map<string, vector<int>> &pac_aux_param_symb_ids,
map<string, expr_t> &pac_expectation_substitution);
+ /* Same as above, but for PAC models which have an associated
+ pac_target_info.
+ Contrary to the above routine, this one will create the growth correction
+ parameters as needed.
+ Those parameter IDs, as well as the IDs for the h parameters, are stored
+ in target_components.
+ The routine also creates the auxiliary variables for the components, and
+ adds the corresponding equations. */
+ void computePacBackwardExpectationSubstitutionWithComponents(const string &name,
+ const vector<int> &lhs,
+ int max_lag,
+ const string &aux_model_type,
+ vector<PacModelTable::target_component_t> &pac_target_components,
+ map<string, expr_t> &pac_expectation_substitution);
+
//! Substitutes pac_expectation operator with expectation based on auxiliary model
void substitutePacExpectation(const map<string, expr_t> &pac_expectation_substitution,
const map<string, string> &pac_eq_name);
+ //! Substitutes the pac_target_nonstationary operator of a given pac_model
+ void substitutePacTargetNonstationary(const string &pac_model_name, expr_t substexpr);
+
//! Table to undiff LHS variables for pac vector z
vector<int> getUndiffLHSForPac(const string &aux_model_name,
const ExprNode::subst_table_t &diff_subst_table) const;
@@ -605,6 +623,10 @@ public:
out otherwise */
void checkNoRemainingPacExpectation() const;
+ /*! Checks that all pac_target_nonstationary operators have been substituted, error
+ out otherwise */
+ void checkNoRemainingPacTargetNonstationary() const;
+
auto
getStaticOnlyEquationsInfo() const
{
diff --git a/src/DynareBison.yy b/src/DynareBison.yy
index 1d2dde48..87c0f9ba 100644
--- a/src/DynareBison.yy
+++ b/src/DynareBison.yy
@@ -174,6 +174,8 @@ class ParsingDriver;
%token RANDOM_FUNCTION_CONVERGENCE_CRITERION RANDOM_PARAMETER_CONVERGENCE_CRITERION NO_INIT_ESTIMATION_CHECK_FIRST_OBS
%token HETEROSKEDASTIC_FILTER TIME_SHIFT STRUCTURAL TERMINAL_STEADY_STATE_AS_GUESS_VALUE
%token SURPRISE OCCBIN_CONSTRAINTS
+%token PAC_TARGET_INFO COMPONENT TARGET AUXNAME AUXNAME_TARGET_NONSTATIONARY PAC_TARGET_NONSTATIONARY
+%token <string> KIND LL DL DD
/* Method of Moments */
%token METHOD_OF_MOMENTS MOM_METHOD
%token BARTLETT_KERNEL_LAG WEIGHTING_MATRIX WEIGHTING_MATRIX_SCALING_FACTOR ANALYTIC_STANDARD_ERRORS ANALYTIC_JACOBIAN PENALIZED_ESTIMATOR VERBOSE
@@ -214,6 +216,7 @@ class ParsingDriver;
%type <vector<tuple<string, BinaryOpNode *, BinaryOpNode *, expr_t, expr_t>>> occbin_constraints_regimes_list
%type <map<string, expr_t>> occbin_constraints_regime_options_list
%type <pair<string, expr_t>> occbin_constraints_regime_option
+%type <PacTargetKind> pac_target_kind
%%
%start statement_list;
@@ -348,6 +351,7 @@ statement : parameters
| model_replace
| model_options
| var_remove
+ | pac_target_info
;
dsample : DSAMPLE INT_NUMBER ';'
@@ -945,6 +949,49 @@ occbin_constraints_regime_option : BIND hand_side ';'
{ $$ = { "error_relax", $2 }; }
;
+pac_target_info : PAC_TARGET_INFO '(' symbol ')' ';'
+ { driver.begin_pac_target_info($3); }
+ pac_target_info_statement_list
+ END ';'
+ { driver.end_pac_target_info(); }
+ ;
+
+pac_target_info_statement_list : pac_target_info_statement
+ | pac_target_info_statement_list pac_target_info_statement
+ ;
+
+pac_target_info_statement : TARGET hand_side ';'
+ { driver.set_pac_target_info_target($2); }
+ | AUXNAME_TARGET_NONSTATIONARY symbol ';'
+ { driver.set_pac_target_info_auxname_target_nonstationary($2); }
+ | pac_target_info_component
+ ;
+
+pac_target_info_component : COMPONENT hand_side ';'
+ pac_target_info_component_list
+ { driver.add_pac_target_info_component($2); }
+ ;
+
+pac_target_info_component_list : pac_target_info_component_elem
+ | pac_target_info_component_list pac_target_info_component_elem
+ ;
+
+pac_target_info_component_elem : GROWTH hand_side ';'
+ { driver.set_pac_target_info_component_growth($2); }
+ | AUXNAME symbol ';'
+ { driver.set_pac_target_info_component_auxname($2); }
+ | KIND pac_target_kind ';'
+ { driver.set_pac_target_info_component_kind($2); }
+ ;
+
+pac_target_kind : LL
+ { $$ = PacTargetKind::ll; }
+ | DL
+ { $$ = PacTargetKind::dl; }
+ | DD
+ { $$ = PacTargetKind::dd; }
+ ;
+
/* The tokens below must be accepted in both DYNARE_STATEMENT and DYNARE_BLOCK
states in the lexer, because of model block and model_options statement */
model_option : BLOCK { driver.block(); }
@@ -1036,6 +1083,8 @@ hand_side : '(' hand_side ')'
{ $$ = driver.add_var_expectation($3); }
| PAC_EXPECTATION '(' symbol ')'
{ $$ = driver.add_pac_expectation($3); }
+ | PAC_TARGET_NONSTATIONARY '(' symbol ')'
+ { $$ = driver.add_pac_target_nonstationary($3); }
| MINUS hand_side %prec UNARY
{ $$ = driver.add_uminus($2); }
| PLUS hand_side %prec UNARY
@@ -4301,6 +4350,10 @@ symbol : NAME
| RELAX
| ERROR_BIND
| ERROR_RELAX
+ | KIND
+ | LL
+ | DL
+ | DD
;
%%
diff --git a/src/DynareFlex.ll b/src/DynareFlex.ll
index 1bc54b01..f9be86f9 100644
--- a/src/DynareFlex.ll
+++ b/src/DynareFlex.ll
@@ -236,6 +236,7 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
<INITIAL>matched_moments {BEGIN DYNARE_BLOCK; return token::MATCHED_MOMENTS;}
<INITIAL>occbin_constraints {BEGIN DYNARE_BLOCK; return token::OCCBIN_CONSTRAINTS;}
<INITIAL>model_replace {BEGIN DYNARE_BLOCK; return token::MODEL_REPLACE;}
+<INITIAL>pac_target_info {BEGIN DYNARE_BLOCK; return token::PAC_TARGET_INFO;}
/* For the semicolon after an "end" keyword */
<INITIAL>; {return Dynare::parser::token_type (yytext[0]);}
@@ -684,7 +685,7 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
<DYNARE_STATEMENT>epilogue {return token::EPILOGUE;}
<DYNARE_STATEMENT>growth_factor {return token::GROWTH_FACTOR;}
<DYNARE_STATEMENT>log_growth_factor {return token::LOG_GROWTH_FACTOR;}
-<DYNARE_STATEMENT>growth {return token::GROWTH;}
+<DYNARE_STATEMENT,DYNARE_BLOCK>growth {return token::GROWTH;}
<DYNARE_STATEMENT>cova_compute {return token::COVA_COMPUTE;}
<DYNARE_STATEMENT>discretionary_tol {return token::DISCRETIONARY_TOL;}
<DYNARE_STATEMENT>analytic_derivation {return token::ANALYTIC_DERIVATION;}
@@ -800,6 +801,27 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
<DYNARE_BLOCK># {return Dynare::parser::token_type (yytext[0]);}
<DYNARE_BLOCK>restriction {return token::RESTRICTION;}
+<DYNARE_BLOCK>component {return token::COMPONENT;}
+<DYNARE_BLOCK>target {return token::TARGET;}
+<DYNARE_BLOCK>auxname {return token::AUXNAME;}
+<DYNARE_BLOCK>auxname_target_nonstationary {return token::AUXNAME_TARGET_NONSTATIONARY;}
+<DYNARE_BLOCK>kind {
+ yylval->build<string>(yytext);
+ return token::KIND;
+}
+<DYNARE_BLOCK>ll {
+ yylval->build<string>(yytext);
+ return token::LL;
+}
+<DYNARE_BLOCK>dl {
+ yylval->build<string>(yytext);
+ return token::DL;
+}
+<DYNARE_BLOCK>dd {
+ yylval->build<string>(yytext);
+ return token::DD;
+}
+
/* Inside Dynare statement */
<DYNARE_STATEMENT>solve_algo {return token::SOLVE_ALGO;}
@@ -945,6 +967,7 @@ DATE -?[0-9]+([ya]|m([1-9]|1[0-2])|q[1-4])
<DYNARE_STATEMENT,DYNARE_BLOCK>expectation {return token::EXPECTATION;}
<DYNARE_BLOCK>var_expectation {return token::VAR_EXPECTATION;}
<DYNARE_BLOCK>pac_expectation {return token::PAC_EXPECTATION;}
+<DYNARE_BLOCK>pac_target_nonstationary {return token::PAC_TARGET_NONSTATIONARY;}
<DYNARE_STATEMENT>discount {return token::DISCOUNT;}
<DYNARE_STATEMENT,DYNARE_BLOCK>varobs {return token::VAROBS;}
<DYNARE_STATEMENT,DYNARE_BLOCK>varexobs {return token::VAREXOBS;}
diff --git a/src/ExprNode.cc b/src/ExprNode.cc
index 6ce01c80..cc7b73fb 100644
--- a/src/ExprNode.cc
+++ b/src/ExprNode.cc
@@ -667,6 +667,12 @@ NumConstNode::substitutePacExpectation(const string &name, expr_t subexpr)
return const_cast<NumConstNode *>(this);
}
+expr_t
+NumConstNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ return const_cast<NumConstNode *>(this);
+}
+
expr_t
NumConstNode::differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
@@ -694,6 +700,12 @@ NumConstNode::containsPacExpectation(const string &pac_model_name) const
return false;
}
+bool
+NumConstNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ return false;
+}
+
expr_t
NumConstNode::replaceTrendVar() const
{
@@ -1588,6 +1600,15 @@ VariableNode::substitutePacExpectation(const string &name, expr_t subexpr)
return const_cast<VariableNode *>(this);
}
+expr_t
+VariableNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ if (get_type() == SymbolType::modelLocalVariable)
+ return datatree.getLocalVariable(symb_id)->substitutePacTargetNonstationary(name, subexpr);
+
+ return const_cast<VariableNode *>(this);
+}
+
expr_t
VariableNode::decreaseLeadsLags(int n) const
{
@@ -1820,6 +1841,15 @@ VariableNode::containsPacExpectation(const string &pac_model_name) const
return false;
}
+bool
+VariableNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ if (get_type() == SymbolType::modelLocalVariable)
+ return datatree.getLocalVariable(symb_id)->containsPacTargetNonstationary(pac_model_name);
+
+ return false;
+}
+
expr_t
VariableNode::replaceTrendVar() const
{
@@ -3538,6 +3568,13 @@ UnaryOpNode::substitutePacExpectation(const string &name, expr_t subexpr)
return buildSimilarUnaryOpNode(argsubst, datatree);
}
+expr_t
+UnaryOpNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ expr_t argsubst = arg->substitutePacTargetNonstationary(name, subexpr);
+ return buildSimilarUnaryOpNode(argsubst, datatree);
+}
+
expr_t
UnaryOpNode::decreaseLeadsLags(int n) const
{
@@ -3665,6 +3702,12 @@ UnaryOpNode::containsPacExpectation(const string &pac_model_name) const
return arg->containsPacExpectation(pac_model_name);
}
+bool
+UnaryOpNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ return arg->containsPacTargetNonstationary(pac_model_name);
+}
+
expr_t
UnaryOpNode::replaceTrendVar() const
{
@@ -5086,6 +5129,14 @@ BinaryOpNode::substitutePacExpectation(const string &name, expr_t subexpr)
return buildSimilarBinaryOpNode(arg1subst, arg2subst, datatree);
}
+expr_t
+BinaryOpNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ expr_t arg1subst = arg1->substitutePacTargetNonstationary(name, subexpr);
+ expr_t arg2subst = arg2->substitutePacTargetNonstationary(name, subexpr);
+ return buildSimilarBinaryOpNode(arg1subst, arg2subst, datatree);
+}
+
expr_t
BinaryOpNode::differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
@@ -5120,6 +5171,12 @@ BinaryOpNode::containsPacExpectation(const string &pac_model_name) const
return arg1->containsPacExpectation(pac_model_name) || arg2->containsPacExpectation(pac_model_name);
}
+bool
+BinaryOpNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ return arg1->containsPacTargetNonstationary(pac_model_name) || arg2->containsPacTargetNonstationary(pac_model_name);
+}
+
expr_t
BinaryOpNode::replaceTrendVar() const
{
@@ -6341,6 +6398,15 @@ TrinaryOpNode::substitutePacExpectation(const string &name, expr_t subexpr)
return buildSimilarTrinaryOpNode(arg1subst, arg2subst, arg3subst, datatree);
}
+expr_t
+TrinaryOpNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ expr_t arg1subst = arg1->substitutePacTargetNonstationary(name, subexpr);
+ expr_t arg2subst = arg2->substitutePacTargetNonstationary(name, subexpr);
+ expr_t arg3subst = arg3->substitutePacTargetNonstationary(name, subexpr);
+ return buildSimilarTrinaryOpNode(arg1subst, arg2subst, arg3subst, datatree);
+}
+
expr_t
TrinaryOpNode::differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
@@ -6368,6 +6434,14 @@ TrinaryOpNode::containsPacExpectation(const string &pac_model_name) const
return (arg1->containsPacExpectation(pac_model_name) || arg2->containsPacExpectation(pac_model_name) || arg3->containsPacExpectation(pac_model_name));
}
+bool
+TrinaryOpNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ return arg1->containsPacTargetNonstationary(pac_model_name)
+ || arg2->containsPacTargetNonstationary(pac_model_name)
+ || arg3->containsPacTargetNonstationary(pac_model_name);
+}
+
expr_t
TrinaryOpNode::replaceTrendVar() const
{
@@ -6741,6 +6815,15 @@ AbstractExternalFunctionNode::substitutePacExpectation(const string &name, expr_
return buildSimilarExternalFunctionNode(arguments_subst, datatree);
}
+expr_t
+AbstractExternalFunctionNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ vector<expr_t> arguments_subst;
+ for (auto argument : arguments)
+ arguments_subst.push_back(argument->substitutePacTargetNonstationary(name, subexpr));
+ return buildSimilarExternalFunctionNode(arguments_subst, datatree);
+}
+
expr_t
AbstractExternalFunctionNode::differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const
{
@@ -6835,6 +6918,15 @@ AbstractExternalFunctionNode::containsPacExpectation(const string &pac_model_nam
return false;
}
+bool
+AbstractExternalFunctionNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ for (auto argument : arguments)
+ if (argument->containsPacTargetNonstationary(pac_model_name))
+ return true;
+ return false;
+}
+
expr_t
AbstractExternalFunctionNode::replaceTrendVar() const
{
@@ -8344,12 +8436,24 @@ VarExpectationNode::substitutePacExpectation(const string &name, expr_t subexpr)
return const_cast<VarExpectationNode *>(this);
}
+expr_t
+VarExpectationNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ return const_cast<VarExpectationNode *>(this);
+}
+
bool
VarExpectationNode::containsPacExpectation(const string &pac_model_name) const
{
return false;
}
+bool
+VarExpectationNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ return false;
+}
+
void
VarExpectationNode::writeJsonAST(ostream &output) const
{
@@ -8417,6 +8521,12 @@ PacExpectationNode::containsPacExpectation(const string &pac_model_name) const
return pac_model_name == model_name;
}
+bool
+PacExpectationNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ return false;
+}
+
void
PacExpectationNode::writeJsonAST(ostream &output) const
{
@@ -8443,6 +8553,99 @@ PacExpectationNode::substitutePacExpectation(const string &name, expr_t subexpr)
return subexpr;
}
+expr_t
+PacExpectationNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ return const_cast<PacExpectationNode *>(this);
+}
+
+PacTargetNonstationaryNode::PacTargetNonstationaryNode(DataTree &datatree_arg,
+ int idx_arg,
+ string model_name_arg) :
+ SubModelNode{datatree_arg, idx_arg, move(model_name_arg)}
+{
+}
+
+expr_t
+PacTargetNonstationaryNode::clone(DataTree &datatree) const
+{
+ return datatree.AddPacTargetNonstationary(model_name);
+}
+
+void
+PacTargetNonstationaryNode::writeOutput(ostream &output, ExprNodeOutputType output_type,
+ const temporary_terms_t &temporary_terms,
+ const temporary_terms_idxs_t &temporary_terms_idxs,
+ const deriv_node_temp_terms_t &tef_terms) const
+{
+ assert(output_type != ExprNodeOutputType::matlabOutsideModel);
+ if (isLatexOutput(output_type))
+ {
+ output << "PAC_TARGET_NONSTATIONARY" << LEFT_PAR(output_type) << model_name << RIGHT_PAR(output_type);
+ return;
+ }
+}
+
+int
+PacTargetNonstationaryNode::maxLagWithDiffsExpanded() const
+{
+ // This node will be replaced by the target lagged by one
+ return 1;
+}
+
+expr_t
+PacTargetNonstationaryNode::substituteVarExpectation(const map<string, expr_t> &subst_table) const
+{
+ return const_cast<PacTargetNonstationaryNode *>(this);
+}
+
+bool
+PacTargetNonstationaryNode::containsPacExpectation(const string &pac_model_name) const
+{
+ return false;
+}
+
+bool
+PacTargetNonstationaryNode::containsPacTargetNonstationary(const string &pac_model_name) const
+{
+ if (pac_model_name.empty())
+ return true;
+ else
+ return pac_model_name == model_name;
+}
+
+void
+PacTargetNonstationaryNode::writeJsonAST(ostream &output) const
+{
+ output << R"({"node_type" : "PacTargetNonstationaryNode", )"
+ << R"("name" : ")" << model_name << R"("})";
+}
+
+void
+PacTargetNonstationaryNode::writeJsonOutput(ostream &output,
+ const temporary_terms_t &temporary_terms,
+ const deriv_node_temp_terms_t &tef_terms,
+ bool isdynamic) const
+{
+ output << "pac_target_nonstationary("
+ << "model_name = " << model_name
+ << ")";
+}
+
+expr_t
+PacTargetNonstationaryNode::substitutePacExpectation(const string &name, expr_t subexpr)
+{
+ return const_cast<PacTargetNonstationaryNode *>(this);
+}
+
+expr_t
+PacTargetNonstationaryNode::substitutePacTargetNonstationary(const string &name, expr_t subexpr)
+{
+ if (model_name != name)
+ return const_cast<PacTargetNonstationaryNode *>(this);
+ return subexpr;
+}
+
void
ExprNode::decomposeAdditiveTerms(vector<pair<expr_t, int>> &terms, int current_sign) const
{
@@ -8649,6 +8852,8 @@ ExprNode::matchParamTimesTargetMinusVariable(int symb_id) const
if (datatree.symbol_table.isAuxiliaryVariable(target->symb_id))
{
auto avi = datatree.symbol_table.getAuxVarInfo(target->symb_id);
+ if (avi.get_type() == AuxVarType::pacTargetNonstationary && target->lag == -1)
+ return true;
return (avi.get_type() == AuxVarType::unaryOp
&& avi.get_unary_op() == "log"
&& avi.get_orig_symb_id() != -1
@@ -8664,3 +8869,59 @@ ExprNode::matchParamTimesTargetMinusVariable(int symb_id) const
else
throw MatchFailureException{"Neither factor is of the form (target-variable) where target is endo or exo (possibly logged), and has one lag"};
}
+
+pair<int, expr_t>
+ExprNode::matchEndogenousTimesConstant() const
+{
+ throw MatchFailureException{"This expression is not of the form endogenous*constant"};
+}
+
+pair<int, expr_t>
+VariableNode::matchEndogenousTimesConstant() const
+{
+ if (get_type() == SymbolType::endogenous)
+ return { symb_id, datatree.One };
+ else
+ throw MatchFailureException{"This expression is not of the form endogenous*constant"};
+}
+
+pair<int, expr_t>
+BinaryOpNode::matchEndogenousTimesConstant() const
+{
+ if (op_code == BinaryOpcode::times)
+ {
+ if (auto varg1 = dynamic_cast<VariableNode *>(arg1);
+ varg1 && varg1->get_type() == SymbolType::endogenous && arg2->isConstant())
+ return { varg1->symb_id, arg2 };
+ if (auto varg2 = dynamic_cast<VariableNode *>(arg2);
+ varg2 && varg2->get_type() == SymbolType::endogenous && arg1->isConstant())
+ return { varg2->symb_id, arg1 };
+ }
+ throw MatchFailureException{"This expression is not of the form endogenous*constant"};
+}
+
+pair<vector<pair<int, expr_t>>, expr_t>
+ExprNode::matchLinearCombinationOfEndogenousWithConstant() const
+{
+ vector<pair<expr_t, int>> all_terms;
+ decomposeAdditiveTerms(all_terms);
+
+ vector<pair<int, expr_t>> endo_terms;
+ expr_t intercept = datatree.Zero;
+ for (auto [term, sign] : all_terms)
+ if (term->isConstant())
+ {
+ if (sign == -1)
+ intercept = datatree.AddMinus(intercept, term);
+ else
+ intercept = datatree.AddPlus(intercept, term);
+ }
+ else
+ {
+ auto [endo_id, constant] = term->matchEndogenousTimesConstant();
+ if (sign == -1)
+ constant = datatree.AddUMinus(constant);
+ endo_terms.emplace_back(endo_id, constant);
+ }
+ return { endo_terms, intercept };
+}
diff --git a/src/ExprNode.hh b/src/ExprNode.hh
index 17bf0405..c4b20990 100644
--- a/src/ExprNode.hh
+++ b/src/ExprNode.hh
@@ -612,6 +612,9 @@ public:
//! Substitute pac_expectation operator
virtual expr_t substitutePacExpectation(const string &name, expr_t subexpr) = 0;
+ //! Substitute pac_target_nonstationary operator
+ virtual expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) = 0;
+
virtual int findTargetVariable(int lhs_symb_id) const = 0;
//! Add ExprNodes to the provided datatree
@@ -626,7 +629,7 @@ public:
//! Substitute auxiliary variables by their expression in static model
virtual expr_t substituteStaticAuxiliaryVariable() const = 0;
- //! Matches a linear combination of variables, where scalars can be constant*parameter
+ //! Matches a linear combination of variables (endo or exo), where scalars can be constant*parameter
/*! Returns a list of (variable_id, lag, param_id, constant)
corresponding to the terms in the expression. When there is no
parameter in a term, param_id == -1.
@@ -638,6 +641,15 @@ public:
*/
vector<tuple<int, int, int, double>> matchLinearCombinationOfVariables(bool variable_obligatory_in_each_term = true) const;
+ /* Matches a linear combination of endogenous, where scalars can be any
+ constant expression (i.e. containing no endogenous, no exogenous and no
+ exogenous deterministic). The linear combination can contain constant
+ terms (intercept).
+ Returns a pair composed of:
+ – the terms of the form endogenous*scalar, as a list of (endo_id, constant expr);
+ – the sum of all constant (intercept) terms */
+ pair<vector<pair<int, expr_t>>, expr_t> matchLinearCombinationOfEndogenousWithConstant() const;
+
pair<int, vector<tuple<int, int, int, double>>> matchParamTimesLinearCombinationOfVariables() const;
/* Matches an expression of the form parameter*(var1-endo2).
@@ -663,6 +675,9 @@ public:
//! Returns true if PacExpectationNode encountered
virtual bool containsPacExpectation(const string &pac_model_name = "") const = 0;
+ //! Returns true if PacTargetNonstationaryNode encountered
+ virtual bool containsPacTargetNonstationary(const string &pac_model_name = "") const = 0;
+
//! Decompose an expression into its additive terms
/*! Returns a list of terms, with their sign (either 1 or -1, depending
on whether the terms appears with a plus or a minus).
@@ -691,7 +706,14 @@ public:
*/
tuple<int, int, int, double> matchVariableTimesConstantTimesParam(bool variable_obligatory = true) const;
- //! Exception thrown by matchVariableTimesConstantTimesParam when matching fails
+ /* Matches an expression of the form endogenous*constant where constant is an
+ expression containing no endogenous, no exogenous and no exogenous deterministic.
+ Returns (endo_id, constant expr).
+ Note that it will also match a simple endogenous (in which case the
+ constant will of course be equal to one). */
+ virtual pair<int, expr_t> matchEndogenousTimesConstant() const;
+
+ //! Exception thrown when matching fails
class MatchFailureException
{
public:
@@ -780,6 +802,7 @@ public:
expr_t substituteDiff(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substituteUnaryOpNodes(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
expr_t decreaseLeadsLagsPredeterminedVariables() const override;
expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
bool isNumConstNodeEqualTo(double value) const override;
@@ -792,6 +815,7 @@ public:
bool isInStaticForm() const override;
expr_t replaceVarsInEquation(map<VariableNode *, NumConstNode *> &table) const override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
bool isParamTimesEndogExpr() const override;
expr_t substituteStaticAuxiliaryVariable() const override;
};
@@ -850,6 +874,7 @@ public:
expr_t substituteDiff(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substituteUnaryOpNodes(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
expr_t decreaseLeadsLagsPredeterminedVariables() const override;
expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
bool isNumConstNodeEqualTo(double value) const override;
@@ -862,10 +887,12 @@ public:
bool isInStaticForm() const override;
expr_t replaceVarsInEquation(map<VariableNode *, NumConstNode *> &table) const override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
bool isParamTimesEndogExpr() const override;
//! Substitute auxiliary variables by their expression in static model
expr_t substituteStaticAuxiliaryVariable() const override;
void matchMatchedMoment(vector<int> &symb_ids, vector<int> &lags, vector<int> &powers) const override;
+ pair<int, expr_t> matchEndogenousTimesConstant() const override;
};
//! Unary operator node
@@ -951,6 +978,7 @@ public:
expr_t substituteDiff(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substituteUnaryOpNodes(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
expr_t decreaseLeadsLagsPredeterminedVariables() const override;
expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
bool isNumConstNodeEqualTo(double value) const override;
@@ -963,6 +991,7 @@ public:
bool isInStaticForm() const override;
expr_t replaceVarsInEquation(map<VariableNode *, NumConstNode *> &table) const override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
bool isParamTimesEndogExpr() const override;
//! Substitute auxiliary variables by their expression in static model
expr_t substituteStaticAuxiliaryVariable() const override;
@@ -1056,6 +1085,7 @@ public:
expr_t substituteDiff(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substituteUnaryOpNodes(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
expr_t decreaseLeadsLagsPredeterminedVariables() const override;
expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
bool isNumConstNodeEqualTo(double value) const override;
@@ -1077,6 +1107,7 @@ public:
void findConstantEquations(map<VariableNode *, NumConstNode *> &table) const;
expr_t replaceVarsInEquation(map<VariableNode *, NumConstNode *> &table) const override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
/*
ec_params_and_vars:
- 1st element = feedback force parameter
@@ -1108,6 +1139,7 @@ public:
void decomposeAdditiveTerms(vector<pair<expr_t, int>> &terms, int current_sign) const override;
void decomposeMultiplicativeFactors(vector<pair<expr_t, int>> &factors, int current_exponent = 1) const override;
void matchMatchedMoment(vector<int> &symb_ids, vector<int> &lags, vector<int> &powers) const override;
+ pair<int, expr_t> matchEndogenousTimesConstant() const override;
};
//! Trinary operator node
@@ -1187,6 +1219,7 @@ public:
expr_t substituteDiff(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substituteUnaryOpNodes(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
expr_t decreaseLeadsLagsPredeterminedVariables() const override;
expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
bool isNumConstNodeEqualTo(double value) const override;
@@ -1199,6 +1232,7 @@ public:
bool isInStaticForm() const override;
expr_t replaceVarsInEquation(map<VariableNode *, NumConstNode *> &table) const override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
bool isParamTimesEndogExpr() const override;
//! Substitute auxiliary variables by their expression in static model
expr_t substituteStaticAuxiliaryVariable() const override;
@@ -1294,6 +1328,7 @@ public:
expr_t substituteDiff(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substituteUnaryOpNodes(const lag_equivalence_table_t &nodes, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
virtual expr_t buildSimilarExternalFunctionNode(vector<expr_t> &alt_args, DataTree &alt_datatree) const = 0;
expr_t decreaseLeadsLagsPredeterminedVariables() const override;
expr_t differentiateForwardVars(const vector<string> &subset, subst_table_t &subst_table, vector<BinaryOpNode *> &neweqs) const override;
@@ -1308,6 +1343,7 @@ public:
bool isInStaticForm() const override;
expr_t replaceVarsInEquation(map<VariableNode *, NumConstNode *> &table) const override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
bool isParamTimesEndogExpr() const override;
//! Substitute auxiliary variables by their expression in static model
expr_t substituteStaticAuxiliaryVariable() const override;
@@ -1498,7 +1534,9 @@ public:
int maxLagWithDiffsExpanded() const override;
expr_t substituteVarExpectation(const map<string, expr_t> &subst_table) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
void writeJsonAST(ostream &output) const override;
void writeJsonOutput(ostream &output, const temporary_terms_t &temporary_terms, const deriv_node_temp_terms_t &tef_terms, bool isdynamic) const override;
};
@@ -1512,7 +1550,25 @@ public:
int maxLagWithDiffsExpanded() const override;
expr_t substituteVarExpectation(const map<string, expr_t> &subst_table) const override;
expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
+ bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
+ void writeJsonAST(ostream &output) const override;
+ void writeJsonOutput(ostream &output, const temporary_terms_t &temporary_terms, const deriv_node_temp_terms_t &tef_terms, bool isdynamic) const override;
+};
+
+class PacTargetNonstationaryNode : public SubModelNode
+{
+public:
+ PacTargetNonstationaryNode(DataTree &datatree_arg, int idx_arg, string model_name);
+ void writeOutput(ostream &output, ExprNodeOutputType output_type, const temporary_terms_t &temporary_terms, const temporary_terms_idxs_t &temporary_terms_idxs, const deriv_node_temp_terms_t &tef_terms) const override;
+ expr_t clone(DataTree &datatree) const override;
+ int maxLagWithDiffsExpanded() const override;
+ expr_t substituteVarExpectation(const map<string, expr_t> &subst_table) const override;
+ expr_t substitutePacExpectation(const string &name, expr_t subexpr) override;
+ expr_t substitutePacTargetNonstationary(const string &name, expr_t subexpr) override;
bool containsPacExpectation(const string &pac_model_name = "") const override;
+ bool containsPacTargetNonstationary(const string &pac_model_name = "") const override;
void writeJsonAST(ostream &output) const override;
void writeJsonOutput(ostream &output, const temporary_terms_t &temporary_terms, const deriv_node_temp_terms_t &tef_terms, bool isdynamic) const override;
};
diff --git a/src/ModFile.cc b/src/ModFile.cc
index 902f0563..916d4717 100644
--- a/src/ModFile.cc
+++ b/src/ModFile.cc
@@ -442,7 +442,9 @@ ModFile::transformPass(bool nostrict, bool stochastic, bool compute_xrefs, bool
original_model.fillVarModelTableFromOrigModel();
// PAC model
- pac_model_table.transformPass(diff_subst_table, dynamic_model, var_model_table,
+ pac_model_table.transformPass(unary_ops_nodes, unary_ops_subst_table,
+ diff_nodes, diff_subst_table,
+ dynamic_model, var_model_table,
trend_component_model_table);
// Create auxiliary vars for Expectation operator
@@ -1119,6 +1121,8 @@ ModFile::writeMOutput(const string &basename, bool clear_all, bool clear_global,
// Create epilogue file
epilogue.writeEpilogueFile(basename);
+
+ pac_model_table.writeTargetCoefficientsFile(basename);
}
cout << "done" << endl;
diff --git a/src/ModelTree.cc b/src/ModelTree.cc
index bbbde0c2..3337eb72 100644
--- a/src/ModelTree.cc
+++ b/src/ModelTree.cc
@@ -2047,3 +2047,12 @@ ModelTree::updateReverseVariableEquationOrderings()
eq_idx_orig2block[eq_idx_block2orig[i]] = i;
}
}
+
+expr_t
+ModelTree::getRHSFromLHS(expr_t lhs) const
+{
+ for (auto eq : equations)
+ if (eq->arg1 == lhs)
+ return eq->arg2;
+ throw ExprNode::MatchFailureException{"Cannot find an equation with the requested LHS"};
+}
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index 373273e6..5421792c 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -444,6 +444,10 @@ public:
//! Helper for writing the Jacobian elements in MATLAB and C
/*! Writes either (i+1,j+1) or [i+j*no_eq] */
void jacobianHelper(ostream &output, int eq_nb, int col_nb, ExprNodeOutputType output_type) const;
+ /* Given an expression, searches for the first equation that has exactly this
+ expression on the LHS, and returns the RHS of that equation.
+ If no such equation can be found, throws an ExprNode::MatchFailureExpression */
+ expr_t getRHSFromLHS(expr_t lhs) const;
//! Returns all the equation tags associated to an equation
inline map<string, string>
diff --git a/src/ParsingDriver.cc b/src/ParsingDriver.cc
index 30948634..215148a4 100644
--- a/src/ParsingDriver.cc
+++ b/src/ParsingDriver.cc
@@ -2615,6 +2615,15 @@ ParsingDriver::add_pac_expectation(const string &model_name)
return data_tree->AddPacExpectation(model_name);
}
+expr_t
+ParsingDriver::add_pac_target_nonstationary(const string &model_name)
+{
+ if (data_tree == occbin_constraints_tree.get())
+ error("The 'pac_target_nonstationary' operator is forbidden in 'occbin_constraints'.");
+
+ return data_tree->AddPacTargetNonstationary(model_name);
+}
+
void
ParsingDriver::begin_pac_growth()
{
@@ -3502,6 +3511,57 @@ ParsingDriver::end_occbin_constraints(const vector<tuple<string, BinaryOpNode *,
reset_data_tree();
}
+void
+ParsingDriver::begin_pac_target_info(string name)
+{
+ pac_target_info_name = move(name);
+ set_current_data_tree(&mod_file->dynamic_model);
+}
+
+void
+ParsingDriver::end_pac_target_info()
+{
+ reset_data_tree();
+}
+
+void
+ParsingDriver::set_pac_target_info_target(expr_t target)
+{
+ mod_file->pac_model_table.setTargetExpr(pac_target_info_name, target);
+}
+
+void
+ParsingDriver::set_pac_target_info_auxname_target_nonstationary(string auxname)
+{
+ mod_file->pac_model_table.setTargetAuxnameNonstationary(pac_target_info_name, move(auxname));
+}
+
+void
+ParsingDriver::add_pac_target_info_component(expr_t component_expr)
+{
+ get<0>(pac_target_info_component) = component_expr;
+ mod_file->pac_model_table.addTargetComponent(pac_target_info_name, pac_target_info_component);
+ pac_target_info_component = {};
+}
+
+void
+ParsingDriver::set_pac_target_info_component_growth(expr_t growth)
+{
+ get<1>(pac_target_info_component) = growth;
+}
+
+void
+ParsingDriver::set_pac_target_info_component_auxname(string auxname)
+{
+ get<2>(pac_target_info_component) = move(auxname);
+}
+
+void
+ParsingDriver::set_pac_target_info_component_kind(PacTargetKind kind)
+{
+ get<3>(pac_target_info_component) = kind;
+}
+
vector<string>
ParsingDriver::strsplit(const string &str, char delim)
{
diff --git a/src/ParsingDriver.hh b/src/ParsingDriver.hh
index b6cd328e..8c9d08fb 100644
--- a/src/ParsingDriver.hh
+++ b/src/ParsingDriver.hh
@@ -249,6 +249,9 @@ private:
SymbolList graph_formats;
//! Temporary storage for equation tags
map<string, string> eq_tags;
+ // Temporary storages for pac_target_info
+ string pac_target_info_name;
+ PacModelTable::target_component_t pac_target_info_component;
//! The mod file representation constructed by this ParsingDriver
unique_ptr<ModFile> mod_file;
@@ -741,6 +744,8 @@ public:
expr_t add_var_expectation(const string &model_name);
//! Writes token "PAC_EXPECTATION(model_name, discount, growth)" to model tree
expr_t add_pac_expectation(const string &model_name);
+ //! Adds a pac_target_nonstationary(model_name, discount, growth) node to model tree
+ expr_t add_pac_target_nonstationary(const string &model_name);
//! Creates pac_model statement
void begin_pac_growth();
void begin_pac_model();
@@ -897,6 +902,14 @@ public:
void begin_model_replace(const vector<pair<string, string>> &listed_eqs_by_tags);
// Add a var_remove statement
void var_remove();
+ void begin_pac_target_info(string name);
+ void end_pac_target_info();
+ void set_pac_target_info_target(expr_t target);
+ void set_pac_target_info_auxname_target_nonstationary(string auxname);
+ void add_pac_target_info_component(expr_t component_expr);
+ void set_pac_target_info_component_growth(expr_t growth);
+ void set_pac_target_info_component_auxname(string auxname);
+ void set_pac_target_info_component_kind(PacTargetKind kind);
// Equivalent of MATLAB’s strsplit. Returns an empty vector given an empty string.
static vector<string> strsplit(const string &str, char delim);
// Returns true iff the string is a legal symbol identifier (see NAME token in lexer)
diff --git a/src/SubModel.cc b/src/SubModel.cc
index 28b09016..73d2ce57 100644
--- a/src/SubModel.cc
+++ b/src/SubModel.cc
@@ -18,6 +18,7 @@
*/
#include <algorithm>
+#include <cassert>
#include "SubModel.hh"
#include "DynamicModel.hh"
@@ -791,7 +792,60 @@ PacModelTable::checkPass(ModFileStructure &mod_file_struct, WarningConsolidation
{
for (auto &[name, gv] : growth)
if (gv)
- gv->collectVariables(SymbolType::exogenous, mod_file_struct.pac_params);
+ {
+ if (target_info.find(name) != target_info.end())
+ {
+ cerr << "ERROR: for PAC model '" << name << "', it is not possible to declare a 'growth' option in the 'pac_model' command when there is also a 'pac_target_info' block" << endl;
+ exit(EXIT_FAILURE);
+ }
+ gv->collectVariables(SymbolType::exogenous, mod_file_struct.pac_params);
+ }
+
+ for (const auto &[name, ti] : target_info)
+ for (auto &[expr, gv, auxname, kind, coeff, growth_neutrality_param, h_indices, original_gv, gv_info] : get<2>(ti))
+ if (gv)
+ gv->collectVariables(SymbolType::exogenous, mod_file_struct.pac_params);
+
+ for (const auto &[name, ti] : target_info)
+ {
+ auto &[target, auxname_target_nonstationary, components] = ti;
+ if (!target)
+ {
+ cerr << "ERROR: the block 'pac_target_info(" << name << ")' is missing the 'target' statement" << endl;
+ exit(EXIT_FAILURE);
+ }
+ if (auxname_target_nonstationary.empty())
+ {
+ cerr << "ERROR: the block 'pac_target_info(" << name << ")' is missing the 'auxname_target_nonstationary' statement" << endl;
+ exit(EXIT_FAILURE);
+ }
+ int nonstationary_nb = 0;
+ for (auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : components)
+ {
+ if (auxname.empty())
+ {
+ cerr << "ERROR: the block 'pac_target_info(" << name << ")' is missing the 'auxname' statement in some 'component'" << endl;
+ exit(EXIT_FAILURE);
+ }
+ if (kind == PacTargetKind::unspecified)
+ {
+ cerr << "ERROR: the block 'pac_target_info(" << name << ")' is missing the 'kind' statement in some 'component'" << endl;
+ exit(EXIT_FAILURE);
+ }
+ if (kind == PacTargetKind::ll && growth_component)
+ {
+ cerr << "ERROR: in the block 'pac_target_info(" << name << ")', a component of 'kind ll' (i.e. stationary) has a 'growth' option. This is not permitted." << endl;
+ exit(EXIT_FAILURE);
+ }
+ if (kind == PacTargetKind::dd || kind == PacTargetKind::dl)
+ nonstationary_nb++;
+ }
+ if (!nonstationary_nb)
+ {
+ cerr << "ERROR: the block 'pac_target_info(" << name << ")' must contain at least one nonstationary component (i.e. of 'kind' equal to either 'dd' or 'dl')." << endl;
+ exit(EXIT_FAILURE);
+ }
+ }
}
void
@@ -800,6 +854,11 @@ PacModelTable::findDiffNodesInGrowth(lag_equivalence_table_t &diff_nodes) const
for (auto &[name, gv] : growth)
if (gv)
gv->findDiffNodes(diff_nodes);
+
+ for (const auto &[name, ti] : target_info)
+ for (auto &[expr, gv, auxname, kind, coeff, growth_neutrality_param, h_indices, original_gv, gv_info] : get<2>(ti))
+ if (gv)
+ gv->findDiffNodes(diff_nodes);
}
void
@@ -808,10 +867,18 @@ PacModelTable::substituteDiffNodesInGrowth(const lag_equivalence_table_t &diff_n
for (auto &[name, gv] : growth)
if (gv)
gv = gv->substituteDiff(diff_nodes, diff_subst_table, neweqs);
+
+ for (auto &[name, ti] : target_info)
+ for (auto &[expr, gv, auxname, kind, coeff, growth_neutrality_param, h_indices, original_gv, gv_info] : get<2>(ti))
+ if (gv)
+ gv = gv->substituteDiff(diff_nodes, diff_subst_table, neweqs);
}
void
-PacModelTable::transformPass(ExprNode::subst_table_t &diff_subst_table,
+PacModelTable::transformPass(const lag_equivalence_table_t &unary_ops_nodes,
+ ExprNode::subst_table_t &unary_ops_subst_table,
+ const lag_equivalence_table_t &diff_nodes,
+ ExprNode::subst_table_t &diff_subst_table,
DynamicModel &dynamic_model, const VarModelTable &var_model_table,
const TrendComponentModelTable &trend_component_model_table)
{
@@ -834,6 +901,123 @@ PacModelTable::transformPass(ExprNode::subst_table_t &diff_subst_table,
exit(EXIT_FAILURE);
}
+ // Perform transformations for the pac_target_info block (if any)
+ if (target_info.find(name) != target_info.end())
+ {
+ // Substitute unary ops and diffs in the target…
+ expr_t &target = get<0>(target_info[name]);
+ vector<BinaryOpNode *> neweqs;
+ target = target->substituteUnaryOpNodes(unary_ops_nodes, unary_ops_subst_table, neweqs);
+ if (neweqs.size() > 0)
+ {
+ cerr << "ERROR: the 'target' expression of 'pac_target_info(" << name << ")' contains a variable with a unary operator that is not present in the model" << endl;
+ exit(EXIT_FAILURE);
+ }
+ target = target->substituteDiff(diff_nodes, diff_subst_table, neweqs);
+ if (neweqs.size() > 0)
+ {
+ cerr << "ERROR: the 'target' expression of 'pac_target_info(" << name << ")' contains a diff'd variable that is not present in the model" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ // …and in component expressions
+ auto &components = get<2>(target_info[name]);
+ for (auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : components)
+ {
+ component = component->substituteUnaryOpNodes(unary_ops_nodes, unary_ops_subst_table, neweqs);
+ if (neweqs.size() > 0)
+ {
+ cerr << "ERROR: a 'component' expression of 'pac_target_info(" << name << ")' contains a variable with a unary operator that is not present in the model" << endl;
+ exit(EXIT_FAILURE);
+ }
+ component = component->substituteDiff(diff_nodes, diff_subst_table, neweqs);
+ if (neweqs.size() > 0)
+ {
+ cerr << "ERROR: a 'component' expression of 'pac_target_info(" << name << ")' contains a diff'd variable that is not present in the model" << endl;
+ exit(EXIT_FAILURE);
+ }
+ }
+
+ /* Fill the growth_info structure.
+ Cannot be done in an earlier pass since growth terms can be
+ transformed by DynamicModel::substituteDiff(). */
+ for (auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : components)
+ {
+ if (growth_component)
+ try
+ {
+ growth_component_info = growth_component->matchLinearCombinationOfVariables(false);
+ }
+ catch (ExprNode::MatchFailureException &e)
+ {
+ cerr << "ERROR: PAC growth must be a linear combination of variables" << endl;
+ exit(EXIT_FAILURE);
+ }
+ }
+
+ // Identify the model equation defining the target
+ expr_t target_expr;
+ try
+ {
+ target_expr = dynamic_model.getRHSFromLHS(target);
+ }
+ catch (ExprNode::MatchFailureException)
+ {
+ cerr << "ERROR: there is no equation whose LHS is equal to the 'target' of 'pac_target_info(" << name << ")'" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ // Parse that model equation
+ vector<pair<int, expr_t>> terms;
+ expr_t constant;
+ try
+ {
+ tie(terms, constant) = target_expr->matchLinearCombinationOfEndogenousWithConstant();
+ }
+ catch (ExprNode::MatchFailureException)
+ {
+ cerr << "ERROR: the model equation defining the 'target' of 'pac_target_info(" << name << ")' is not of the right form (should be a linear combination of endogenous variables)" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ // Associate the coefficients of the linear combination with the right components
+ for (auto [var, coeff] : terms)
+ if (auto it = find_if(components.begin(), components.end(),
+ [&](const auto &v) { return get<0>(v) == dynamic_model.AddVariable(var); });
+ it != components.end())
+ get<4>(*it) = coeff;
+ else
+ {
+ cerr << "ERROR: the model equation defining the 'target' of 'pac_target_info(" << name << ")' contains a variable (" << symbol_table.getName(var) << ") that is not declared as a 'component'" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ // Verify that all declared components appear in that equation
+ for (const auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : components)
+ if (!coeff)
+ {
+ cerr << "ERROR: a 'component' of 'pac_target_info(" << name << ")' does not appear in the model equation defining the 'target'" << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ /* Add the variable and equation defining the stationary part of the
+ target. Note that it includes the constant. */
+ expr_t yns = constant;
+ for (const auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : components)
+ if (kind != PacTargetKind::ll)
+ yns = dynamic_model.AddPlus(yns, dynamic_model.AddTimes(coeff, component));
+ int target_nonstationary_id = symbol_table.addPacTargetNonstationaryAuxiliaryVar(get<1>(target_info[name]), yns);
+ expr_t neweq = dynamic_model.AddEqual(dynamic_model.AddVariable(target_nonstationary_id), yns);
+ dynamic_model.addEquation(neweq, -1);
+ dynamic_model.addAuxEquation(neweq);
+
+ /* Perform the substitution of the pac_target_nonstationary operator.
+ This needs to be done here, otherwise
+ DynamicModel::analyzePacEquationStructure() will not be able to
+ identify the error-correction part */
+ dynamic_model.substitutePacTargetNonstationary(name, dynamic_model.AddVariable(target_nonstationary_id, -1));
+ }
+
// Collect some information about PAC models
int max_lag;
if (trend_component_model_table.isExistingTrendComponentModelName(aux_model_name[name]))
@@ -878,31 +1062,106 @@ PacModelTable::transformPass(ExprNode::subst_table_t &diff_subst_table,
if (growth[name])
growth_correction_term = dynamic_model.AddTimes(growth[name], dynamic_model.AddVariable(growth_neutrality_params[name]));
if (aux_model_name[name].empty())
- dynamic_model.computePacModelConsistentExpectationSubstitution(name,
- symbol_table.getID(discount[name]),
- pacEquationMaxLag(name),
- growth_correction_term,
- diff_subst_table,
- aux_var_symb_ids,
- aux_param_symb_ids,
- pac_expectation_substitution);
+ {
+ if (target_info.find(name) != target_info.end())
+ {
+ cerr << "ERROR: the block 'pac_target_info(" << name << ")' is not supported in the context of a PAC model with model-consistent expectations (MCE)." << endl;
+ exit(EXIT_FAILURE);
+ }
+ else
+ dynamic_model.computePacModelConsistentExpectationSubstitution(name,
+ symbol_table.getID(discount[name]),
+ pacEquationMaxLag(name),
+ growth_correction_term,
+ diff_subst_table,
+ aux_var_symb_ids,
+ aux_param_symb_ids,
+ pac_expectation_substitution);
+ }
else
- dynamic_model.computePacBackwardExpectationSubstitution(name, lhs[name], max_lag,
- aux_model_type[name],
- growth_correction_term,
- aux_var_symb_ids,
- aux_param_symb_ids,
- pac_expectation_substitution);
+ {
+ if (target_info.find(name) != target_info.end())
+ {
+ assert(growth_correction_term == dynamic_model.Zero);
+ dynamic_model.computePacBackwardExpectationSubstitutionWithComponents(name, lhs[name],
+ max_lag,
+ aux_model_type[name],
+ get<2>(target_info[name]),
+ pac_expectation_substitution);
+ }
+ else
+ dynamic_model.computePacBackwardExpectationSubstitution(name, lhs[name], max_lag,
+ aux_model_type[name],
+ growth_correction_term,
+ aux_var_symb_ids,
+ aux_param_symb_ids,
+ pac_expectation_substitution);
+ }
}
// Actually perform the substitution of pac_expectation
dynamic_model.substitutePacExpectation(pac_expectation_substitution, eq_name);
dynamic_model.checkNoRemainingPacExpectation();
+
+ // Check that there is no remaining pac_target_nonstationary operator
+ dynamic_model.checkNoRemainingPacTargetNonstationary();
}
void
PacModelTable::writeOutput(const string &basename, ostream &output) const
{
+ // Helper to print the “growth_info” structure (linear decomposition of growth)
+ auto growth_info_helper = [&](const string &fieldname, const growth_info_t &gi)
+ {
+ int i = 1;
+ for (auto [growth_symb_id, growth_lag, param_id, constant] : gi)
+ {
+ string structname = fieldname + "(" + to_string(i++) + ").";
+ if (growth_symb_id >= 0)
+ {
+ string var_field = "endo_id";
+ if (symbol_table.getType(growth_symb_id) == SymbolType::exogenous)
+ {
+ var_field = "exo_id";
+ output << structname << "endo_id = 0;" << endl;
+ }
+ else
+ output << structname << "exo_id = 0;" << endl;
+ try
+ {
+ // case when this is not the highest lag of the growth variable
+ int aux_symb_id = symbol_table.searchAuxiliaryVars(growth_symb_id, growth_lag);
+ output << structname << var_field << " = " << symbol_table.getTypeSpecificID(aux_symb_id) + 1 << ";" << endl
+ << structname << "lag = 0;" << endl;
+ }
+ catch (...)
+ {
+ try
+ {
+ // case when this is the highest lag of the growth variable
+ int tmp_growth_lag = growth_lag + 1;
+ int aux_symb_id = symbol_table.searchAuxiliaryVars(growth_symb_id, tmp_growth_lag);
+ output << structname << var_field << " = " << symbol_table.getTypeSpecificID(aux_symb_id) + 1 << ";" << endl
+ << structname << "lag = -1;" << endl;
+ }
+ catch (...)
+ {
+ // case when there is no aux var for the variable
+ output << structname << var_field << " = "<< symbol_table.getTypeSpecificID(growth_symb_id) + 1 << ";" << endl
+ << structname << "lag = " << growth_lag << ";" << endl;
+ }
+ }
+ }
+ else
+ output << structname << "endo_id = 0;" << endl
+ << structname << "exo_id = 0;" << endl
+ << structname << "lag = 0;" << endl;
+ output << structname << "param_id = "
+ << (param_id == -1 ? 0 : symbol_table.getTypeSpecificID(param_id) + 1) << ";" << endl
+ << structname << "constant = " << constant << ";" << endl;
+ }
+ };
+
for (const auto &name : names)
{
output << "M_.pac." << name << ".auxiliary_model_name = '" << aux_model_name.at(name) << "';" << endl
@@ -913,53 +1172,7 @@ PacModelTable::writeOutput(const string &basename, ostream &output) const
output << "M_.pac." << name << ".growth_str = '";
original_growth.at(name)->writeJsonOutput(output, {}, {}, true);
output << "';" << endl;
- int i = 0;
- for (auto [growth_symb_id, growth_lag, param_id, constant] : growth_info.at(name))
- {
- string structname = "M_.pac." + name + ".growth_linear_comb(" + to_string(++i) + ").";
- if (growth_symb_id >= 0)
- {
- string var_field = "endo_id";
- if (symbol_table.getType(growth_symb_id) == SymbolType::exogenous)
- {
- var_field = "exo_id";
- output << structname << "endo_id = 0;" << endl;
- }
- else
- output << structname << "exo_id = 0;" << endl;
- try
- {
- // case when this is not the highest lag of the growth variable
- int aux_symb_id = symbol_table.searchAuxiliaryVars(growth_symb_id, growth_lag);
- output << structname << var_field << " = " << symbol_table.getTypeSpecificID(aux_symb_id) + 1 << ";" << endl
- << structname << "lag = 0;" << endl;
- }
- catch (...)
- {
- try
- {
- // case when this is the highest lag of the growth variable
- int tmp_growth_lag = growth_lag + 1;
- int aux_symb_id = symbol_table.searchAuxiliaryVars(growth_symb_id, tmp_growth_lag);
- output << structname << var_field << " = " << symbol_table.getTypeSpecificID(aux_symb_id) + 1 << ";" << endl
- << structname << "lag = -1;" << endl;
- }
- catch (...)
- {
- // case when there is no aux var for the variable
- output << structname << var_field << " = "<< symbol_table.getTypeSpecificID(growth_symb_id) + 1 << ";" << endl
- << structname << "lag = " << growth_lag << ";" << endl;
- }
- }
- }
- else
- output << structname << "endo_id = 0;" << endl
- << structname << "exo_id = 0;" << endl
- << structname << "lag = 0;" << endl;
- output << structname << "param_id = "
- << (param_id == -1 ? 0 : symbol_table.getTypeSpecificID(param_id) + 1) << ";" << endl
- << structname << "constant = " << constant << ";" << endl;
- }
+ growth_info_helper("M_.pac." + name + ".growth_linear_comb", growth_info.at(name));
}
}
@@ -1160,6 +1373,34 @@ PacModelTable::writeOutput(const string &basename, ostream &output) const
output << "];" << endl;
}
}
+
+ for (auto &[name, val] : target_info)
+ {
+ int component_idx = 1;
+ for (auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : get<2>(val))
+ {
+ string fieldname = "M_.pac." + name + ".components(" + to_string(component_idx) + ")";
+ output << fieldname << ".aux_id = " << symbol_table.getTypeSpecificID(auxname) + 1 << ";" << endl
+ << fieldname << ".endo_var = " << symbol_table.getTypeSpecificID(dynamic_cast<VariableNode *>(component)->symb_id) + 1 << ";" << endl
+ << fieldname << ".kind = '" << kindToString(kind) << "';" << endl
+ << fieldname << ".h_param_indices = [";
+ for (int id : h_indices)
+ output << symbol_table.getTypeSpecificID(id) + 1 << " ";
+ output << "];" << endl
+ << fieldname << ".coeff_str = '";
+ coeff->writeJsonOutput(output, {}, {}, true);
+ output << "';" << endl;
+ if (growth_component)
+ {
+ output << fieldname << ".growth_neutrality_param_index = " << symbol_table.getTypeSpecificID(growth_neutrality_param) + 1 << ";" << endl
+ << fieldname << ".growth_str = '";
+ original_growth_component->writeJsonOutput(output, {}, {}, true);
+ output << "';" << endl;
+ growth_info_helper(fieldname + ".growth_linear_comb", growth_component_info);
+ }
+ component_idx++;
+ }
+ }
}
void
@@ -1167,6 +1408,8 @@ PacModelTable::writeJsonOutput(ostream &output) const
{
for (const auto &name : names)
{
+ /* The calling method has already added a comma, so don’t output one for
+ the first statement */
if (name != *names.begin())
output << ", ";
output << R"({"statementName": "pac_model",)"
@@ -1181,6 +1424,31 @@ PacModelTable::writeJsonOutput(ostream &output) const
}
output << "}" << endl;
}
+
+ for (auto &[name, val] : target_info)
+ {
+ output << R"(, {"statementName": "pac_target_info", "model_name": ")" << name
+ << R"(", "target": ")";
+ get<0>(val)->writeJsonOutput(output, {}, {}, true);
+ output << R"(", "auxname_target_nonstationary": ")" << get<1>(val)
+ << R"(", "components": [)";
+ for (auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : get<2>(val))
+ {
+ if (component != get<0>(get<2>(val).front()))
+ output << ", ";
+ output << R"({"component": ")";
+ component->writeJsonOutput(output, {}, {}, true);
+ output << R"(", "auxname": ")" << auxname
+ << R"(", "kind": ")" << kindToString(kind);
+ if (growth_component)
+ {
+ output << R"(", "growth_str": ")";
+ original_growth_component->writeJsonOutput(output, {}, {}, true);
+ }
+ output << R"("})";
+ }
+ output << "]}" << endl;
+ }
}
int
@@ -1188,3 +1456,75 @@ PacModelTable::pacEquationMaxLag(const string &name_arg) const
{
return get<2>(equation_info.at(name_arg)).size();
}
+
+string
+PacModelTable::kindToString(PacTargetKind kind)
+{
+ switch (kind)
+ {
+ case PacTargetKind::unspecified:
+ cerr << "Internal error: kind should not be unspecified" << endl;
+ exit(EXIT_FAILURE);
+ case PacTargetKind::ll:
+ return "ll";
+ case PacTargetKind::dl:
+ return "dl";
+ case PacTargetKind::dd:
+ return "dd";
+ }
+ // Silent GCC warning
+ assert(false);
+}
+
+void
+PacModelTable::setTargetExpr(const string &name_arg, expr_t target)
+{
+ get<0>(target_info[name_arg]) = target;
+}
+
+void
+PacModelTable::setTargetAuxnameNonstationary(const string &name_arg, string auxname)
+{
+ get<1>(target_info[name_arg]) = move(auxname);
+}
+
+void
+PacModelTable::addTargetComponent(const string &name_arg, target_component_t component)
+{
+ get<7>(component) = get<1>(component); // original_growth = growth
+ get<2>(target_info[name_arg]).emplace_back(move(component));
+}
+
+void
+PacModelTable::writeTargetCoefficientsFile(const string &basename) const
+{
+ if (target_info.empty())
+ return;
+
+ string filename = DataTree::packageDir(basename) + "/pac_target_coefficients.m";
+ ofstream output;
+ output.open(filename, ios::out | ios::binary);
+ if (!output.is_open())
+ {
+ cerr << "ERROR: Can't open file " << filename << " for writing" << endl;
+ exit(EXIT_FAILURE);
+ }
+ output << "function coeffs = pac_target_coefficients(model_name, params)" << endl;
+ for (auto &[model_name, val] : target_info)
+ {
+ output << " if strcmp(model_name, '" << model_name << "')" << endl
+ << " coeffs = NaN(" << get<2>(val).size() << ",1);" << endl;
+ int i = 1;
+ for (auto &[component, growth_component, auxname, kind, coeff, growth_neutrality_param, h_indices, original_growth_component, growth_component_info] : get<2>(val))
+ {
+ output << " coeffs(" << i++ << ") = ";
+ coeff->writeOutput(output, ExprNodeOutputType::matlabDynamicModel);
+ output << ";" << endl;
+ }
+ output << " return" << endl
+ << " end" << endl;
+ }
+ output << " error([ 'Unknown PAC model: ' model_name ])" << endl
+ << "end" << endl;
+ output.close();
+}
diff --git a/src/SubModel.hh b/src/SubModel.hh
index d45f29c7..74d924c7 100644
--- a/src/SubModel.hh
+++ b/src/SubModel.hh
@@ -197,12 +197,15 @@ private:
set<string> names;
map<string, string> aux_model_name;
map<string, string> discount;
- // The growth expressions belong to the main dynamic_model from the ModFile instance
+ /* The growth expressions belong to the main dynamic_model from the ModFile
+ instance. The growth expression is necessarily nullptr for a model with a
+ pac_target_info block. */
map<string, expr_t> growth, original_growth;
/* Information about the structure of growth expressions (which must be a
linear combination of variables).
Each tuple represents a term: (endo_id, lag, param_id, constant) */
- map<string, vector<tuple<int, int, int, double>>> growth_info;
+ using growth_info_t = vector<tuple<int, int, int, double>>;
+ map<string, growth_info_t> growth_info;
/* Stores the name of the PAC equation associated to the model.
pac_model_name → eq_name */
@@ -214,6 +217,8 @@ private:
pac_expectation value
- in the MCE case, this auxiliary represents Z₁ (i.e. without the growth
correction term)
+ Note that this structure is not used in the presence of the
+ pac_target_info block.
pac_model_name → symb_id */
map<string, int> aux_var_symb_ids;
/* Stores symb_ids for auxiliary parameters created for the expression
@@ -221,10 +226,13 @@ private:
neutrality correction):
- in the backward case, contains the “h” parameters
- in the MCE case, contains the “α” parameters
+ Note that this structure is not used in the presence of the
+ pac_target_info block.
pac_model_name → symb_ids */
map<string, vector<int>> aux_param_symb_ids;
/* Stores indices for growth neutrality parameters
- pac_model_name → growth_neutrality_param_index */
+ pac_model_name → growth_neutrality_param_index.
+ This map is not used for PAC models with a pac_target_info block. */
map<string, int> growth_neutrality_params;
// Stores LHS vars (only for backward PAC models)
@@ -243,8 +251,21 @@ public:
private:
equation_info_t equation_info;
+public:
+ /* (component variable/expr, growth, auxname, kind, coeff. in the linear
+ combination, growth_param ID possibly equal to -1, vector of h parameters,
+ original_growth, growth_info) */
+ using target_component_t = tuple<expr_t, expr_t, string, PacTargetKind, expr_t, int, vector<int>, expr_t, growth_info_t>;
+
+private:
+ // pac_model_name → (target variable/expr, auxname_target_nonstationary, target components)
+ map<string, tuple<expr_t, string, vector<target_component_t>>> target_info;
+
int pacEquationMaxLag(const string &name_arg) const;
+ // Return a text representation of a kind (but fails on “unspecified” kind value)
+ static string kindToString(PacTargetKind kind);
+
public:
explicit PacModelTable(SymbolTable &symbol_table_arg);
void addPacModel(string name_arg, string aux_model_name_arg, string discount_arg, expr_t growth_arg);
@@ -255,11 +276,20 @@ public:
// Called by DynamicModel::substituteDiff()
void substituteDiffNodesInGrowth(const lag_equivalence_table_t &diff_nodes, ExprNode::subst_table_t &diff_subst_table, vector<BinaryOpNode *> &neweqs);
// Must be called after substituteDiffNodesInGrowth()
- void transformPass(ExprNode::subst_table_t &diff_subst_table,
+ void transformPass(const lag_equivalence_table_t &unary_ops_nodes,
+ ExprNode::subst_table_t &unary_ops_subst_table,
+ const lag_equivalence_table_t &diff_nodes,
+ ExprNode::subst_table_t &diff_subst_table,
DynamicModel &dynamic_model, const VarModelTable &var_model_table,
const TrendComponentModelTable &trend_component_model_table);
void writeOutput(const string &basename, ostream &output) const;
void writeJsonOutput(ostream &output) const;
+ void setTargetExpr(const string &name_arg, expr_t target);
+ void setTargetAuxnameNonstationary(const string &name_arg, string auxname);
+ /* Only the first four elements of the tuple are expected to be set by the
+ caller. The other ones will be filled by this class. */
+ void addTargetComponent(const string &name_arg, target_component_t component);
+ void writeTargetCoefficientsFile(const string &basename) const;
};
diff --git a/src/SymbolTable.cc b/src/SymbolTable.cc
index 23c028ea..32ea5670 100644
--- a/src/SymbolTable.cc
+++ b/src/SymbolTable.cc
@@ -377,6 +377,7 @@ SymbolTable::writeOutput(ostream &output) const noexcept(false)
break;
case AuxVarType::expectation:
case AuxVarType::pacExpectation:
+ case AuxVarType::pacTargetNonstationary:
break;
case AuxVarType::diff:
case AuxVarType::diffLag:
@@ -687,6 +688,24 @@ SymbolTable::addPacExpectationAuxiliaryVar(const string &name, expr_t expr_arg)
return symb_id;
}
+int
+SymbolTable::addPacTargetNonstationaryAuxiliaryVar(const string &name, expr_t expr_arg)
+{
+ int symb_id;
+ try
+ {
+ symb_id = addSymbol(name, SymbolType::endogenous);
+ }
+ catch (AlreadyDeclaredException &e)
+ {
+ cerr << "ERROR: the variable/parameter '" << name << "' conflicts with a variable that will be generated for a 'pac_target_nonstationary' expression. Please rename it." << endl;
+ exit(EXIT_FAILURE);
+ }
+
+ aux_vars.emplace_back(symb_id, AuxVarType::pacTargetNonstationary, 0, 0, 0, 0, expr_arg, "");
+ return symb_id;
+}
+
int
SymbolTable::searchAuxiliaryVars(int orig_symb_id, int orig_lead_lag) const noexcept(false)
{
@@ -1004,6 +1023,7 @@ SymbolTable::writeJsonOutput(ostream &output) const
break;
case AuxVarType::expectation:
case AuxVarType::pacExpectation:
+ case AuxVarType::pacTargetNonstationary:
break;
case AuxVarType::diff:
case AuxVarType::diffLag:
diff --git a/src/SymbolTable.hh b/src/SymbolTable.hh
index b3fd7be4..9d3579ac 100644
--- a/src/SymbolTable.hh
+++ b/src/SymbolTable.hh
@@ -53,7 +53,8 @@ enum class AuxVarType
diffLag = 9, //!< Variable for timing between Diff operators (lag)
unaryOp = 10, //!< Variable for allowing the undiff operator to work when diff was taken of unary op, eg diff(log(x))
diffLead = 11, //!< Variable for timing between Diff operators (lead)
- pacExpectation = 12 //!< Variable created for the substitution of the pac_expectation operator
+ pacExpectation = 12, //!< Variable created for the substitution of the pac_expectation operator
+ pacTargetNonstationary = 13 //!< Variable created for the substitution of the pac_target_nonstationary operator
};
//! Information on some auxiliary variables
@@ -337,6 +338,8 @@ public:
int addUnaryOpAuxiliaryVar(int index, expr_t expr_arg, string unary_op, int orig_symb_id = -1, int orig_lag = 0) noexcept(false);
//! An auxiliary variable for a pac_expectation operator
int addPacExpectationAuxiliaryVar(const string &name, expr_t expr_arg);
+ //! An auxiliary variable for a pac_target_nonstationary operator
+ int addPacTargetNonstationaryAuxiliaryVar(const string &name, expr_t expr_arg);
//! Returns the number of auxiliary variables
int
AuxVarsSize() const
--
GitLab